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         assert(approx);
 185         if (methods[i].sign == BV_APPROX_SIGN_LOWER)
 186             assert(isl_val_le(approx, exact));
 187         if (methods[i].sign == BV_APPROX_SIGN_UPPER)
 188             assert(isl_val_ge(approx, exact));
 189         approx = isl_val_sub(approx, isl_val_copy(exact));
 190         if (isl_val_is_zero(exact))
 191             error = fabs(isl_val_get_d(approx));
 192         else
 193             error = fabs(isl_val_get_d(approx)) / isl_val_get_d(exact);
 194         isl_val_free(approx);
 195         ta_data->result->RE_sum[i] += error;
 196     }
 197
 198     if (!ta_data->vpd.options->print_all &&
 199         (ta_data->vpd.n % ta_data->vpd.s) == 0) {
 200             printf("o");
 201             fflush(stdout);
 202     }
 203
 204     isl_val_free(exact);
 205     isl_point_free(pnt);
 206
 207     return (ta_data->vpd.n >= 1) ? isl_stat_ok : isl_stat_error;
 208 }
 209
 210 static int test(__isl_keep isl_set *context,
 211         __isl_keep isl_pw_qpolynomial **pwqp, struct result_data *result,
 212         struct verify_options *options)
 213 {
 214         struct test_approx_data data = { { options } };
 215         int r;
 216
 217         r = verify_point_data_init(&data.vpd, context);
 218
 219         data.pwqp = pwqp;
 220         data.result = result;
 221         if (r == 0)
 222                 isl_set_foreach_point(context, &test_approx, &data);
 223         if (data.vpd.error)
 224                 r = -1;
 225
 226         verify_point_data_fini(&data.vpd);
 227
 228         return r;
 229 }
 230
 231 /* Turn the set dimensions of "context" into parameters and return
 232  * the corresponding parameter domain.
 233  */
 234 static __isl_give isl_set *to_parameter_domain(__isl_take isl_set *context)
 235 {
 236         context = isl_set_move_dims(context, isl_dim_param, 0, isl_dim_set, 0,
 237                                        isl_set_dim(context, isl_dim_set));
 238         return isl_set_params(context);
 239 }
 240
 241 static int handle(isl_ctx *ctx, FILE *in, struct result_data *result,
 242         struct verify_options *options)
 243 {
 244     int i;
 245     int r;
 246     int nparam;
 247     isl_space *space;
 248     isl_set *set;
 249     isl_set *context;
 250     isl_pw_qpolynomial *pwqp[nr_methods];
 251
 252     set = isl_set_read_from_file(ctx, in);
 253     context = isl_set_read_from_file(ctx, in);
 254
 255     context = to_parameter_domain(context);
 256     nparam = isl_set_dim(context, isl_dim_param);
 257     if (nparam != isl_set_dim(set, isl_dim_param)) {
 258             int dim = isl_set_dim(set, isl_dim_set);
 259             set = isl_set_move_dims(set, isl_dim_param, 0,
 260                                         isl_dim_set, dim - nparam, nparam);
 261     }
 262
 263     set = isl_set_intersect_params(set, context);
 264     context = isl_set_params(isl_set_copy(set));
 265     space = isl_set_get_space(context);
 266
 267     context = verify_context_set_bounds(context, options);
 268
 269     for (i = 0; i < nr_methods; ++i) {
 270         struct tms st_cpu;
 271         struct tms en_cpu;
 272         options->barvinok->approx->approximation = methods[i].sign;
 273         options->barvinok->approx->method = methods[i].method;
 274         if (methods[i].method == BV_APPROX_SCALE)
 275             options->barvinok->approx->scale_flags = methods[i].flags;
 276         else if (methods[i].method == BV_APPROX_VOLUME)
 277             options->barvinok->approx->volume_triangulate = methods[i].flags;
 278
 279         times(&st_cpu);
 280         pwqp[i] = isl_set_card(isl_set_copy(set));
 281         times(&en_cpu);
 282         result->ticks[i] = time_diff(&en_cpu, &st_cpu);
 283     }
 284     for (i = 0; i < nr_methods; ++i)
 285         result->size[i] = 0;
 286     r = test(context, pwqp, result, options);
 287     for (i = 0; i < nr_methods; ++i)
 288         isl_pw_qpolynomial_free(pwqp[i]);
 289
 290     isl_space_free(space);
 291     isl_set_free(context);
 292     isl_set_free(set);
 293
 294     return r;
 295 }
 296
 297 int main(int argc, char **argv)
 298 {
 299     isl_ctx *ctx;
 300     char path[PATH_MAX+1];
 301     struct result_data all_result;
 302     int n = 0;
 303     int r = EXIT_SUCCESS;
 304     struct options *options = options_new_with_defaults();
 305
 306     argc = options_parse(options, argc, argv, ISL_ARG_ALL);
 307     ctx = isl_ctx_alloc_with_options(&options_args, options);
 308
 309     if (options->verify->M == INT_MIN)
 310         options->verify->M = 100;
 311     if (options->verify->m == INT_MAX)
 312         options->verify->m = -100;
 313
 314     result_data_init(&all_result);
 315
 316     while (fgets(path, sizeof(path), stdin)) {
 317         struct result_data result;
 318         FILE *in;
 319         int i;
 320
 321         ++n;
 322         result_data_init(&result);
 323         fprintf(stderr, "%s", path);
 324         *strchr(path, '\n') = '\0';
 325         in = fopen(path, "r");
 326         assert(in);
 327         if (handle(ctx, in, &result, options->verify) < 0)
 328             r = EXIT_FAILURE;
 329         fclose(in);
 330
 331         if (!options->quiet)
 332             result_data_print(&result, 1);
 333
 334         value_addto(all_result.n, all_result.n, result.n);
 335         for (i = 0; i < nr_methods; ++i) {
 336             all_result.RE_sum[i] += result.RE_sum[i];
 337             all_result.ticks[i] += result.ticks[i];
 338             all_result.size[i] += result.size[i];
 339         }
 340
 341         result_data_clear(&result);
 342
 343         if (!options->quiet) {
 344             fprintf(stderr, "average\n");
 345             result_data_print(&all_result, n);
 346         }
 347     }
 348
 349     result_data_clear(&all_result);
 350
 351     isl_ctx_free(ctx);
 352
 353     return r;
 354 }