omega/polyfunc.cc

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