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