verify.c

   1 #include <assert.h>
   2 #include <stdlib.h>
   3 #include <barvinok/options.h>
   4 #include <barvinok/util.h>
   5 #include "verify.h"
   6
   7 /* RANGE : normal range for evalutations (-RANGE -> RANGE) */
   8 #define RANGE 50
   9
  10 /* SRANGE : small range for evalutations */
  11 #define SRANGE 15
  12
  13 /* if dimension >= BIDDIM, use SRANGE */
  14 #define BIGDIM 5
  15
  16 /* VSRANGE : very small range for evalutations */
  17 #define VSRANGE 5
  18
  19 /* if dimension >= VBIDDIM, use VSRANGE */
  20 #define VBIGDIM 8
  21
  22 static struct argp_option argp_options[] = {
  23     { "verify",             'T',    0,      0 },
  24     { "exact",              'E',    0,      0 },
  25     { "print-all",          'A',    0,      0 },
  26     { "continue-on-error",  'C',    0,      0 },
  27     { "min",                'm',    "int",  0 },
  28     { "max",                'M',    "int",  0 },
  29     { "range",              'r',    "int",  0 },
  30     { 0 }
  31 };
  32
  33 static error_t parse_opt(int key, char *arg, struct argp_state *state)
  34 {
  35     struct verify_options *options = state->input;
  36
  37     switch (key) {
  38     case ARGP_KEY_INIT:
  39         options->verify = 0;
  40         options->exact = 0;
  41         options->print_all = 0;
  42         options->continue_on_error = 0;
  43         options->m = INT_MAX;
  44         options->M = INT_MIN;
  45         options->r = -1;
  46         break;
  47     case ARGP_KEY_FINI:
  48         break;
  49     case 'T':
  50         options->verify = 1;
  51         break;
  52     case 'E':
  53         options->exact = 1;
  54         break;
  55     case 'A':
  56         options->print_all = 1;
  57         break;
  58     case 'C':
  59         options->continue_on_error = 1;
  60         break;
  61     case 'm':
  62         options->m = atoi(arg);
  63         options->verify = 1;
  64         break;
  65     case 'M':
  66         options->M = atoi(arg);
  67         options->verify = 1;
  68         break;
  69     case 'r':
  70         options->r = atoi(arg);
  71         options->M = atoi(arg);
  72         options->m = -options->M;
  73         options->verify = 1;
  74         break;
  75     default:
  76         return ARGP_ERR_UNKNOWN;
  77     }
  78     return 0;
  79 }
  80
  81 void verify_options_set_range(struct verify_options *options, int dim)
  82 {
  83     int r;
  84
  85     if (dim >= VBIGDIM)
  86         r = VSRANGE;
  87     else if (dim >= BIGDIM)
  88         r = SRANGE;
  89     else
  90         r = RANGE;
  91     /* If the user didn't set m or M, then we try to adjust the window
  92      * to the context in check_poly_context_scan.
  93      */
  94     if (options->m == INT_MAX && options->M == INT_MIN)
  95         options->r = r;
  96     if (options->M == INT_MIN)
  97         options->M = r;
  98     if (options->m == INT_MAX)
  99         options->m = -r;
 100
 101     if (options->verify && options->m > options->M) {
 102         fprintf(stderr,"Nothing to do: min > max !\n");
 103         exit(0);
 104     }
 105 }
 106
 107 struct argp verify_argp = {
 108     argp_options, parse_opt, 0, 0
 109 };
 110
 111 static Polyhedron *project_on(Polyhedron *P, int i)
 112 {
 113     unsigned dim = P->Dimension;
 114     Matrix *T;
 115     Polyhedron *I;
 116
 117     if (dim == 1)
 118         return Polyhedron_Copy(P);
 119
 120     T = Matrix_Alloc(2, dim+1);
 121     value_set_si(T->p[0][i], 1);
 122     value_set_si(T->p[1][dim], 1);
 123     I = Polyhedron_Image(P, T, P->NbConstraints);
 124     Matrix_Free(T);
 125     return I;
 126 }
 127
 128 static void set_bounds(Polyhedron *C, Value **rows, int i, unsigned nparam,
 129                         const struct verify_options *options)
 130 {
 131     Value min;
 132     Value max;
 133
 134     value_init(min);
 135     value_init(max);
 136     value_set_si(min, options->m);
 137     value_set_si(max, options->M);
 138
 139     if (options->r > 0) {
 140         Polyhedron *I = project_on(C, i);
 141         line_minmax(I, &min, &max);
 142         if (value_cmp_si(min, options->M) >= 0)
 143             value_add_int(max, min, options->r);
 144         else if (value_cmp_si(max, options->m) <= 0)
 145             value_sub_int(min, max, options->r);
 146         else {
 147             value_set_si(min, options->m);
 148             value_set_si(max, options->M);
 149         }
 150     }
 151
 152     value_set_si(rows[0][0], 1);
 153     value_set_si(rows[0][1+i], 1);
 154     value_oppose(rows[0][1+nparam], min);
 155     value_set_si(rows[1][0], 1);
 156     value_set_si(rows[1][1+i], -1);
 157     value_assign(rows[1][1+nparam], max);
 158
 159     value_clear(min);
 160     value_clear(max);
 161 }
 162
 163 Polyhedron *check_poly_context_scan(Polyhedron *P, Polyhedron **C,
 164                                     unsigned nparam,
 165                                     const struct verify_options *options)
 166 {
 167     int i;
 168     Matrix *MM;
 169     Polyhedron *CC, *CC2, *CS, *U;
 170     unsigned MaxRays = options->barvinok->MaxRays;
 171
 172     if (nparam <= 0)
 173         return NULL;
 174
 175     if (!P)
 176         CC = *C;
 177     else {
 178         CC = Polyhedron_Project(P, nparam);
 179         if (*C) {
 180             CC2 = DomainIntersection(*C, CC, MaxRays);
 181             Domain_Free(CC);
 182             CC = CC2;
 183         }
 184     }
 185
 186     /* Intersect context with range */
 187     MM = Matrix_Alloc(2*nparam, nparam+2);
 188     for (i = 0; i < nparam; ++i)
 189         set_bounds(CC, &MM->p[2*i], i, nparam, options);
 190     CC2 = AddConstraints(MM->p[0], 2*nparam, CC, options->barvinok->MaxRays);
 191     if (CC != *C)
 192         Domain_Free(CC);
 193     CC = CC2;
 194     U = Universe_Polyhedron(0);
 195     CS = Polyhedron_Scan(CC, U, MaxRays & POL_NO_DUAL ? 0 : MaxRays);
 196     Polyhedron_Free(U);
 197     *C = CC;
 198     Matrix_Free(MM);
 199     return CS;
 200 }
 201
 202 void check_poly_init(Polyhedron *C, struct verify_options *options)
 203 {
 204     int d, i;
 205
 206     if (options->print_all)
 207         return;
 208     if (C->Dimension <= 0)
 209         return;
 210
 211     d = options->M - options->m;
 212     if (d > 80)
 213         options->st = 1+d/80;
 214     else
 215         options->st = 1;
 216     for (i = options->m; i <= options->M; i += options->st)
 217         printf(".");
 218     printf( "\r" );
 219     fflush(stdout);
 220 }
 221
 222 static void print_rational(FILE *out, Value n, Value d)
 223 {
 224     value_print(out, VALUE_FMT, n);
 225     if (value_notone_p(d)) {
 226         fprintf(out, "/");
 227         value_print(out, VALUE_FMT, d);
 228     }
 229 }
 230
 231 void check_poly_print(int ok, int nparam, Value *z,
 232                       Value want_n, Value want_d,
 233                       Value got_n, Value got_d,
 234                       const char *op, const char *check,
 235                       const char *long_op,
 236                       const struct verify_options *options)
 237 {
 238     int k;
 239
 240     if (options->print_all) {
 241         printf("%s(", op);
 242         value_print(stdout, VALUE_FMT, z[0]);
 243         for (k = 1; k < nparam; ++k) {
 244             printf(", ");
 245             value_print(stdout, VALUE_FMT, z[k]);
 246         }
 247         printf(") = ");
 248         print_rational(stdout, got_n, got_d);
 249         printf(", %s = ", check);
 250         print_rational(stdout, want_n, want_d);
 251         printf(". ");
 252     }
 253
 254     if (!ok) {
 255         printf("\n");
 256         fflush(stdout);
 257         fprintf(stderr, "Error !\n");
 258         fprintf(stderr, "%s(", op);
 259         value_print(stderr, VALUE_FMT, z[0]);
 260         for (k = 1; k < nparam; ++k) {
 261             fprintf(stderr,", ");
 262             value_print(stderr, VALUE_FMT, z[k]);
 263         }
 264         fprintf(stderr, ") should be ");
 265         print_rational(stderr, want_n, want_d);
 266         fprintf(stderr, ", while %s gives ", long_op);
 267         print_rational(stderr, got_n, got_d);
 268         fprintf(stderr, ".\n");
 269     } else if (options->print_all)
 270         printf("OK.\n");
 271 }
 272
 273 /****************************************************/
 274 /* function check_poly :                            */
 275 /* scans the parameter space from Min to Max (all   */
 276 /* directions). Computes the number of points in    */
 277 /* the polytope using both methods, and compare them*/
 278 /* returns 1 on success                             */
 279 /****************************************************/
 280
 281 int check_poly(Polyhedron *CS, const struct check_poly_data *data,
 282                int nparam, int pos, Value *z,
 283                const struct verify_options *options)
 284 {
 285     if (pos == nparam) {
 286         if (!data->check(data, nparam, z, options) && !options->continue_on_error)
 287             return 0;
 288     } else {
 289         Value LB, UB;
 290         int ok;
 291         value_init(LB);
 292         value_init(UB);
 293         ok = !(lower_upper_bounds(1+pos, CS, z-1, &LB, &UB));
 294         assert(ok);
 295         for (; value_le(LB, UB); value_increment(LB, LB)) {
 296             if (!options->print_all) {
 297                 int k = VALUE_TO_INT(LB);
 298                 if (!pos && !(k % options->st)) {
 299                     printf("o");
 300                     fflush(stdout);
 301                 }
 302             }
 303
 304             value_assign(z[pos], LB);
 305             if (!check_poly(CS->next, data, nparam, pos+1, z, options)) {
 306                 value_clear(LB);
 307                 value_clear(UB);
 308                 return 0;
 309             }
 310         }
 311         value_set_si(z[pos], 0);
 312         value_clear(LB);
 313         value_clear(UB);
 314     }
 315     return 1;
 316 } /* check_poly */