summaryrefslogtreecommitdiff
path: root/omegalib/omega/src/hull.cc
diff options
context:
space:
mode:
Diffstat (limited to 'omegalib/omega/src/hull.cc')
-rw-r--r--omegalib/omega/src/hull.cc1489
1 files changed, 0 insertions, 1489 deletions
diff --git a/omegalib/omega/src/hull.cc b/omegalib/omega/src/hull.cc
deleted file mode 100644
index f1b0601..0000000
--- a/omegalib/omega/src/hull.cc
+++ /dev/null
@@ -1,1489 +0,0 @@
-/*****************************************************************************
- Copyright (C) 1994-2000 University of Maryland
- Copyright (C) 2008 University of Southern California
- Copyright (C) 2009-2010 University of Utah
- All Rights Reserved.
-
- Purpose:
- Various hull calculations.
-
- Notes:
-
- History:
- 06/15/09 ConvexRepresentation, Chun Chen
- 11/25/09 RectHull, Chun Chen
-*****************************************************************************/
-
-#include <omega.h>
-#include <omega/farkas.h>
-#include <omega/hull.h>
-#include <basic/Bag.h>
-#include <basic/Map.h>
-#include <basic/omega_error.h>
-#include <list>
-#include <vector>
-#include <set>
-
-namespace omega {
-
-int hull_debug = 0;
-
-Relation ConvexHull(NOT_CONST Relation &R) {
- Relation S = Approximate(consume_and_regurgitate(R));
- if (!S.is_upper_bound_satisfiable())
- return S;
- if (S.has_single_conjunct())
- return S;
- return Farkas(Farkas(S,Basic_Farkas), Convex_Combination_Farkas);
-}
-
-Relation DecoupledConvexHull(NOT_CONST Relation &R) {
- Relation S = Approximate(consume_and_regurgitate(R));
- if (!S.is_upper_bound_satisfiable())
- return S;
- if (S.has_single_conjunct())
- return S;
- return Farkas(Farkas(S,Decoupled_Farkas), Convex_Combination_Farkas);
-}
-
-Relation AffineHull(NOT_CONST Relation &R) {
- Relation S = Approximate(consume_and_regurgitate(R));
- if (!S.is_upper_bound_satisfiable())
- return S;
- return Farkas(Farkas(S,Basic_Farkas), Affine_Combination_Farkas);
-}
-
-Relation LinearHull(NOT_CONST Relation &R) {
- Relation S = Approximate(consume_and_regurgitate(R));
- if (!S.is_upper_bound_satisfiable())
- return S;
- return Farkas(Farkas(S,Basic_Farkas), Linear_Combination_Farkas);
-}
-
-Relation ConicHull(NOT_CONST Relation &R) {
- Relation S = Approximate(consume_and_regurgitate(R));
- if (!S.is_upper_bound_satisfiable())
- return S;
- return Farkas(Farkas(S,Basic_Farkas), Positive_Combination_Farkas);
-}
-
-
-Relation FastTightHull(NOT_CONST Relation &input_R, NOT_CONST Relation &input_H) {
- Relation R = Approximate(consume_and_regurgitate(input_R));
- Relation H = Approximate(consume_and_regurgitate(input_H));
-
- if (hull_debug) {
- fprintf(DebugFile,"[ Computing FastTightHull of:\n");
- R.prefix_print(DebugFile);
- fprintf(DebugFile,"given known hull of:\n");
- H.prefix_print(DebugFile);
- }
-
- if (!H.has_single_conjunct()) {
- if (hull_debug)
- fprintf(DebugFile, "] bailing out of FastTightHull, known hull not convex\n");
- return H;
- }
-
- if (!H.is_obvious_tautology()) {
- R = Gist(R,copy(H));
- R.simplify(1,0);
- }
-
- if (R.has_single_conjunct()) {
- R = Intersection(R,H);
- if (hull_debug) {
- fprintf(DebugFile, "] quick easy answer to FastTightHull\n");
- R.prefix_print(DebugFile);
- }
- return R;
- }
- if (R.has_local(coefficient_of_constant_term)) {
- if (hull_debug) {
- fprintf(DebugFile, "] Can't handle recursive application of Farkas lemma\n");
- }
- return H;
- }
-
- if (hull_debug) {
- fprintf(DebugFile,"Gist of R given H is:\n");
- R.prefix_print(DebugFile);
- }
-
- if (1) {
- Set<Variable_ID> vars;
- int conjuncts = 0;
- for (DNF_Iterator s(R.query_DNF()); s.live(); s.next()) {
- conjuncts++;
- for (Variable_ID_Iterator v(*((*s)->variables())); v.live(); v++) {
- bool found = false;
- for (EQ_Iterator eq = (*s)->EQs(); eq.live(); eq.next())
- if ((*eq).get_coef(*v) != 0) {
- if (!found) vars.insert(*v);
- found = true;
- break;
- }
- if (!found)
- for (GEQ_Iterator geq = (*s)->GEQs(); geq.live(); geq.next())
- if ((*geq).get_coef(*v) != 0) {
- if (!found) vars.insert(*v);
- found = true;
- break;
- }
- }
- }
-
-
- // We now know which variables appear in R
- if (hull_debug) {
- fprintf(DebugFile,"Variables we need a better hull on are: ");
- foreach(v,Variable_ID,vars,
- fprintf(DebugFile," %s",v->char_name()));
- fprintf(DebugFile,"\n");
- }
- Conjunct *c = H.single_conjunct();
- int total=0;
- int copied = 0;
- for (EQ_Iterator eq = c->EQs(); eq.live(); eq.next()) {
- total++;
- foreach(v,Variable_ID,vars,
- if ((*eq).get_coef(v) != 0) {
- R.and_with_EQ(*eq);
- copied++;
- break; // out of variable loop
- }
- );
- }
- for (GEQ_Iterator geq = c->GEQs(); geq.live(); geq.next()) {
- total++;
- foreach(v,Variable_ID,vars,
- if ((*geq).get_coef(v) != 0) {
- R.and_with_GEQ(*geq);
- copied++;
- break; // out of variable loop
- }
- );
- }
- if (copied < total) {
- R = Approximate(R);
-
- if (hull_debug) {
- fprintf(DebugFile,"Decomposed relation, copied only %d of %d constraints\n",copied,total);
- fprintf(DebugFile,"Original R:\n");
- R.prefix_print(DebugFile);
- fprintf(DebugFile,"Known hull:\n");
- H.prefix_print(DebugFile);
- fprintf(DebugFile,"New R:\n");
- R.prefix_print(DebugFile);
- }
- }
- }
-
- Relation F = Farkas(copy(R), Basic_Farkas, true);
- if (hull_debug)
- fprintf(DebugFile,"Farkas Difficulty = " coef_fmt "\n", farkasDifficulty);
- if (farkasDifficulty > 260) {
- if (hull_debug) {
- fprintf(DebugFile, "] bailing out, farkas is way too complex\n");
- fprintf(DebugFile,"Farkas:\n");
- F.prefix_print(DebugFile);
- }
- return H;
- }
- else if (farkasDifficulty > 130) {
- // Bail out
- if (hull_debug) {
- fprintf(DebugFile, coef_fmt " non-zeros in original farkas\n", farkasDifficulty);
- }
- Relation tmp = Farkas(R, Decoupled_Farkas, true);
-
- if (hull_debug) {
- fprintf(DebugFile, coef_fmt " non-zeros in decoupled farkas\n", farkasDifficulty);
- }
- if (farkasDifficulty > 260) {
- if (hull_debug) {
- fprintf(DebugFile, "] bailing out, farkas is way too complex\n");
- fprintf(DebugFile,"Farkas:\n");
- F.prefix_print(DebugFile);
- }
- return H;
- }
- else {
- if (farkasDifficulty > 130)
- R = Intersection(H, Farkas(tmp, Affine_Combination_Farkas, true));
- else R = Intersection(H,
- Intersection(Farkas(tmp, Convex_Combination_Farkas, true),
- Farkas(F, Affine_Combination_Farkas, true)));
- if (hull_debug) {
- fprintf(DebugFile, "] bailing out, farkas is too complex, using affine hull\n");
- fprintf(DebugFile,"Farkas:\n");
- F.prefix_print(DebugFile);
- fprintf(DebugFile,"Affine hull:\n");
- R.prefix_print(DebugFile);
- }
- return R;
- }
- }
-
- R = Intersection(H, Farkas(F, Convex_Combination_Farkas, true));
- if (hull_debug) {
- fprintf(DebugFile, "] Result of FastTightHull:\n");
- R.prefix_print(DebugFile);
- }
- return R;
-}
-
-
-
-namespace {
- bool parallel(const GEQ_Handle &g1, const GEQ_Handle &g2) {
- for(Constr_Vars_Iter cvi(g1, false); cvi; cvi++) {
- coef_t c1 = (*cvi).coef;
- coef_t c2 = g2.get_coef((*cvi).var);
- if (c1 != c2) return false;
- }
- {
- for(Constr_Vars_Iter cvi(g2, false); cvi; cvi++) {
- coef_t c1 = g1.get_coef((*cvi).var);
- coef_t c2 = (*cvi).coef;
- if (c1 != c2) return false;
- }
- }
- return true;
- }
-
-
- bool hull(const EQ_Handle &e, const GEQ_Handle &g, coef_t &hull) {
- int sign = 0;
- for(Constr_Vars_Iter cvi(e, false); cvi; cvi++) {
- coef_t c1 = (*cvi).coef;
- coef_t c2 = g.get_coef((*cvi).var);
- if (sign == 0) sign = (c1*c2>=0?1:-1);
- if (sign*c1 != c2) return false;
- }
- assert(sign != 0);
- {
- for(Constr_Vars_Iter cvi(g, false); cvi; cvi++) {
- coef_t c1 = e.get_coef((*cvi).var);
- coef_t c2 = (*cvi).coef;
- if (sign*c1 != c2) return false;
- }
- }
- hull = max(sign * e.get_const(), g.get_const());
- if (hull_debug) {
- fprintf(DebugFile,"Hull of:\n %s\n", e.print_to_string().c_str());
- fprintf(DebugFile," %s\n", g.print_to_string().c_str());
- fprintf(DebugFile,"is " coef_fmt "\n\n",hull);
- }
- return true;
- }
-
- bool eq(const EQ_Handle &e1, const EQ_Handle &e2) {
- int sign = 0;
- for(Constr_Vars_Iter cvi(e1, false); cvi; cvi++) {
- coef_t c1 = (*cvi).coef;
- coef_t c2 = e2.get_coef((*cvi).var);
- if (sign == 0) sign = (c1*c2>=0?1:-1);
- if (sign*c1 != c2) return false;
- }
- assert(sign != 0);
- {
- for(Constr_Vars_Iter cvi(e2, false); cvi; cvi++) {
- coef_t c1 = e1.get_coef((*cvi).var);
- coef_t c2 = (*cvi).coef;
- if (sign*c1 != c2) return false;
- }
- }
- return sign * e1.get_const() == e2.get_const();
- }
-}
-
-
-// This function is deprecated!!!
-Relation QuickHull(Relation &R) {
- Tuple<Relation> Rs(1);
- Rs[1] = R;
- return QuickHull(Rs);
-}
-
-
-// This function is deprecated!!!
-Relation QuickHull(Tuple<Relation> &Rs) {
- assert(!Rs.empty());
-
- // if (Rs.size() == 1) return Rs[1];
-
- Tuple<Relation> l_Rs;
- for (int i = 1; i <= Rs.size(); i++)
- for (DNF_Iterator c(Rs[i].query_DNF()); c; c++) {
- Relation r = Relation(Rs[i], c.curr());
- l_Rs.append(Approximate(r));
- }
-
- if (l_Rs.size() == 1)
- return l_Rs[1];
-
- Relation result = Relation::True(Rs[1]);
- result.copy_names(Rs[1]);
-
- use_ugly_names++;
-
- if (hull_debug > 1)
- for (int i = 1; i <= l_Rs.size(); i++) {
- fprintf(DebugFile,"#%d \n",i);
- l_Rs[i].prefix_print(DebugFile);
- }
-
-
-// Relation R = copy(Rs[1]);
-// for (int i = 2; i <= Rs.size(); i++)
-// R = Union(R,copy(Rs[i]));
-
-// #if 0
-// if (!R.is_set()) {
-// if (R.n_inp() == R.n_out()) {
-// Relation AC = DeltasToRelation(Hull(Deltas(copy(R),
-// min(R.n_inp(),R.n_out()))),
-// R.n_inp(),R.n_out());
-// Relation dH = Hull(Domain(copy(R)),false);
-// Relation rH = Hull(Range(copy(R)),false);
-// result = Intersection(AC,Cross_Product(dH,rH));
-// }
-// else {
-// Relation dH = Hull(Domain(copy(R)),false);
-// Relation rH = Hull(Range(copy(R)),false);
-// result = Cross_Product(dH,rH);
-// assert(Must_Be_Subset(copy(R),copy(result)));
-// }
-// }
-
-// #endif
-
- Conjunct *first;
- l_Rs[1] = EQs_to_GEQs(l_Rs[1]);
- first = l_Rs[1].single_conjunct();
- for (GEQ_Iterator candidate(first->GEQs()); candidate.live(); candidate.next()) {
- coef_t maxConstantTerm = (*candidate).get_const();
- bool found = true;
- if (hull_debug > 1) {
- fprintf(DebugFile,"searching for bound on:\n %s\n", (*candidate).print_to_string().c_str());
- }
- for (int i = 2; i <= l_Rs.size(); i++) {
- Conjunct *other = l_Rs[i].single_conjunct();
- bool found_for_i = false;
- for (GEQ_Iterator target(other->GEQs()); target.live(); target.next()) {
- if (hull_debug > 2) {
- fprintf(DebugFile,"candidate:\n %s\n", (*candidate).print_to_string().c_str());
- fprintf(DebugFile,"target:\n %s\n", (*target).print_to_string().c_str());
- }
- if (parallel(*candidate,*target)) {
- if (hull_debug > 1)
- fprintf(DebugFile,"Found bound:\n %s\n", (*target).print_to_string().c_str());
- maxConstantTerm = max(maxConstantTerm,(*target).get_const());
- found_for_i = true;
- break;
- }
- }
- if (!found_for_i) {
- for (EQ_Iterator target_e(other->EQs()); target_e.live(); target_e.next()) {
- coef_t h;
- if (hull(*target_e,*candidate,h)) {
- if (hull_debug > 1)
- fprintf(DebugFile,"Found bound of " coef_fmt ":\n %s\n", h, (*target_e).print_to_string().c_str());
- maxConstantTerm = max(maxConstantTerm,h);
- found_for_i = true;
- break;
- }
- };
- if (!found_for_i) {
- if (hull_debug > 1) {
- fprintf(DebugFile,"No bound found in:\n");
- fprintf(DebugFile, "%s", l_Rs[i].print_with_subs_to_string().c_str());
- }
- //if nothing found
- found = false;
- break;
- }
- }
- }
-
- if (found) {
- GEQ_Handle h = result.and_with_GEQ();
- copy_constraint(h,*candidate);
- if (hull_debug > 1)
- fprintf(DebugFile,"Setting constant term to " coef_fmt " in\n %s\n", maxConstantTerm, h.print_to_string().c_str());
- h.update_const(maxConstantTerm - (*candidate).get_const());
- if (hull_debug > 1)
- fprintf(DebugFile,"Updated constraint is\n %s\n", h.print_to_string().c_str());
- }
- }
-
-
- for (EQ_Iterator candidate_eq(first->EQs()); candidate_eq.live(); candidate_eq.next()) {
- bool found = true;
- for (int i = 2; i <= l_Rs.size(); i++) {
- Conjunct *C = l_Rs[i].single_conjunct();
- bool found_for_i = false;
-
- for (EQ_Iterator target(C->EQs()); target.live(); target.next()) {
- if (eq(*candidate_eq,*target)) {
- found_for_i = true;
- break;
- }
- }
- if (!found_for_i) {
- //if nothing found
- found = false;
- break;
- }
- }
-
- if (found) {
- EQ_Handle h = result.and_with_EQ();
- copy_constraint(h,*candidate_eq);
- if (hull_debug > 1)
- fprintf(DebugFile,"Adding eq constraint: %s\n", h.print_to_string().c_str());
- }
- }
-
- use_ugly_names--;
- if (hull_debug > 1) {
- fprintf(DebugFile,"quick hull is of:");
- result.print_with_subs(DebugFile);
- }
- return result;
-}
-
-
-// Relation Hull2(Tuple<Relation> &Rs, Tuple<int> &active) {
-// assert(Rs.size() == active.size() && Rs.size() > 0);
-
-// Tuple<Relation> l_Rs;
-// for (int i = 1; i <= Rs.size(); i++)
-// if (active[i])
-// l_Rs.append(copy(Rs[i]));
-
-// if (l_Rs.size() == 0)
-// return Relation::False(Rs[1]);
-
-// try {
-// Relation r = l_Rs[1];
-// for (int i = 2; i <= l_Rs.size(); i++) {
-// r = Union(r, copy(l_Rs[i]));
-// r.simplify();
-// }
-
-// // Relation F = Farkas(r, Basic_Farkas, true);
-// // if (farkasDifficulty >= 500)
-// // throw std::overflow_error("loop convex hull too complicated.");
-// // F = Farkas(F, Convex_Combination_Farkas, true);
-// return Farkas(Farkas(r, Basic_Farkas, true), Convex_Combination_Farkas, true);
-// }
-// catch (std::overflow_error) {
-// return QuickHull(l_Rs);
-// }
-// }
-
-
-namespace {
- void printRs(Tuple<Relation> &Rs) {
- fprintf(DebugFile,"Rs:\n");
- for (int i = 1; i <= Rs.size(); i++)
- fprintf(DebugFile,"#%d : %s\n",i,
- Rs[i].print_with_subs_to_string().c_str());
- }
-}
-
-Relation BetterHull(Tuple<Relation> &Rs, bool stridesAllowed, bool checkSubsets,
- NOT_CONST Relation &input_knownHull = Relation::Null()) {
- Relation knownHull = consume_and_regurgitate(input_knownHull);
- static int OMEGA_WHINGE = -1;
- if (OMEGA_WHINGE < 0) {
- OMEGA_WHINGE = getenv("OMEGA_WHINGE") ? atoi(getenv("OMEGA_WHINGE")) : 0;
- }
- assert(!Rs.empty());
- if (Rs.size() == 1) {
- if (stridesAllowed) return Rs[1];
- else return Approximate(Rs[1]);
- }
-
- if (checkSubsets) {
- Tuple<bool> live(Rs.size());
- if (hull_debug) {
- fprintf(DebugFile,"Checking subsets in hull computation:\n");
- printRs(Rs);
- }
- int i;
- for(i=1;i <=Rs.size(); i++) live[i] = true;
- for(i=1;i <=Rs.size(); i++)
- for(int j=1;j <=Rs.size(); j++) if (i != j && live[j]) {
- if (hull_debug) fprintf(DebugFile,"checking %d Is_Obvious_Subset %d\n",i,j);
- if (Is_Obvious_Subset(copy(Rs[i]),copy(Rs[j]))) {
- if (hull_debug) fprintf(DebugFile,"yes...\n");
- live[i] = false;
- break;
- }
- }
- for(i=1;i <=Rs.size(); i++) if (!live[i]) {
- if (i < Rs.size()) {
- Rs[i] = Rs[Rs.size()];
- live[i] = live[Rs.size()];
- }
- Rs[Rs.size()] = Relation();
- Rs.delete_last();
- i--;
- }
- }
- Relation hull;
- if (hull_debug) {
- fprintf(DebugFile,"Better Hull:\n");
- printRs(Rs);
- fprintf(DebugFile,"known hull: %s\n", knownHull.print_with_subs_to_string().c_str());
- }
- if (knownHull.is_null()) hull = QuickHull(Rs);
- else hull = Intersection(QuickHull(Rs),knownHull);
- // for (int i = 1; i <= Rs.size(); i++)
- // hull = RectHull(Union(hull, copy(Rs[i])));
- // hull = Intersection(hull, knownHull);
- hull.simplify();
- if (hull_debug) {
- fprintf(DebugFile,"quick hull: %s\n", hull.print_with_subs_to_string().c_str());
- }
-
- Relation orig = Relation::False(Rs[1]);
- int i;
- for (i = 1; i <= Rs.size(); i++)
- orig = Union(orig,copy(Rs[i]));
-
- orig.simplify();
-
- for (i = 1; i <= Rs.size(); i++) {
- if (!hull.is_obvious_tautology()) Rs[i] = Gist(Rs[i],copy(hull));
- Rs[i].simplify();
- if (Rs[i].is_obvious_tautology()) return hull;
- if (Rs[i].has_single_conjunct()) {
- Rs[i] = EQs_to_GEQs(Rs[i]);
- if (hull_debug) {
- fprintf(DebugFile,"Checking for hull constraints in:\n %s\n", Rs[i].print_with_subs_to_string().c_str());
- }
- Conjunct *c = Rs[i].single_conjunct();
- for (GEQ_Iterator g(c->GEQs()); g.live(); g.next()) {
- Relation tmp = Relation::True(Rs[i]);
- tmp.and_with_GEQ(*g);
- if (!Difference(copy(orig),tmp).is_upper_bound_satisfiable())
- hull.and_with_GEQ(*g);
- }
- for (EQ_Iterator e(c->EQs()); e.live(); e.next()) {
- Relation tmp = Relation::True(Rs[i]);
- tmp.and_with_EQ(*e);
- if (!Difference(copy(orig),tmp).is_upper_bound_satisfiable())
- hull.and_with_EQ(*e);
- }
- }
- }
-
- hull = FastTightHull(orig,hull);
- assert(hull.has_single_conjunct());
-
- if (stridesAllowed) return hull;
- else return Approximate(hull);
-
-}
-
-
-
-Relation Hull(NOT_CONST Relation &S,
- bool stridesAllowed,
- int effort,
- NOT_CONST Relation &knownHull) {
- Relation R = consume_and_regurgitate(S);
- R.simplify(1,0);
- if (!R.is_upper_bound_satisfiable()) return R;
- Tuple<Relation> Rs;
- for (DNF_Iterator c(R.query_DNF()); c.live(); ) {
- Rs.append(Relation(R,c.curr()));
- c.next();
- }
- if (effort == 1)
- return BetterHull(Rs,stridesAllowed,false,knownHull);
- else
- return QuickHull(Rs);
-}
-
-
-
-Relation Hull(Tuple<Relation> &Rs,
- Tuple<int> &validMask,
- int effort,
- bool stridesAllowed,
- NOT_CONST Relation &knownHull) {
- // Use relation of index i only when validMask[i] != 0
- Tuple<Relation> Rs2;
- for(int i = 1; i <= Rs.size(); i++) {
- if (validMask[i]) {
- Rs[i].simplify();
- for (DNF_Iterator c(Rs[i].query_DNF()); c.live(); ) {
- Rs2.append(Relation(Rs[i],c.curr()));
- c.next();
- }
- }
- }
- assert(effort == 0 || effort == 1);
- if (effort == 1)
- return BetterHull(Rs2,stridesAllowed,true,knownHull);
- else
- return QuickHull(Rs2);
-}
-
-
-// This function is deprecated!!!
-Relation CheckForConvexRepresentation(NOT_CONST Relation &R_In) {
- Relation R = consume_and_regurgitate(R_In);
- Relation h = Hull(copy(R));
- if (!Difference(copy(h),copy(R)).is_upper_bound_satisfiable())
- return h;
- else
- return R;
-}
-
-// This function is deprecated!!!
-Relation CheckForConvexPairs(NOT_CONST Relation &S) {
- Relation R = consume_and_regurgitate(S);
- Relation hull = FastTightHull(copy(R),Relation::True(R));
- R.simplify(1,0);
- if (!R.is_upper_bound_satisfiable() || R.number_of_conjuncts() < 2) return R;
- Tuple<Relation> Rs;
- for (DNF_Iterator c(R.query_DNF()); c.live(); ) {
- Rs.append(Relation(R,c.curr()));
- c.next();
- }
-
- bool *dead = new bool[Rs.size()+1];
- int i;
- for(i = 1; i<=Rs.size();i++) dead[i] = false;
-
- for(i = 1; i<=Rs.size();i++)
- if (!dead[i])
- for(int j = i+1; j<=Rs.size();j++) if (!dead[j]) {
- if (hull_debug) {
- fprintf(DebugFile,"Comparing #%d and %d\n",i,j);
- }
- Relation U = Union(copy(Rs[i]),copy(Rs[j]));
- Relation H_ij = FastTightHull(copy(U),copy(hull));
- if (!Difference(copy(H_ij),U).is_upper_bound_satisfiable()) {
- Rs[i] = H_ij;
- dead[j] = true;
- if (hull_debug) {
- fprintf(DebugFile,"Combined them\n");
- }
- }
- }
- i = 1;
- while(i<=Rs.size() && dead[i]) i++;
- assert(i<=Rs.size());
- R = Rs[i];
- i++;
- for(; i<=Rs.size();i++)
- if (!dead[i])
- R = Union(R,Rs[i]);
- delete []dead;
- return R;
-}
-
-//
-// Supporting functions for ConvexRepresentation
-//
-namespace {
-struct Interval {
- std::list<std::pair<Relation, Relation> >::iterator pos;
- coef_t lb;
- coef_t ub;
- bool change;
- coef_t modulo;
- Interval(std::list<std::pair<Relation, Relation> >::iterator pos_, coef_t lb_, coef_t ub_):
- pos(pos_), lb(lb_), ub(ub_) {}
- friend bool operator<(const Interval &a, const Interval &b);
-};
-
-bool operator<(const Interval &a, const Interval &b) {
- return a.lb < b.lb;
-}
-
-struct Modulo_Interval {
- coef_t modulo;
- coef_t start;
- coef_t size;
- Modulo_Interval(coef_t modulo_, coef_t start_, coef_t size_):
- modulo(modulo_), start(start_), size(size_) {}
- friend bool operator<(const Interval &a, const Interval &b);
-};
-
-bool operator<(const Modulo_Interval &a, const Modulo_Interval &b) {
- if (a.modulo == b.modulo) {
- if (a.start == b.start)
- return a.size < b.size;
- else
- return a.start < b.start;
- }
- else
- return a.modulo < b.modulo;
-}
-
-void merge_intervals(std::list<Interval> &intervals, coef_t modulo, std::list<std::pair<Relation, Relation> > &Rs, std::list<std::pair<Relation, Relation> >::iterator orig) {
- // normalize intervals
- for (std::list<Interval>::iterator i = intervals.begin(); i != intervals.end(); i++) {
- (*i).modulo = modulo;
- (*i).change = false;
- if ((*i).ub - (*i).lb + 1>= modulo) {
- (*i).lb = 0;
- (*i).ub = modulo - 1;
- }
- else if ((*i).ub < 0 || (*i).lb >= modulo) {
- coef_t range = (*i).ub - (*i).lb;
- (*i).lb = int_mod((*i).lb, modulo);
- (*i).ub = (*i).lb + range;
- }
- }
-
- intervals.sort();
-
- // merge neighboring intervals
- std::list<Interval>::iterator p = intervals.begin();
- while (p != intervals.end()) {
- std::list<Interval>::iterator q = p;
- q++;
- while (q != intervals.end()) {
- if ((*p).ub + 1 >= (*q).lb) {
- Relation hull = ConvexHull(Union(copy((*(*p).pos).first), copy((*(*q).pos).first)));
- Relation remainder = Difference(Difference(copy(hull), copy((*(*p).pos).first)), copy((*(*q).pos).first));
- if (!remainder.is_upper_bound_satisfiable()) {
- if ((*q).pos == orig)
- std::swap((*p).pos, (*q).pos);
- (*(*p).pos).first = hull;
- (*p).ub = max((*p).ub, (*q).ub);
- (*p).change = true;
- Rs.erase((*q).pos);
- q = intervals.erase(q);
- }
- else
- break;
- }
- else
- break;
- }
-
- bool p_moved = false;
- q = p;
- q++;
- while (q != intervals.end()) {
- if ((*q).ub >= modulo && int_mod((*q).ub, modulo) + 1 >= (*p).lb) {
- Relation hull = ConvexHull(Union(copy((*(*p).pos).first), copy((*(*q).pos).first)));
- Relation remainder = Difference(Difference(copy(hull), copy((*(*p).pos).first)), copy((*(*q).pos).first));
- if (!remainder.is_upper_bound_satisfiable()) {
- if ((*p).pos == orig)
- std::swap((*p).pos, (*q).pos);
- (*(*q).pos).first = hull;
- coef_t t = (*p).ub - int_mod((*q).ub, modulo);
- if (t > 0)
- (*q).ub = (*q).ub + t;
- (*q).change = true;
- Rs.erase((*p).pos);
- p = intervals.erase(p);
- p_moved = true;
- break;
- }
- else
- q++;
- }
- else
- q++;
- }
-
- if (!p_moved)
- p++;
- }
-
- // merge by reducing the strengh of modulo
- std::list<Modulo_Interval> modulo_intervals;
- coef_t max_distance = modulo/2;
- for (std::list<Interval>::iterator p = intervals.begin(); p != intervals.end(); p++) {
- if ((*p).lb >= max_distance)
- break;
-
- coef_t size = (*p).ub - (*p).lb;
-
- std::list<Interval>::iterator q = p;
- q++;
- while (q != intervals.end()) {
- coef_t distance = (*q).lb - (*p).lb;
- if (distance > max_distance)
- break;
-
- if ((*q).ub - (*q).lb != size || int_mod(modulo, distance) != 0) {
- q++;
- continue;
- }
-
- int num_reduced = 0;
- coef_t looking_for = int_mod((*p).lb, distance);
- for (std::list<Interval>::iterator k = intervals.begin(); k != intervals.end(); k++) {
- if ((*k).lb == looking_for && (*k).ub - (*k).lb == size) {
- num_reduced++;
- looking_for += distance;
- if (looking_for >= modulo)
- break;
- }
- else if ((*k).lb <= looking_for && (*k).ub >= looking_for + size) {
- looking_for += distance;
- if (looking_for >= modulo)
- break;
- }
- else if ((*k).lb > looking_for)
- break;
- }
-
- if (looking_for >= modulo && num_reduced > 1)
- modulo_intervals.push_back(Modulo_Interval(distance, int_mod((*p).lb, distance), size));
-
- q++;
- }
- }
-
- modulo_intervals.sort();
-
- // remove redundant reduced-strength intervals
- std::list<Modulo_Interval>::iterator p2 = modulo_intervals.begin();
- while (p2 != modulo_intervals.end()) {
- std::list<Modulo_Interval>::iterator q2 = p2;
- q2++;
- while (q2 != modulo_intervals.end()) {
- if ((*p2).modulo == (*q2).modulo && (*p2).start == (*q2).start)
- q2 = modulo_intervals.erase(q2);
- else if (int_mod((*q2).modulo, (*p2).modulo) == 0 &&
- (*p2).start == int_mod((*q2).start, (*p2).modulo) &&
- (*p2).size >= (*q2).size)
- q2 = modulo_intervals.erase(q2);
- else
- q2++;
- }
- p2++;
- }
-
- // replace original intervals with new reduced-strength ones
- for (std::list<Modulo_Interval>::iterator i = modulo_intervals.begin(); i != modulo_intervals.end(); i++) {
- std::vector<Relation *> candidates;
- int num_replaced = 0;
- for (std::list<Interval>::iterator j = intervals.begin(); j != intervals.end(); j++)
- if (int_mod((*j).modulo, (*i).modulo) == 0 &&
- (*j).ub - (*j).lb >= (*i).size &&
- (int_mod((*j).lb, (*i).modulo) == (*i).start ||
- int_mod((*j).ub, (*i).modulo) == (*i).start + (*i).size)) {
- candidates.push_back(&((*(*j).pos).first));
- if (int_mod((*j).lb, (*i).modulo) == (*i).start &&
- (*j).ub - (*j).lb == (*i).size)
- num_replaced++;
- }
- if (num_replaced <= 1)
- continue;
-
- Relation R = copy(*candidates[0]);
- for (size_t k = 1; k < candidates.size(); k++)
- R = Union(R, copy(*candidates[k]));
- Relation hull = ConvexHull(copy(R));
- Relation remainder = Difference(copy(hull), copy(R));
- if (!remainder.is_upper_bound_satisfiable()) {
- std::list<Interval>::iterator replaced_one = intervals.end();
- for (std::list<Interval>::iterator j = intervals.begin(); j != intervals.end();)
- if (int_mod((*j).modulo, (*i).modulo) == 0 &&
- (*j).ub - (*j).lb >= (*i).size &&
- (int_mod((*j).lb, (*i).modulo) == (*i).start ||
- int_mod((*j).ub, (*i).modulo) == (*i).start + (*i).size)) {
- if (int_mod((*j).lb, (*i).modulo) == (*i).start &&
- (*j).ub - (*j).lb == (*i).size) {
- if (replaced_one == intervals.end()) {
- (*(*j).pos).first = hull;
- (*j).lb = int_mod((*j).lb, (*i).modulo);
- (*j).ub = int_mod((*j).ub, (*i).modulo);
- (*j).modulo = (*i).modulo;
- (*j).change = true;
- replaced_one = j;
- j++;
- }
- else {
- if ((*j).pos == orig) {
- std::swap((*replaced_one).pos, (*j).pos);
- (*(*replaced_one).pos).first = (*(*j).pos).first;
- }
- Rs.erase((*j).pos);
- j = intervals.erase(j);
- }
- }
- else {
- if (int_mod((*j).lb, (*i).modulo) == (*i).start)
- (*j).lb = (*j).lb + (*i).size + 1;
- else
- (*j).ub = (*j).ub - (*i).size - 1;
- (*j).change = true;
- j++;
- }
- }
- else
- j++;
- }
- }
-}
-} // namespace
-
-
-//
-// Simplify a union of sets/relations to a minimal (may not be
-// optimal) number of convex regions. It intends to replace
-// CheckForConvexRepresentation and CheckForConvexPairs functions.
-//
-Relation ConvexRepresentation(NOT_CONST Relation &R) {
- Relation l_R = copy(R);
- if (!l_R.is_upper_bound_satisfiable() || l_R.number_of_conjuncts() < 2)
- return R;
-
- // separate each conjunct into smooth convex region and holes
- std::list<std::pair<Relation, Relation> > Rs; // pair(smooth region, hole condition)
- for (DNF_Iterator c(l_R.query_DNF()); c.live(); c++) {
- Relation r1 = Relation(l_R, c.curr());
- Relation r2 = Approximate(copy(r1));
- r1 = Gist(r1, copy(r2));
- Rs.push_back(std::make_pair(r2, r1));
- }
-
- try {
- bool change = true;
- while (change) {
- change = false;
-
- std::list<std::pair<Relation, Relation> >::iterator i = Rs.begin();
- while (i != Rs.end()) {
- // find regions with identical hole conditions to merge
- {
- std::list<std::pair<Relation, Relation> >::iterator j = i;
- j++;
- while (j != Rs.end()) {
- if (!Difference(copy((*i).second), copy((*j).second)).is_upper_bound_satisfiable() &&
- !Difference(copy((*j).second), copy((*i).second)).is_upper_bound_satisfiable()) {
- if (Must_Be_Subset(copy((*j).first), copy((*i).first))) {
- j = Rs.erase(j);
- }
- else if (Must_Be_Subset(copy((*i).first), copy((*j).first))) {
- (*i).first = (*j).first;
- j = Rs.erase(j);
- change = true;
- }
- else {
- Relation r;
- bool already_use_recthull = false;
- try {
- // chun's debug
- // throw std::runtime_error("dfdf");
-
- r = ConvexHull(Union(copy((*i).first), copy((*j).first)));
- }
- catch (const std::overflow_error &e) {
- r = RectHull(Union(copy((*i).first), copy((*j).first)));
- already_use_recthull = true;
- }
- retry_recthull:
- Relation r2 = Difference(Difference(copy(r), copy((*i).first)), copy((*j).first));
- if (!r2.is_upper_bound_satisfiable()) { // convex hull is tight
- (*i).first = r;
- j = Rs.erase(j);
- change = true;
- }
- else {
- if (!already_use_recthull) {
- r = RectHull(Union(copy((*i).first), copy((*j).first)));
- already_use_recthull = true;
- goto retry_recthull;
- }
- else
- j++;
- }
- }
- }
- else
- j++;
- }
- }
-
- // find identical smooth regions as candidates for hole merge
- std::list<std::list<std::pair<Relation, Relation> >::iterator> s;
- for (std::list<std::pair<Relation, Relation> >::iterator j = Rs.begin(); j != Rs.end(); j++)
- if (j != i) {
- if (!Intersection(Difference(copy((*i).first), copy((*j).first)), copy((*j).second)).is_upper_bound_satisfiable() &&
- !Intersection(Difference(copy((*j).first), copy((*i).first)), copy((*i).second)).is_upper_bound_satisfiable())
- s.push_back(j);
- }
-
- if (s.size() != 0) {
- // convert hole condition c1*x1+c2*x2+... = c*alpha+d to a pair of inequalities
- (*i).second = EQs_to_GEQs((*i).second, false);
-
- // find potential wildcards that can be used for hole conditions
- std::set<Variable_ID> nonsingle_wild;
- for (EQ_Iterator ei((*i).second.single_conjunct()); ei; ei++)
- if ((*ei).has_wildcards())
- for (Constr_Vars_Iter cvi(*ei, true); cvi; cvi++)
- nonsingle_wild.insert(cvi.curr_var());
- for (GEQ_Iterator gei((*i).second.single_conjunct()); gei; gei++)
- if ((*gei).has_wildcards()) {
- Constr_Vars_Iter cvi(*gei, true);
- Constr_Vars_Iter cvi2 = cvi;
- cvi2++;
- if (cvi2) {
- nonsingle_wild.insert(cvi.curr_var());
- for (; cvi2; cvi2++)
- nonsingle_wild.insert(cvi2.curr_var());
- }
- }
-
- // find hole condition in c*alpha+d1<=c1*x1+c2*x2+...<=c*alpha+d2 format
- for (GEQ_Iterator gei((*i).second.single_conjunct()); gei; gei++)
- if ((*gei).has_wildcards()) {
- coef_t c;
- Variable_ID v;
- {
- Constr_Vars_Iter cvi(*gei, true);
- v = cvi.curr_var();
- c = cvi.curr_coef();
- if (c < 0 || nonsingle_wild.find(v) != nonsingle_wild.end())
- continue;
- }
-
- coef_t lb = posInfinity;
- for (GEQ_Iterator gei2((*i).second.single_conjunct()); gei2; gei2++) {
- if (!(*gei2 == *gei) && (*gei2).get_coef(v) != 0) {
- if (lb != posInfinity) {
- nonsingle_wild.insert(v);
- break;
- }
-
- bool match = true;
- for (Constr_Vars_Iter cvi2(*gei); cvi2; cvi2++)
- if (cvi2.curr_coef() != -((*gei2).get_coef(cvi2.curr_var()))) {
- match = false;
- break;
- }
- if (match)
- for (Constr_Vars_Iter cvi2(*gei2); cvi2; cvi2++)
- if (cvi2.curr_coef() != -((*gei).get_coef(cvi2.curr_var()))) {
- match = false;
- break;
- }
- if (!match) {
- nonsingle_wild.insert(v);
- break;
- }
-
- lb = -(*gei2).get_const();
- }
- }
-
- if (nonsingle_wild.find(v) != nonsingle_wild.end())
- continue;
-
- Relation stride_cond = Relation::True((*i).second);
- F_Exists *f_exists = stride_cond.and_with_and()->add_exists();
- Variable_ID e = f_exists->declare();
- F_And *f_root = f_exists->add_and();
- GEQ_Handle h1 = f_root->add_GEQ();
- GEQ_Handle h2 = f_root->add_GEQ();
- for (Constr_Vars_Iter cvi2(*gei); cvi2; cvi2++) {
- Variable_ID v = cvi2.curr_var();
- switch (v->kind()) {
- case Wildcard_Var:
- h1.update_coef(e, cvi2.curr_coef());
- h2.update_coef(e, -cvi2.curr_coef());
- break;
- case Global_Var: {
- Global_Var_ID g = v->get_global_var();
- Variable_ID v2;
- if (g->arity() == 0)
- v2 = stride_cond.get_local(g);
- else
- v2 = stride_cond.get_local(g, v->function_of());
- h1.update_coef(v2, cvi2.curr_coef());
- h2.update_coef(v2, -cvi2.curr_coef());
- break;
- }
- default:
- h1.update_coef(v, cvi2.curr_coef());
- h2.update_coef(v, -cvi2.curr_coef());
- }
- }
- h1.update_const((*gei).get_const());
- h2.update_const(-lb);
-
- stride_cond.simplify();
- Relation other_cond = Gist(copy((*i).second), copy(stride_cond));
-
- // find regions with potential mergeable stride condition with this one
- std::list<Interval> intervals;
- intervals.push_back(Interval(i, lb, (*gei).get_const()));
-
- for (std::list<std::list<std::pair<Relation, Relation> >::iterator>::iterator j = s.begin(); j != s.end(); j++)
- if (Must_Be_Subset(copy((**j).second), copy(other_cond))) {
- Relation stride_cond2 = Gist(copy((**j).second), copy(other_cond));
-
- // interval can be removed
- if (stride_cond2.is_obvious_tautology()) {
- intervals.push_back(Interval(*j, 0, c-1));
- continue;
- }
-
- stride_cond2 = EQs_to_GEQs(stride_cond2, false);
- coef_t lb, ub;
- GEQ_Iterator gei2(stride_cond2.single_conjunct());
- coef_t sign = 0;
- for (Constr_Vars_Iter cvi(*gei2, true); cvi; cvi++)
- if (sign != 0) {
- sign = 0;
- break;
- }
- else if (cvi.curr_coef() == c)
- sign = 1;
- else if (cvi.curr_coef() == -c)
- sign = -1;
- else {
- sign = 0;
- break;
- }
- if (sign == 0)
- continue;
-
- bool match = true;
- for (Constr_Vars_Iter cvi(*gei2); cvi; cvi++) {
- Variable_ID v = cvi.curr_var();
- if (v->kind() == Wildcard_Var)
- continue;
- else if (v->kind() == Global_Var) {
- Global_Var_ID g = v->get_global_var();
- if (g->arity() == 0)
- v = (*i).second.get_local(g);
- else
- v = (*i).second.get_local(g, v->function_of());
- }
-
- if (cvi.curr_coef() != sign * (*gei).get_coef(v)) {
- match = false;
- break;
- }
- }
- if (!match)
- continue;
-
- for (Constr_Vars_Iter cvi(*gei); cvi; cvi++) {
- Variable_ID v = cvi.curr_var();
- if (v->kind() == Wildcard_Var)
- continue;
- else if (v->kind() == Global_Var) {
- Global_Var_ID g = v->get_global_var();
- if (g->arity() == 0)
- v = stride_cond2.get_local(g);
- else
- v = stride_cond2.get_local(g, v->function_of());
- }
-
- if (cvi.curr_coef() != sign * (*gei2).get_coef(v)) {
- match = false;
- break;
- }
- }
- if (!match)
- continue;
- if (sign > 0)
- ub = (*gei2).get_const();
- else
- lb = -(*gei2).get_const();
-
- gei2++;
- if (!gei2)
- continue;
-
- coef_t sign2 = 0;
- for (Constr_Vars_Iter cvi(*gei2, true); cvi; cvi++)
- if (sign2 != 0) {
- sign2 = 0;
- break;
- }
- else if (cvi.curr_coef() == c)
- sign2 = 1;
- else if (cvi.curr_coef() == -c)
- sign2 = -1;
- else {
- sign2 = 0;
- break;
- }
- if (sign2 != -sign)
- continue;
-
- for (Constr_Vars_Iter cvi(*gei2); cvi; cvi++) {
- Variable_ID v = cvi.curr_var();
- if (v->kind() == Wildcard_Var)
- continue;
- else if (v->kind() == Global_Var) {
- Global_Var_ID g = v->get_global_var();
- if (g->arity() == 0)
- v = (*i).second.get_local(g);
- else
- v = (*i).second.get_local(g, v->function_of());
- }
-
- if (cvi.curr_coef() != sign2 * (*gei).get_coef(v)) {
- match = false;
- break;
- }
- }
- if (!match)
- continue;
-
- for (Constr_Vars_Iter cvi(*gei); cvi; cvi++) {
- Variable_ID v = cvi.curr_var();
- if (v->kind() == Wildcard_Var)
- continue;
- else if (v->kind() == Global_Var) {
- Global_Var_ID g = v->get_global_var();
- if (g->arity() == 0)
- v = stride_cond2.get_local(g);
- else
- v = stride_cond2.get_local(g, v->function_of());
- }
-
- if (cvi.curr_coef() != sign2 * (*gei2).get_coef(v)) {
- match = false;
- break;
- }
- }
- if (!match)
- continue;
- if (sign2 > 0)
- ub = (*gei2).get_const();
- else
- lb = -(*gei2).get_const();
-
- gei2++;
- if (gei2)
- continue;
-
- intervals.push_back(Interval(*j, lb, ub));
- }
-
- merge_intervals(intervals, c, Rs, i);
-
- // make current region the last one being updated
- bool invalid = false;
- for (std::list<Interval>::iterator ii = intervals.begin(); ii != intervals.end(); ii++)
- if ((*ii).change && (*ii).pos == i) {
- invalid = true;
- intervals.push_back(*ii);
- intervals.erase(ii);
- break;
- }
-
- // update hole condition for each region
- for (std::list<Interval>::iterator ii = intervals.begin(); ii != intervals.end(); ii++)
- if ((*ii).change) {
- change = true;
-
- if ((*ii).ub - (*ii).lb + 1 >= (*ii).modulo)
- (*(*ii).pos).second = copy(other_cond);
- else {
- Relation stride_cond = Relation::True((*i).second);
- F_Exists *f_exists = stride_cond.and_with_and()->add_exists();
- Variable_ID e = f_exists->declare();
- F_And *f_root = f_exists->add_and();
- GEQ_Handle h1 = f_root->add_GEQ();
- GEQ_Handle h2 = f_root->add_GEQ();
- for (Constr_Vars_Iter cvi2(*gei); cvi2; cvi2++) {
- Variable_ID v = cvi2.curr_var();
- switch (v->kind()) {
- case Wildcard_Var:
- h1.update_coef(e, (*ii).modulo);
- h2.update_coef(e, -(*ii).modulo);
- break;
- case Global_Var: {
- Global_Var_ID g = v->get_global_var();
- Variable_ID v2;
- if (g->arity() == 0)
- v2 = stride_cond.get_local(g);
- else
- v2 = stride_cond.get_local(g, v->function_of());
- h1.update_coef(v2, cvi2.curr_coef());
- h2.update_coef(v2, -cvi2.curr_coef());
- break;
- }
- default:
- h1.update_coef(v, cvi2.curr_coef());
- h2.update_coef(v, -cvi2.curr_coef());
- }
- }
- h1.update_const((*ii).ub);
- h2.update_const(-(*ii).lb);
-
- (*(*ii).pos).second = Intersection(copy(other_cond), stride_cond);
- (*(*ii).pos).second.simplify();
- }
- }
-
- if (invalid)
- break;
- }
- }
- i++;
- }
- }
- }
- catch (const presburger_error &e) {
- throw e;
- }
-
- Relation R2 = Relation::False(l_R);
- for (std::list<std::pair<Relation, Relation> >::iterator i = Rs.begin(); i != Rs.end(); i++)
- R2 = Union(R2, Intersection((*i).first, (*i).second));
- R2.simplify(0, 1);
-
- return R2;
-}
-
-
-//
-// Use gist and value range to calculate a quick rectangular hull. It
-// intends to replace all hull calculations (QuickHull, BetterHull,
-// FastTightHull) beyond the method of ConvexHull (dual
-// representations). In the future, it will support max(...)-like
-// upper bound. So RectHull complements ConvexHull in two ways: first
-// for relations that ConvexHull gets too complicated, second for
-// relations where different conjuncts have different symbolic upper
-// bounds.
-//
-Relation RectHull(NOT_CONST Relation &Rel) {
- Relation R = Approximate(consume_and_regurgitate(Rel));
- if (!R.is_upper_bound_satisfiable())
- return R;
- if (R.has_single_conjunct())
- return R;
-
- std::vector<std::string> input_names(R.n_inp());
- for (int i = 1; i <= R.n_inp(); i++)
- input_names[i-1] = R.input_var(i)->name();
- std::vector<std::string> output_names(R.n_out());
- for (int i = 1; i <= R.n_out(); i++)
- output_names[i-1] = R.output_var(i)->name();
-
- DNF_Iterator c(R.query_DNF());
- Relation r = Relation(R, c.curr());
- c++;
- std::vector<std::pair<coef_t, coef_t> > bounds1(R.n_inp());
- std::vector<std::pair<coef_t, coef_t> > bounds2(R.n_out());
- {
- Relation t = Project_Sym(copy(r));
- t.simplify();
- for (int i = 1; i <= R.n_inp(); i++) {
- Tuple<Variable_ID> v;
- for (int j = 1; j <= R.n_inp(); j++)
- if (j != i)
- v.append(r.input_var(j));
- for (int j = 1; j <= R.n_out(); j++)
- v.append(r.output_var(j));
- Relation t2 = Project(copy(t), v);
- t2.query_variable_bounds(t2.input_var(i), bounds1[i-1].first, bounds1[i-1].second);
- }
- for (int i = 1; i <= R.n_out(); i++) {
- Tuple<Variable_ID> v;
- for (int j = 1; j <= R.n_out(); j++)
- if (j != i)
- v.append(r.output_var(j));
- for (int j = 1; j <= R.n_inp(); j++)
- v.append(r.input_var(j));
- Relation t2 = Project(copy(t), v);
- t2.query_variable_bounds(t2.output_var(i), bounds2[i-1].first, bounds2[i-1].second);
- }
- }
-
- while (c.live()) {
- Relation r2 = Relation(R, c.curr());
- c++;
- Relation x = Gist(copy(r), Gist(copy(r), copy(r2), 1), 1);
- if (Difference(copy(r2), copy(x)).is_upper_bound_satisfiable())
- x = Relation::True(R);
- Relation y = Gist(copy(r2), Gist(copy(r2), copy(r), 1), 1);
- if (Difference(copy(r), copy(y)).is_upper_bound_satisfiable())
- y = Relation::True(R);
- r = Intersection(x, y);
-
- {
- Relation t = Project_Sym(copy(r2));
- t.simplify();
- for (int i = 1; i <= R.n_inp(); i++) {
- Tuple<Variable_ID> v;
- for (int j = 1; j <= R.n_inp(); j++)
- if (j != i)
- v.append(r2.input_var(j));
- for (int j = 1; j <= R.n_out(); j++)
- v.append(r2.output_var(j));
- Relation t2 = Project(copy(t), v);
- coef_t lbound, ubound;
- t2.query_variable_bounds(t2.input_var(i), lbound, ubound);
- bounds1[i-1].first = min(bounds1[i-1].first, lbound);
- bounds1[i-1].second = max(bounds1[i-1].second, ubound);
- }
- for (int i = 1; i <= R.n_out(); i++) {
- Tuple<Variable_ID> v;
- for (int j = 1; j <= R.n_out(); j++)
- if (j != i)
- v.append(r2.output_var(j));
- for (int j = 1; j <= R.n_inp(); j++)
- v.append(r2.input_var(j));
- Relation t2 = Project(copy(t), v);
- coef_t lbound, ubound;
- t2.query_variable_bounds(t2.output_var(i), lbound, ubound);
- bounds2[i-1].first = min(bounds2[i-1].first, lbound);
- bounds2[i-1].second = max(bounds2[i-1].second, ubound);
- }
- }
-
- Relation r3(R.n_inp(), R.n_out());
- F_And *f_root = r3.add_and();
- for (int i = 1; i <= R.n_inp(); i++) {
- if (bounds1[i-1].first != -posInfinity) {
- GEQ_Handle h = f_root->add_GEQ();
- h.update_coef(r3.input_var(i), 1);
- h.update_const(-bounds1[i-1].first);
- }
- if (bounds1[i-1].second != posInfinity) {
- GEQ_Handle h = f_root->add_GEQ();
- h.update_coef(r3.input_var(i), -1);
- h.update_const(bounds1[i-1].second);
- }
- }
- for (int i = 1; i <= R.n_out(); i++) {
- if (bounds2[i-1].first != -posInfinity) {
- GEQ_Handle h = f_root->add_GEQ();
- h.update_coef(r3.output_var(i), 1);
- h.update_const(-bounds2[i-1].first);
- }
- if (bounds2[i-1].second != posInfinity) {
- GEQ_Handle h = f_root->add_GEQ();
- h.update_coef(r3.output_var(i), -1);
- h.update_const(bounds2[i-1].second);
- }
- }
- r = Intersection(r, r3);
- r.simplify();
- }
-
- for (int i = 1; i <= r.n_inp(); i++)
- r.name_input_var(i, input_names[i-1]);
- for (int i = 1; i <= r.n_out(); i++)
- r.name_output_var(i, output_names[i-1]);
- r.setup_names();
- return r;
-}
-
-} // namespace