omega/polyfunc.cc

   1 #include <bernstein/bernstein.h>
   2 #include <bernstein/piecewise_lst.h>
   3 #include <barvinok/options.h>
   4 #include <barvinok/bernstein.h>
   5 #include "polyfunc.h"
   6 #include "param_util.h"
   7 #include "omega/convert.h"
   8
   9 using namespace GiNaC;
  10 using namespace bernstein;
  11 using namespace std;
  12
  13 void maximize(PolyFunc *polyfunc, Map<Variable_Ref *, GiNaC::ex>& variableMap)
  14 {
  15     varvector vv;
  16     varvector params;
  17     Param_Polyhedron *PP;
  18     exvector exvars;
  19     exvector exparams;
  20     struct barvinok_options *options = barvinok_options_new_with_defaults();
  21
  22     cout << "maximize " << polyfunc->poly << " over ";
  23     polyfunc->domain.simplify();
  24     polyfunc->domain.print(stdout);
  25
  26     Polyhedron *D = relation2Domain(polyfunc->domain, vv, params);
  27     assert(!D->next);
  28     assert(polyfunc->domain.is_set());
  29     int dim = polyfunc->domain.n_set();
  30     assert(D->Dimension == dim + params.size());
  31     Polyhedron *ctx = Universe_Polyhedron(params.size());
  32
  33     for (int i = 0; i < dim; ++i) {
  34         ex var = 0;
  35         foreach_map(vr,Variable_Ref *,s,ex,variableMap, {
  36             if (vr->vid == vv[i]) {
  37                 var = s;
  38                 break;
  39             }
  40         });
  41         assert(var != 0);
  42         exvars.push_back(var);
  43     }
  44     for (int i = dim; i < vv.size(); ++i) {
  45         Global_Var_ID global = vv[i]->get_global_var();
  46         ex var = 0;
  47         foreach_map(vr,Variable_Ref *,s,ex,variableMap, {
  48             if (vr->g == global) {
  49                 var = s;
  50                 break;
  51             }
  52         });
  53         if (var != 0)
  54             exparams.push_back(var);
  55         else
  56             exparams.push_back(symbol(vv[i]->char_name()));
  57     }
  58
  59     PP = Polyhedron2Param_Polyhedron(D, ctx, options);
  60     piecewise_lst *pl = new piecewise_lst(exparams);
  61     for (Param_Domain *Q = PP->D; Q; Q = Q->next) {
  62         GiNaC::matrix VM = domainVertices(PP, Q, exparams);
  63         lst coeffs = bernsteinExpansion(VM, polyfunc->poly, exvars, exparams);
  64         pl->list.push_back(guarded_lst(Polyhedron_Copy(Q->Domain), coeffs));
  65     }
  66     cout << "coefficients: " << *pl << endl;
  67     pl->maximize();
  68     cout << "maximum: " << *pl << endl;
  69     delete pl;
  70     Param_Polyhedron_Free(PP);
  71
  72     Polyhedron_Free(ctx);
  73     Domain_Free(D);
  74     barvinok_options_free(options);
  75 }