summate.cc

   1 #include <assert.h>
   2 #include <iostream>
   3 #include <barvinok/barvinok.h>
   4 #include <barvinok/options.h>
   5 #include <barvinok/util.h>
   6 #include "argp.h"
   7 #include "progname.h"
   8 #include "evalue_convert.h"
   9 #include "evalue_read.h"
  10 #include "verify.h"
  11
  12 using std::cout;
  13 using std::cerr;
  14 using std::endl;
  15
  16 #define ALLOCN(type,n) (type*)malloc((n) * sizeof(type))
  17
  18 #define OPT_VARS            (BV_OPT_LAST+1)
  19
  20 struct argp_option argp_options[] = {
  21     { "variables",          OPT_VARS,   "list", 0,
  22         "comma separated list of variables over which to sum" },
  23     { "verbose",            'v',        0,      0, },
  24     { 0 }
  25 };
  26
  27 struct options {
  28     struct convert_options   convert;
  29     struct verify_options    verify;
  30     char* var_list;
  31     int verbose;
  32 };
  33
  34 static error_t parse_opt(int key, char *arg, struct argp_state *state)
  35 {
  36     struct options *options = (struct options*) state->input;
  37
  38     switch (key) {
  39     case ARGP_KEY_INIT:
  40         state->child_inputs[0] = &options->convert;
  41         state->child_inputs[1] = &options->verify;
  42         state->child_inputs[2] = options->verify.barvinok;
  43         options->var_list = NULL;
  44         options->verbose = 0;
  45         break;
  46     case 'v':
  47         options->verbose = 1;
  48         break;
  49     case OPT_VARS:
  50         options->var_list = strdup(arg);
  51         break;
  52     default:
  53         return ARGP_ERR_UNKNOWN;
  54     }
  55     return 0;
  56 }
  57
  58 static int check_poly_sum(const struct check_poly_data *data,
  59                           int nparam, Value *z,
  60                           const struct verify_options *options);
  61
  62 struct check_poly_sum_data : public check_poly_data {
  63     Polyhedron          **S;
  64     evalue               *EP;
  65     evalue               *sum;
  66
  67     check_poly_sum_data(Value *z, evalue *EP, evalue *sum) :
  68                         EP(EP), sum(sum) {
  69         this->z = z;
  70         this->check = check_poly_sum;
  71     }
  72 };
  73
  74 static void sum(Polyhedron *S, int pos, const check_poly_sum_data *data,
  75                 evalue *s, const struct verify_options *options)
  76 {
  77     if (!S) {
  78         evalue *e = evalue_eval(data->EP, data->z+1);
  79         eadd(e, s);
  80         free_evalue_refs(e);
  81         free(e);
  82     } else {
  83         Value LB, UB;
  84         int ok;
  85         value_init(LB);
  86         value_init(UB);
  87         ok = !(lower_upper_bounds(1+pos, S, data->z, &LB, &UB));
  88         assert(ok);
  89         for (; value_le(LB, UB); value_increment(LB, LB)) {
  90             value_assign(data->z[1+pos], LB);
  91             sum(S->next, pos+1, data, s, options);
  92         }
  93         value_set_si(data->z[1+pos], 0);
  94         value_clear(LB);
  95         value_clear(UB);
  96     }
  97 }
  98
  99 static evalue *sum(const check_poly_sum_data *data,
 100                           const struct verify_options *options)
 101 {
 102     evalue *s = evalue_zero();
 103     for (int i = 0; i < data->EP->x.p->size/2; ++i)
 104         if (!emptyQ2(data->S[i]))
 105             sum(data->S[i], 0, data, s, options);
 106     return s;
 107 }
 108
 109 static int check_poly_sum(const struct check_poly_data *data,
 110                           int nparam, Value *z,
 111                           const struct verify_options *options)
 112 {
 113     const check_poly_sum_data *sum_data;
 114     sum_data = static_cast<const check_poly_sum_data *>(data);
 115     evalue *e, *s;
 116     int k;
 117     int ok;
 118
 119     e = evalue_eval(sum_data->sum, z);
 120     if (options->print_all) {
 121         printf("sum(");
 122         value_print(stdout, VALUE_FMT, z[0]);
 123         for (k = 1; k < nparam; ++k) {
 124             printf(", ");
 125             value_print(stdout, VALUE_FMT, z[k]);
 126         }
 127         printf(") = ");
 128         value_print(stdout, VALUE_FMT, e->x.n);
 129         if (value_notone_p(e->d)) {
 130             printf("/");
 131             value_print(stdout, VALUE_FMT, e->d);
 132         }
 133     }
 134
 135     s = sum(sum_data, options);
 136
 137     if (options->print_all) {
 138         printf(", sum(EP) = ");
 139         value_print(stdout, VALUE_FMT, s->x.n);
 140         if (value_notone_p(s->d)) {
 141             printf("/");
 142             value_print(stdout, VALUE_FMT, s->d);
 143         }
 144         printf(". ");
 145     }
 146
 147     ok = eequal(e, s);
 148
 149     if (!ok) {
 150         printf("\n");
 151         fflush(stdout);
 152         fprintf(stderr,"Error !\n");
 153         fprintf(stderr,"sum( ");
 154         value_print(stderr, VALUE_FMT, z[0]);
 155         for (k = 1; k < nparam; ++k) {
 156             fprintf(stderr, ", ");
 157             value_print(stderr, VALUE_FMT, z[k]);
 158         }
 159         fprintf(stderr," ) should be ");
 160         value_print(stderr, VALUE_FMT, s->x.n);
 161         if (value_notone_p(s->d)) {
 162             fprintf(stderr, "/");
 163             value_print(stderr, VALUE_FMT, s->d);
 164         }
 165         fprintf(stderr,", while summation gives ");
 166         value_print(stderr, VALUE_FMT, e->x.n);
 167         if (value_notone_p(e->d)) {
 168             fprintf(stderr, "/");
 169             value_print(stderr, VALUE_FMT, e->d);
 170         }
 171         fprintf(stderr, ".\n");
 172     } else if (options->print_all)
 173         printf("OK.\n");
 174
 175     free_evalue_refs(s);
 176     free(s);
 177     free_evalue_refs(e);
 178     free(e);
 179
 180     return ok;
 181 }
 182
 183 static int verify(Polyhedron *P, evalue *sum, evalue *EP,
 184                   unsigned nvar, unsigned nparam, Vector *p,
 185                   struct verify_options *options)
 186 {
 187     Polyhedron *CS;
 188     unsigned MaxRays = options->barvinok->MaxRays;
 189     int error = 0;
 190
 191     CS = check_poly_context_scan(NULL, &P, P->Dimension, options);
 192
 193     check_poly_init(P, options);
 194
 195     if (!(CS && emptyQ2(CS))) {
 196         check_poly_sum_data data(p->p, EP, sum);
 197         data.S = ALLOCN(Polyhedron *, EP->x.p->size/2);
 198         for (int i = 0; i < EP->x.p->size/2; ++i) {
 199             Polyhedron *A = EVALUE_DOMAIN(EP->x.p->arr[2*i]);
 200             data.S[i] = Polyhedron_Scan(A, P, MaxRays & POL_NO_DUAL ? 0 : MaxRays);
 201         }
 202         error = !check_poly(CS, &data, nparam, 0, p->p+1+nvar, options);
 203         for (int i = 0; i < EP->x.p->size/2; ++i)
 204             Domain_Free(data.S[i]);
 205         free(data.S);
 206     }
 207
 208     if (!options->print_all)
 209         printf("\n");
 210
 211     if (CS) {
 212         Domain_Free(CS);
 213         Domain_Free(P);
 214     }
 215
 216     return error;
 217 }
 218
 219 /*
 220  * Project on final dim dimensions
 221  */
 222 Polyhedron *DomainProject(Polyhedron *D, unsigned dim, unsigned MaxRays)
 223 {
 224     Polyhedron *P;
 225     Polyhedron *R;
 226
 227     R = Polyhedron_Project(D, dim);
 228     for (P = D->next; P; P = P->next) {
 229         Polyhedron *R2 = Polyhedron_Project(P, dim);
 230         Polyhedron *R3 = DomainUnion(R, R2, MaxRays);
 231         Polyhedron_Free(R2);
 232         Domain_Free(R);
 233         R = R3;
 234     }
 235     return R;
 236 }
 237
 238 static int verify(evalue *EP, evalue *sum, unsigned nvar, unsigned nparam,
 239                   struct verify_options *options)
 240 {
 241     Vector *p;
 242
 243     p = Vector_Alloc(nvar+nparam+2);
 244     value_set_si(p->p[nvar+nparam+1], 1);
 245
 246     assert(value_zero_p(EP->d));
 247     assert(EP->x.p->type == partition);
 248
 249     Polyhedron *EP_D = EVALUE_DOMAIN(EP->x.p->arr[0]);
 250     Polyhedron *D = Polyhedron_Project(EP_D, nparam);
 251
 252     for (int i = 1; i < EP->x.p->size/2; ++i) {
 253         Polyhedron *D2 = D;
 254         EP_D = DomainProject(EVALUE_DOMAIN(EP->x.p->arr[2*i]), nparam,
 255                              options->barvinok->MaxRays);
 256         D = DomainUnion(EP_D, D, options->barvinok->MaxRays);
 257         Domain_Free(D2);
 258     }
 259
 260     int error = 0;
 261
 262     for (Polyhedron *P = D; P; P = P->next) {
 263         error = verify(P, sum, EP, nvar, nparam, p, options);
 264         if (error && !options->continue_on_error)
 265             break;
 266     }
 267
 268     Domain_Free(D);
 269     Vector_Free(p);
 270
 271     return error;
 272 }
 273
 274 int main(int argc, char **argv)
 275 {
 276     evalue *EP;
 277     char **all_vars = NULL;
 278     unsigned nvar;
 279     unsigned nparam;
 280     struct options options;
 281     struct barvinok_options *bv_options = barvinok_options_new_with_defaults();
 282     static struct argp_child argp_children[] = {
 283         { &convert_argp,        0,      "input conversion",     1 },
 284         { &verify_argp,         0,      "verification",         2 },
 285         { &barvinok_argp,       0,      "barvinok options",     3 },
 286         { 0 }
 287     };
 288     static struct argp argp = { argp_options, parse_opt, 0, 0, argp_children };
 289     int result = 0;
 290
 291     options.verify.barvinok = bv_options;
 292     set_program_name(argv[0]);
 293     argp_parse(&argp, argc, argv, 0, 0, &options);
 294
 295     EP = evalue_read_from_file(stdin, options.var_list, &all_vars,
 296                                &nvar, &nparam, bv_options->MaxRays);
 297     assert(EP);
 298
 299     if (options.verify.verify)
 300         verify_options_set_range(&options.verify, nvar+nparam);
 301
 302     evalue_convert(EP, &options.convert, options.verbose, nparam, all_vars);
 303
 304     if (EVALUE_IS_ZERO(*EP))
 305         print_evalue(stdout, EP, all_vars);
 306     else {
 307         evalue *sum = barvinok_summate(EP, nvar, bv_options);
 308         if (options.verify.verify)
 309             result = verify(EP, sum, nvar, nparam, &options.verify);
 310         else
 311             print_evalue(stdout, sum, all_vars+nvar);
 312         free_evalue_refs(sum);
 313         free(sum);
 314     }
 315
 316     free_evalue_refs(EP);
 317     free(EP);
 318
 319     if (options.var_list)
 320         free(options.var_list);
 321     Free_ParamNames(all_vars, nvar+nparam);
 322     barvinok_options_free(bv_options);
 323     return result;
 324 }