verify_lexsmaller.c

   1 #include <assert.h>
   2 #include <barvinok/options.h>
   3 #include <barvinok/util.h>
   4 #include "argp.h"
   5 #include "progname.h"
   6
   7 struct argp_option argp_options[] = {
   8     { "continue",  'k' },
   9     { "max",       'M' },
  10     { "live",      'l' },
  11     { 0 }
  12 };
  13
  14 struct arguments {
  15     int live;
  16     int print_max;
  17     int keep_going;
  18     Value max;
  19     struct barvinok_options     *barvinok;
  20 };
  21
  22 static error_t parse_opt(int key, char *arg, struct argp_state *state)
  23 {
  24     struct arguments *options = (struct arguments*) state->input;
  25
  26     switch (key) {
  27     case ARGP_KEY_INIT:
  28         state->child_inputs[0] = options->barvinok;
  29         options->keep_going = 0;
  30         options->print_max = 0;
  31         options->live = 0;
  32         break;
  33     case 'k':
  34         options->keep_going = 1;
  35         break;
  36     case 'M':
  37         options->print_max = 1;
  38         break;
  39     case 'l':
  40         options->live = 1;
  41         break;
  42     default:
  43         return ARGP_ERR_UNKNOWN;
  44     }
  45     return 0;
  46 }
  47
  48 #define LS_OK       1
  49 #define LS_P        2       /* continue searching P */
  50 #define LS_D        4       /* continue searching D */
  51
  52 static int check_lexsmaller(Polyhedron *SP, Polyhedron *SD, Enumeration *en,
  53                          int pos, int nvar, Value *zP, Value *zD, Value *zE,
  54                          Value *count, struct arguments *options)
  55 {
  56     int i;
  57     int ok;
  58     Value PLB, PUB, DLB, DUB, LB, UB, tmp, c;
  59
  60     if (!SP && !SD)
  61         return LS_OK | LS_P | LS_D;
  62
  63     value_init(PLB); value_init(PUB);
  64     value_init(DLB); value_init(DUB);
  65     value_init(LB); value_init(UB);
  66     value_init(tmp);
  67
  68     if (SP) {
  69         ok = !(lower_upper_bounds(1+pos, SP, zP, &PLB, &PUB));
  70         assert(ok);
  71     }
  72     if (SD) {
  73         ok = !(lower_upper_bounds(1+pos, SD, zD, &DLB, &DUB));
  74         assert(ok);
  75     }
  76     if (!SD || (SP && value_lt(PLB, DLB)))
  77         value_assign(LB, PLB);
  78     else
  79         value_assign(LB, DLB);
  80     if (!SD || (SP && value_gt(PUB, DUB)))
  81         value_assign(UB, PUB);
  82     else
  83         value_assign(UB, DUB);
  84
  85     if (SD && !SD->next)
  86         value_init(c);
  87
  88     ok = LS_OK | LS_P | LS_D;
  89
  90     for(value_assign(tmp,LB); value_le(tmp,UB); value_increment(tmp,tmp)) {
  91         int inP = SP && value_ge(tmp, PLB) && value_le(tmp, PUB);
  92         int inD = SD && value_ge(tmp, DLB) && value_le(tmp, DUB);
  93         if (!inP && !inD)
  94             continue;
  95
  96         if (inP)
  97             value_assign(zP[pos+1], tmp);
  98         if (inD)
  99             value_assign(zD[pos+1], tmp);
 100         if (inD && pos < nvar)
 101             value_assign(zE[pos], tmp);
 102
 103         if (inD && !SD->next) {
 104             Value *ctmp;
 105
 106             value_assign(c,*(ctmp=compute_poly(en, zE)));
 107             value_clear(*ctmp);
 108             free(ctmp);
 109
 110             if (options->barvinok->verbose >= 2) {
 111                 printf("EP( ");
 112                 value_print(stdout, VALUE_FMT, zE[0]);
 113                 for (i = 1; i < nvar; ++i) {
 114                     printf(", ");
 115                     value_print(stdout, VALUE_FMT, zE[i]);
 116                 }
 117                 printf(" ) = ");
 118                 value_print(stdout, VALUE_FMT, c);
 119                 printf("; count = ");
 120                 value_print(stdout, VALUE_FMT, *count);
 121                 printf("\n");
 122             }
 123
 124             if (value_ne(c, *count)) {
 125                 printf("\n");
 126                 fflush(stdout);
 127                 fprintf(stderr,"Error !\n");
 128                 fprintf(stderr, "EP( ");
 129                 value_print(stderr, VALUE_FMT, zE[0]);
 130                 for (i = 1; i < nvar; ++i) {
 131                     fprintf(stderr, ", ");
 132                     value_print(stderr, VALUE_FMT, zE[i]);
 133                 }
 134                 fprintf(stderr, " ) = ");
 135                 value_print(stderr, VALUE_FMT, c);
 136                 fprintf(stderr, " but count = ");
 137                 value_print(stderr, VALUE_FMT, *count);
 138                 printf("\n");
 139                 ok = 0;
 140             }
 141
 142             if (options->live)
 143                 value_decrement(*count, *count);
 144
 145             ok &= ~LS_D;
 146         }
 147
 148         if (pos < nvar-1)
 149             ok &= check_lexsmaller(inP ? SP->next : NULL,
 150                                    inD ? SD->next : NULL,
 151                                    en, pos+1, nvar, zP, zD, zE, count,
 152                                    options);
 153         else {
 154             ok &= check_lexsmaller(NULL, inD ? SD->next : NULL,
 155                                    en, pos+1, nvar, zP, zD, zE, count,
 156                                    options)
 157                &  check_lexsmaller(inP ? SP->next : NULL, NULL,
 158                                    en, pos+1, nvar, zP, zD, zE, count,
 159                                    options);
 160             if (pos >= nvar && !(ok & LS_D))
 161                 break;
 162             if (pos >= nvar && !(ok & LS_P))
 163                 break;
 164         }
 165
 166         if (!ok && !options->keep_going)
 167             goto end;
 168
 169         if (inP && !SP->next) {
 170             value_increment(*count, *count);
 171             if (value_gt(*count, options->max))
 172                 value_assign(options->max, *count);
 173             ok &= ~LS_P;
 174         }
 175     }
 176     if (SP)
 177         value_set_si(zP[pos+1], 0);
 178     if (SD)
 179         value_set_si(zD[pos+1], 0);
 180
 181 end:
 182     if (SD && !SD->next)
 183         value_clear(c);
 184
 185     value_clear(PLB); value_clear(PUB);
 186     value_clear(DLB); value_clear(DUB);
 187     value_clear(LB); value_clear(UB);
 188     value_clear(tmp);
 189
 190     return ok;
 191 }
 192
 193 int main(int argc,char *argv[])
 194 {
 195     Matrix *M;
 196     Polyhedron *P, *D, *C;
 197     Polyhedron *SP, *SD;
 198     int nb_parms;
 199     const char **param_name = NULL;
 200     evalue *EP;
 201     Enumeration *en;
 202     Vector *zP, *zD, *zE;
 203     Value count;
 204     char s[128];
 205     unsigned dim;
 206     struct arguments options;
 207     static struct argp_child argp_children[] = {
 208         { &barvinok_argp,       0,      0,              0 },
 209         { 0 }
 210     };
 211     static struct argp argp = { argp_options, parse_opt, 0, 0, argp_children };
 212     struct barvinok_options *bv_options = barvinok_options_new_with_defaults();
 213
 214     options.barvinok = bv_options;
 215     set_program_name(argv[0]);
 216     argp_parse(&argp, argc, argv, 0, 0, &options);
 217
 218     M = Matrix_Read();
 219     P = Constraints2Polyhedron(M, bv_options->MaxRays);
 220     assert(P != NULL);
 221     Matrix_Free(M);
 222     M = Matrix_Read();
 223     D = Constraints2Polyhedron(M, bv_options->MaxRays);
 224     assert(D != NULL);
 225     Matrix_Free(M);
 226
 227     fgets(s, 128, stdin);
 228     while ((*s=='#') || (sscanf(s, "D %u", &dim)<1))
 229         fgets(s, 128, stdin);
 230
 231     M = Matrix_Read();
 232     C = Constraints2Polyhedron(M, bv_options->MaxRays);
 233     assert(C != NULL);
 234     Matrix_Free(M);
 235
 236     nb_parms = D->Dimension;
 237     param_name = Read_ParamNames(stdin, nb_parms);
 238
 239     EP = barvinok_lexsmaller_ev(P, D, dim, C, bv_options->MaxRays);
 240     if (options.live) {
 241         evalue mone;
 242         evalue *EC = barvinok_lexsmaller_ev(D, D, dim, C, bv_options->MaxRays);
 243         if (options.barvinok->verbose >= 2) {
 244             puts("EP");
 245             print_evalue(stdout, EP, (const char **)param_name);
 246             puts("EC");
 247             print_evalue(stdout, EC, (const char **)param_name);
 248         }
 249         value_init(mone.d);
 250         evalue_set_si(&mone, -1, 1);
 251         emul(&mone, EC);
 252         eadd(EC, EP);
 253         free_evalue_refs(&mone);
 254         evalue_free(EC);
 255         reduce_evalue(EP);
 256     }
 257     if (bv_options->verbose >= 1) {
 258         puts("Enumeration");
 259         print_evalue(stdout, EP, (const char **)param_name);
 260     }
 261     en = partition2enumeration(EP);
 262
 263     assert(C->Dimension == 0); /* for now */
 264
 265     /* S = scanning list of polyhedra */
 266     SP = Polyhedron_Scan(P, C, bv_options->MaxRays);
 267     SD = Polyhedron_Scan(D, C, bv_options->MaxRays);
 268
 269     zP = Vector_Alloc(1+P->Dimension+1);
 270     value_set_si(zP->p[1+P->Dimension], 1);
 271     zD = Vector_Alloc(1+D->Dimension+1);
 272     value_set_si(zD->p[1+D->Dimension], 1);
 273     zE = Vector_Alloc(dim+C->Dimension);
 274
 275     if (options.print_max)
 276         value_init(options.max);
 277     value_init(count);
 278     check_lexsmaller(SP, SD, en, 0, dim, zP->p, zD->p, zE->p, &count, &options);
 279     value_clear(count);
 280
 281     if (options.print_max) {
 282         printf("max = ");
 283         value_print(stdout, VALUE_FMT, options.max);
 284         printf("\n");
 285         value_clear(options.max);
 286     }
 287
 288     Enumeration_Free(en);
 289     Free_ParamNames(param_name, nb_parms);
 290     Polyhedron_Free(P);
 291     Polyhedron_Free(D);
 292     Polyhedron_Free(C);
 293     Domain_Free(SP);
 294     Domain_Free(SD);
 295     Vector_Free(zP);
 296     Vector_Free(zD);
 297     Vector_Free(zE);
 298     barvinok_options_free(bv_options);
 299     return 0;
 300 }