verify_lexsmaller.c

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