verify_lexsmaller.c

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