test_approx.c

   1 #include <assert.h>
   2 #include <ctype.h>
   3 #include <limits.h>
   4 #include <math.h>
   5 #include <stdlib.h>
   6 #include <string.h>
   7 #include <isl/ctx.h>
   8 #include <isl/val.h>
   9 #include <isl/space.h>
  10 #include <isl/set.h>
  11 #include <isl/polynomial.h>
  12 #include <barvinok/polylib.h>
  13 #include <barvinok/isl.h>
  14 #include <barvinok/options.h>
  15 #include <barvinok/barvinok.h>
  16 #include "verify.h"
  17 #include "config.h"
  18
  19 #ifdef HAVE_SYS_TIMES_H
  20
  21 #include <sys/times.h>
  22
  23 typedef clock_t         my_clock_t;
  24
  25 static my_clock_t time_diff(struct tms *before, struct tms *after)
  26 {
  27         return after->tms_utime - before->tms_utime;
  28 }
  29
  30 #else
  31
  32 typedef int             my_clock_t;
  33
  34 struct tms { int dummy; };
  35 static void times(struct tms* time)
  36 {
  37 }
  38 static my_clock_t time_diff(struct tms *before, struct tms *after)
  39 {
  40         return 0;
  41 }
  42
  43 #endif
  44
  45 struct {
  46     int     sign;
  47     int     method;
  48     int     flags;
  49 } methods[] = {
  50 { BV_APPROX_SIGN_NONE,      BV_APPROX_NONE,     0 },
  51 { BV_APPROX_SIGN_APPROX,    BV_APPROX_BERNOULLI,        0 },
  52 { BV_APPROX_SIGN_APPROX,    BV_APPROX_DROP,     0 },
  53 { BV_APPROX_SIGN_APPROX,    BV_APPROX_VOLUME,   BV_VOL_LIFT },
  54 { BV_APPROX_SIGN_APPROX,    BV_APPROX_VOLUME,   BV_VOL_VERTEX },
  55 { BV_APPROX_SIGN_APPROX,    BV_APPROX_VOLUME,   BV_VOL_BARYCENTER },
  56 { BV_APPROX_SIGN_APPROX,    BV_APPROX_SCALE,    0 },
  57 { BV_APPROX_SIGN_APPROX,    BV_APPROX_SCALE,    BV_APPROX_SCALE_CHAMBER },
  58 { BV_APPROX_SIGN_APPROX,    BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST },
  59 { BV_APPROX_SIGN_APPROX,    BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
  60 { BV_APPROX_SIGN_LOWER,     BV_APPROX_DROP,     0 },
  61 { BV_APPROX_SIGN_LOWER,     BV_APPROX_VOLUME,   BV_VOL_LIFT },
  62 { BV_APPROX_SIGN_LOWER,     BV_APPROX_VOLUME,   BV_VOL_VERTEX },
  63 { BV_APPROX_SIGN_LOWER,     BV_APPROX_VOLUME,   BV_VOL_BARYCENTER },
  64 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    0 },
  65 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_CHAMBER },
  66 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST },
  67 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
  68 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW },
  69 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW  | BV_APPROX_SCALE_CHAMBER},
  70 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW },
  71 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
  72 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW2 },
  73 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
  74 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 },
  75 { BV_APPROX_SIGN_LOWER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
  76 { BV_APPROX_SIGN_UPPER,     BV_APPROX_DROP,     0 },
  77 { BV_APPROX_SIGN_UPPER,     BV_APPROX_VOLUME,   BV_VOL_LIFT },
  78 { BV_APPROX_SIGN_UPPER,     BV_APPROX_VOLUME,   BV_VOL_VERTEX },
  79 { BV_APPROX_SIGN_UPPER,     BV_APPROX_VOLUME,   BV_VOL_BARYCENTER },
  80 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    0 },
  81 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_CHAMBER },
  82 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST },
  83 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_CHAMBER },
  84 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW },
  85 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
  86 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW },
  87 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW | BV_APPROX_SCALE_CHAMBER },
  88 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW2 },
  89 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
  90 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 },
  91 { BV_APPROX_SIGN_UPPER,     BV_APPROX_SCALE,    BV_APPROX_SCALE_FAST | BV_APPROX_SCALE_NARROW2 | BV_APPROX_SCALE_CHAMBER },
  92 };
  93
  94 #define nr_methods (sizeof(methods)/sizeof(*methods))
  95
  96 struct options {
  97     int quiet;
  98     struct verify_options    *verify;
  99 };
 100
 101 ISL_ARGS_START(struct options, options_args)
 102 ISL_ARG_CHILD(struct options, verify, NULL, &verify_options_args, NULL)
 103 ISL_ARG_BOOL(struct options, quiet, 'q', "quiet", 0, NULL)
 104 ISL_ARGS_END
 105
 106 ISL_ARG_DEF(options, struct options, options_args)
 107
 108 struct result_data {
 109     Value                   n;
 110     double                  RE_sum[nr_methods];
 111
 112     my_clock_t              ticks[nr_methods];
 113     size_t                  size[nr_methods];
 114 };
 115
 116 void result_data_init(struct result_data *result)
 117 {
 118     int i;
 119     for (i = 0; i < nr_methods; ++i) {
 120         result->RE_sum[i] = 0;
 121         result->ticks[i] = 0;
 122         result->size[i] = 0;
 123     }
 124     value_init(result->n);
 125 }
 126
 127 void result_data_clear(struct result_data *result)
 128 {
 129     value_clear(result->n);
 130 }
 131
 132 void result_data_print(struct result_data *result, int n)
 133 {
 134     int i;
 135
 136     fprintf(stderr, "%g", (double)result->ticks[0]/n);
 137     for (i = 1; i < nr_methods; ++i)
 138         fprintf(stderr, ", %g", (double)result->ticks[i]/n);
 139     fprintf(stderr, "\n");
 140
 141     fprintf(stderr, "%zd", result->size[0]/n);
 142     for (i = 1; i < nr_methods; ++i)
 143         fprintf(stderr, ", %zd", result->size[i]/n);
 144     fprintf(stderr, "\n");
 145
 146     fprintf(stderr, "%g\n", VALUE_TO_DOUBLE(result->n));
 147     fprintf(stderr, "%g", result->RE_sum[0]/VALUE_TO_DOUBLE(result->n));
 148     for (i = 1; i < nr_methods; ++i)
 149         fprintf(stderr, ", %g", result->RE_sum[i]/VALUE_TO_DOUBLE(result->n));
 150     fprintf(stderr, "\n");
 151 }
 152
 153 struct test_approx_data {
 154     struct verify_point_data vpd;
 155     isl_pw_qpolynomial      **pwqp;
 156     struct result_data      *result;
 157 };
 158
 159 static __isl_give isl_val *eval(__isl_keep isl_pw_qpolynomial *pwqp,
 160         __isl_keep isl_point *pnt, int sign)
 161 {
 162     isl_val *res;
 163
 164     pwqp = isl_pw_qpolynomial_copy(pwqp);
 165     res = isl_pw_qpolynomial_eval(pwqp, isl_point_copy(pnt));
 166     if (sign == BV_APPROX_SIGN_LOWER)
 167         res = isl_val_ceil(res);
 168     else if (sign == BV_APPROX_SIGN_UPPER)
 169         res = isl_val_floor(res);
 170     else if (sign == BV_APPROX_SIGN_APPROX)
 171         res = isl_val_trunc(res);
 172     else
 173         assert(isl_val_is_int(res));
 174
 175     return res;
 176 }
 177
 178 static isl_stat test_approx(__isl_take isl_point *pnt, void *user)
 179 {
 180     struct test_approx_data *ta_data = (struct test_approx_data *) user;
 181     isl_val *exact, *approx;
 182     int i;
 183
 184     ta_data->vpd.n--;
 185
 186     exact = eval(ta_data->pwqp[0], pnt, BV_APPROX_SIGN_NONE);
 187
 188     value_increment(ta_data->result->n, ta_data->result->n);
 189     for (i = 1; i < nr_methods; ++i) {
 190         double error;
 191         approx = eval(ta_data->pwqp[i], pnt, methods[i].sign);
 192         assert(approx);
 193         if (methods[i].sign == BV_APPROX_SIGN_LOWER)
 194             assert(isl_val_le(approx, exact));
 195         if (methods[i].sign == BV_APPROX_SIGN_UPPER)
 196             assert(isl_val_ge(approx, exact));
 197         approx = isl_val_sub(approx, isl_val_copy(exact));
 198         if (isl_val_is_zero(exact))
 199             error = fabs(isl_val_get_d(approx));
 200         else
 201             error = fabs(isl_val_get_d(approx)) / isl_val_get_d(exact);
 202         isl_val_free(approx);
 203         ta_data->result->RE_sum[i] += error;
 204     }
 205
 206     if (!ta_data->vpd.options->print_all &&
 207         (ta_data->vpd.n % ta_data->vpd.s) == 0) {
 208             printf("o");
 209             fflush(stdout);
 210     }
 211
 212     isl_val_free(exact);
 213     isl_point_free(pnt);
 214
 215     return (ta_data->vpd.n >= 1) ? isl_stat_ok : isl_stat_error;
 216 }
 217
 218 static int test(__isl_keep isl_set *context,
 219         __isl_keep isl_pw_qpolynomial **pwqp, struct result_data *result,
 220         struct verify_options *options)
 221 {
 222         struct test_approx_data data = { { options } };
 223         int r;
 224
 225         r = verify_point_data_init(&data.vpd, context);
 226
 227         data.pwqp = pwqp;
 228         data.result = result;
 229         if (r == 0)
 230                 isl_set_foreach_point(context, &test_approx, &data);
 231         if (data.vpd.error)
 232                 r = -1;
 233
 234         verify_point_data_fini(&data.vpd);
 235
 236         return r;
 237 }
 238
 239 /* Turn the set dimensions of "context" into parameters and return
 240  * the corresponding parameter domain.
 241  */
 242 static __isl_give isl_set *to_parameter_domain(__isl_take isl_set *context)
 243 {
 244         context = isl_set_move_dims(context, isl_dim_param, 0, isl_dim_set, 0,
 245                                        isl_set_dim(context, isl_dim_set));
 246         return isl_set_params(context);
 247 }
 248
 249 static int handle(isl_ctx *ctx, FILE *in, struct result_data *result,
 250         struct verify_options *options)
 251 {
 252     int i;
 253     int r;
 254     int nparam;
 255     isl_space *space;
 256     isl_set *set;
 257     isl_set *context;
 258     isl_pw_qpolynomial *pwqp[nr_methods];
 259
 260     set = isl_set_read_from_file(ctx, in);
 261     context = isl_set_read_from_file(ctx, in);
 262
 263     context = to_parameter_domain(context);
 264     nparam = isl_set_dim(context, isl_dim_param);
 265     if (nparam != isl_set_dim(set, isl_dim_param)) {
 266             int dim = isl_set_dim(set, isl_dim_set);
 267             set = isl_set_move_dims(set, isl_dim_param, 0,
 268                                         isl_dim_set, dim - nparam, nparam);
 269     }
 270
 271     set = isl_set_intersect_params(set, context);
 272     context = isl_set_params(isl_set_copy(set));
 273     space = isl_set_get_space(context);
 274
 275     context = verify_context_set_bounds(context, options);
 276
 277     for (i = 0; i < nr_methods; ++i) {
 278         struct tms st_cpu;
 279         struct tms en_cpu;
 280         options->barvinok->approx->approximation = methods[i].sign;
 281         options->barvinok->approx->method = methods[i].method;
 282         if (methods[i].method == BV_APPROX_SCALE)
 283             options->barvinok->approx->scale_flags = methods[i].flags;
 284         else if (methods[i].method == BV_APPROX_VOLUME)
 285             options->barvinok->approx->volume_triangulate = methods[i].flags;
 286
 287         times(&st_cpu);
 288         pwqp[i] = isl_set_card(isl_set_copy(set));
 289         times(&en_cpu);
 290         result->ticks[i] = time_diff(&en_cpu, &st_cpu);
 291     }
 292     for (i = 0; i < nr_methods; ++i)
 293         result->size[i] = 0;
 294     r = test(context, pwqp, result, options);
 295     for (i = 0; i < nr_methods; ++i)
 296         isl_pw_qpolynomial_free(pwqp[i]);
 297
 298     isl_space_free(space);
 299     isl_set_free(context);
 300     isl_set_free(set);
 301
 302     return r;
 303 }
 304
 305 int main(int argc, char **argv)
 306 {
 307     isl_ctx *ctx;
 308     char path[PATH_MAX+1];
 309     struct result_data all_result;
 310     int n = 0;
 311     int r = EXIT_SUCCESS;
 312     struct options *options = options_new_with_defaults();
 313
 314     argc = options_parse(options, argc, argv, ISL_ARG_ALL);
 315     ctx = isl_ctx_alloc_with_options(&options_args, options);
 316
 317     if (options->verify->M == INT_MIN)
 318         options->verify->M = 100;
 319     if (options->verify->m == INT_MAX)
 320         options->verify->m = -100;
 321
 322     result_data_init(&all_result);
 323
 324     while (fgets(path, sizeof(path), stdin)) {
 325         struct result_data result;
 326         FILE *in;
 327         int i;
 328
 329         ++n;
 330         result_data_init(&result);
 331         fprintf(stderr, "%s", path);
 332         *strchr(path, '\n') = '\0';
 333         in = fopen(path, "r");
 334         assert(in);
 335         if (handle(ctx, in, &result, options->verify) < 0)
 336             r = EXIT_FAILURE;
 337         fclose(in);
 338
 339         if (!options->quiet)
 340             result_data_print(&result, 1);
 341
 342         value_addto(all_result.n, all_result.n, result.n);
 343         for (i = 0; i < nr_methods; ++i) {
 344             all_result.RE_sum[i] += result.RE_sum[i];
 345             all_result.ticks[i] += result.ticks[i];
 346             all_result.size[i] += result.size[i];
 347         }
 348
 349         result_data_clear(&result);
 350
 351         if (!options->quiet) {
 352             fprintf(stderr, "average\n");
 353             result_data_print(&all_result, n);
 354         }
 355     }
 356
 357     result_data_clear(&all_result);
 358
 359     isl_ctx_free(ctx);
 360
 361     return r;
 362 }