omega_interface/count.cc

   1 #include <vector>
   2 #include <assert.h>
   3 #include <barvinok/barvinok.h>
   4 #include <barvinok/util.h>
   5 #include <omega.h>
   6 #include "omega_interface/convert.h"
   7 #include "normalization.h"
   8 #include "count.h"
   9 #include "config.h"
  10
  11 #include <iostream>
  12 using namespace std;
  13
  14 #define ALLOC(t,p) p = (t*)malloc(sizeof(*p))
  15
  16 #ifdef USE_PARKER
  17
  18 #include "parker/count_solutions.h"
  19
  20 /*
  21  * Use parker's method to compute the number of integer points in D.
  22  * Since this method assumes all variables are non-negative,
  23  * we have to transform the input polytope first.
  24  */
  25 evalue *barvinok_enumerate_parker(Polyhedron *D,
  26                                         unsigned nvar, unsigned nparam,
  27                                         unsigned MaxRays)
  28 {
  29     Polyhedron *R;
  30     evalue *res;
  31
  32     if (nparam != 0) {
  33         fprintf(stderr, "parker method doesn't support parameters\n");
  34         return NULL;
  35     }
  36     R = skew_to_positive_orthant(D, nvar, MaxRays);
  37     Relation r = Domain2relation(R, nvar, 0, NULL);
  38     Polyhedron_Free(R);
  39     double d = count_solutions(r);
  40     ALLOC(evalue, res);
  41     value_init(res->d);
  42     value_set_si(res->d, 0);
  43     res->x.p = new_enode(::partition, 2, 0);
  44     EVALUE_SET_DOMAIN(res->x.p->arr[0], Universe_Polyhedron(0));
  45     value_set_si(res->x.p->arr[1].d, 1);
  46     value_init(res->x.p->arr[1].x.n);
  47     value_set_double(res->x.p->arr[1].x.n, d);
  48     return res;
  49 }
  50
  51 #else
  52
  53 evalue *barvinok_enumerate_parker(Polyhedron *D,
  54                                         unsigned nvar, unsigned nparam,
  55                                         unsigned MaxRays)
  56 {
  57     fprintf(stderr, "support for parker method not compiled in\n");
  58     return NULL;
  59 }
  60
  61 #endif
  62
  63 static evalue *count_relation_barvinok(Polyhedron *D, unsigned nvar,
  64                                         unsigned nparam, unsigned MaxRays)
  65 {
  66     int d = 0;
  67     evalue *EP = NULL;
  68     for (Polyhedron *P = D; P; P = P->next, ++d) {
  69         assert(d == 0);
  70         int exist = P->Dimension - nparam - nvar;
  71         EP = barvinok_enumerate_e(P, exist, nparam, MaxRays);
  72     }
  73     reduce_evalue(EP);
  74     Domain_Free(D);
  75     return EP;
  76 }
  77
  78 evalue *count_relation(Relation& r, int counting_method)
  79 {
  80     varvector vv;
  81     varvector params;
  82     struct barvinok_options *options = barvinok_options_new_with_defaults();
  83     Polyhedron *D = relation2Domain(r, vv, params, options->MaxRays);
  84
  85     evalue *EP = NULL;
  86     if (counting_method == COUNT_RELATION_BARVINOK)
  87         EP = count_relation_barvinok(D, vv.size(), params.size(),
  88                                         options->MaxRays);
  89     else if (counting_method == COUNT_RELATION_PARKER)
  90         EP = barvinok_enumerate_parker(D, vv.size(), params.size(),
  91                                         options->MaxRays);
  92     barvinok_options_free(options);
  93     return EP;
  94 }
  95
  96 evalue *rank_relation(Relation& r)
  97 {
  98     varvector vv;
  99     varvector params;
 100     struct barvinok_options *options = barvinok_options_new_with_defaults();
 101     Polyhedron *D = relation2Domain(r, vv, params, options->MaxRays);
 102     int dim = r.is_set() ? r.n_set() : r.n_inp() + r.n_out();
 103
 104     evalue *EP = NULL;
 105     if (D) {
 106         assert(D->next == NULL);
 107         Polyhedron *C = Universe_Polyhedron(params.size());
 108         EP = barvinok_lexsmaller_ev(D, D, dim, C, options->MaxRays);
 109         Polyhedron_Free(C);
 110     }
 111     Domain_Free(D);
 112     barvinok_options_free(options);
 113     return EP;
 114 }
 115
 116 evalue *count_lexsmaller(Relation& r, Relation& domain)
 117 {
 118     varvector P_vv;
 119     varvector P_params;
 120     varvector D_vv;
 121     varvector D_params;
 122     struct barvinok_options *options = barvinok_options_new_with_defaults();
 123     Polyhedron *P = relation2Domain(r, P_vv, P_params, options->MaxRays);
 124     int P_dim = r.is_set() ? r.n_set() : r.n_inp() + r.n_out();
 125     Polyhedron *D = relation2Domain(domain, D_vv, D_params, options->MaxRays);
 126     int D_dim = r.is_set() ? r.n_set() : r.n_inp() + r.n_out();
 127     assert(P_dim == D_dim);
 128
 129     evalue *EP = NULL;
 130     if (P && D) {
 131         assert(P->next == NULL);
 132         assert(D->next == NULL);
 133         Polyhedron *C = Universe_Polyhedron(D_params.size());
 134         EP = barvinok_lexsmaller_ev(P, D, D_dim, C, options->MaxRays);
 135         Polyhedron_Free(C);
 136     }
 137     Domain_Free(P);
 138     Domain_Free(D);
 139     barvinok_options_free(options);
 140     return EP;
 141 }