summaryrefslogtreecommitdiff
path: root/omegalib/omega/src/omega_core
diff options
context:
space:
mode:
Diffstat (limited to 'omegalib/omega/src/omega_core')
-rw-r--r--omegalib/omega/src/omega_core/oc.cc754
-rw-r--r--omegalib/omega/src/omega_core/oc_eq.cc653
-rw-r--r--omegalib/omega/src/omega_core/oc_exp_kill.cc297
-rw-r--r--omegalib/omega/src/omega_core/oc_global.cc45
-rw-r--r--omegalib/omega/src/omega_core/oc_print.cc686
-rw-r--r--omegalib/omega/src/omega_core/oc_problems.cc198
-rw-r--r--omegalib/omega/src/omega_core/oc_query.cc478
-rw-r--r--omegalib/omega/src/omega_core/oc_quick_kill.cc775
-rw-r--r--omegalib/omega/src/omega_core/oc_simple.cc1373
-rw-r--r--omegalib/omega/src/omega_core/oc_solve.cc1378
-rw-r--r--omegalib/omega/src/omega_core/oc_util.cc327
11 files changed, 6964 insertions, 0 deletions
diff --git a/omegalib/omega/src/omega_core/oc.cc b/omegalib/omega/src/omega_core/oc.cc
new file mode 100644
index 0000000..0dc9b49
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc.cc
@@ -0,0 +1,754 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+ Simplify a problem.
+
+ Notes:
+
+ History:
+ 12/10/06 Improved gist function, by Chun Chen.
+*****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+
+namespace omega {
+
+eqn SUBs[maxVars+1];
+Memory redMemory[maxVars+1];
+
+int Problem::reduceProblem() {
+ int result;
+ checkVars(nVars+1);
+ assert(omegaInitialized);
+ if (nVars > nEQs + 3 * safeVars)
+ freeEliminations(safeVars);
+
+ check();
+ if (!mayBeRed && nSUBs == 0 && safeVars == 0) {
+ result = solve(OC_SOLVE_UNKNOWN);
+ nGEQs = 0;
+ nEQs = 0;
+ nSUBs = 0;
+ nMemories = 0;
+ if (!result) {
+ int e = newEQ();
+ assert(e == 0);
+ eqnnzero(&EQs[0], nVars);
+ EQs[0].color = EQ_BLACK;
+ EQs[0].coef[0] = 1;
+ }
+ check();
+ return result;
+ }
+ return solve(OC_SOLVE_SIMPLIFY);
+}
+
+
+int Problem::simplifyProblem(int verify, int subs, int redundantElimination) {
+ checkVars(nVars+1);
+ assert(omegaInitialized);
+ setInternals();
+ check();
+ if (!reduceProblem()) goto returnFalse;
+ if (verify) {
+ addingOuterEqualities++;
+ int r = verifyProblem();
+ addingOuterEqualities--;
+ if (!r) goto returnFalse;
+ if (nEQs) { // found some equality constraints during verification
+ int numRed = 0;
+ if (mayBeRed)
+ for (int e = nGEQs - 1; e >= 0; e--) if (GEQs[e].color == EQ_RED) numRed++;
+ if (mayBeRed && nVars == safeVars && numRed == 1)
+ nEQs = 0; // discard them
+ else if (!reduceProblem()) {
+ assert(0 && "Added equality constraint to verified problem generates false");
+ }
+ }
+ }
+ if (redundantElimination) {
+ if (redundantElimination > 1) {
+ if (!expensiveEqualityCheck()) goto returnFalse;
+ }
+ if (!quickKill(0)) goto returnFalse;
+ if (redundantElimination > 1) {
+ if (!expensiveKill()) goto returnFalse;
+ }
+ }
+ resurrectSubs();
+ if (redundantElimination)
+ simplifyStrideConstraints();
+ if (redundantElimination > 2 && safeVars < nVars) {
+ if (!quickKill(0)) goto returnFalse;
+ return simplifyProblem(verify, subs, redundantElimination-2);
+ }
+ setExternals();
+ assert(nMemories == 0);
+ return (1);
+
+returnFalse:
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nGEQs = 0;
+ nEQs = 0;
+ int neweq = newEQ();
+ assert(neweq == 0);
+ eqnnzero(&EQs[neweq], nVars);
+ EQs[neweq].color = EQ_BLACK;
+ EQs[neweq].coef[0] = 1;
+ nMemories = 0;
+ return 0;
+}
+
+int Problem::simplifyApproximate(bool strides_allowed) {
+ int result;
+ checkVars(nVars+1);
+ assert(inApproximateMode == 0);
+
+ inApproximateMode = 1;
+ inStridesAllowedMode = strides_allowed;
+ if (TRACE)
+ fprintf(outputFile, "Entering Approximate Mode [\n");
+
+ assert(omegaInitialized);
+ result = simplifyProblem(0,0,0);
+
+ while (result && !strides_allowed && nVars > safeVars) {
+ int e;
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[nVars]) deleteGEQ(e);
+ for (e = nEQs - 1; e >= 0; e--)
+ if (EQs[e].coef[nVars]) deleteEQ(e);
+ nVars--;
+ result = simplifyProblem(0,0,0);
+ }
+
+ if (TRACE)
+ fprintf(outputFile, "] Leaving Approximate Mode\n");
+
+ assert(inApproximateMode == 1);
+ inApproximateMode=0;
+ inStridesAllowedMode = 0;
+
+ assert(nMemories == 0);
+ return (result);
+}
+
+
+
+
+/*
+ * Return 1 if red equations constrain the set of possible
+ * solutions. We assume that there are solutions to the black
+ * equations by themselves, so if there is no solution to the combined
+ * problem, we return 1.
+ */
+
+#ifdef GIST_CHECK
+Problem full_answer, context;
+Problem redProblem;
+#endif
+
+redCheck Problem::redSimplifyProblem(int effort, int computeGist) {
+ int result;
+ int e;
+
+ checkVars(nVars+1);
+ assert(mayBeRed >= 0);
+ mayBeRed++;
+
+ assert(omegaInitialized);
+ if (TRACE) {
+ fprintf(outputFile, "Checking for red equations:\n");
+ printProblem();
+ }
+ setInternals();
+
+#ifdef GIST_CHECK
+ int r1,r2;
+ if (TRACE)
+ fprintf(outputFile,"Set-up for gist invariant checking[\n");
+ redProblem = *this;
+ redProblem.check();
+ full_answer = *this;
+ full_answer.check();
+ full_answer.turnRedBlack();
+ full_answer.check();
+ r1 = full_answer.simplifyProblem(1,0,1);
+ full_answer.check();
+ if (DBUG) fprintf(outputFile,"Simplifying context [\n");
+ context = *this;
+ context.check();
+ context.deleteRed();
+ context.check();
+ r2 = context.simplifyProblem(1,0,1);
+ context.check();
+ if (DBUG) fprintf(outputFile,"] Simplifying context\n");
+
+ if (!r2 && TRACE) fprintf(outputFile, "WARNING: Gist context is false!\n");
+ if (TRACE)
+ fprintf(outputFile,"] Set-up for gist invariant checking done\n");
+#endif
+
+ // Save known integer modular equations, -- by chun 12/10/2006
+ eqn ModularEQs[nEQs];
+ int nModularEQs = 0;
+ int old_nVars = nVars;
+ for (int e = 0; e < nEQs; e++)
+ if (EQs[e].color != EQ_RED)
+ for (int i = safeVars+1; i <= nVars; i++)
+ if (EQs[e].coef[i] != 0) {
+ eqnncpy(&(ModularEQs[nModularEQs++]), &(EQs[e]), nVars);
+ break;
+ }
+
+
+ if (solveEQ() == false) {
+ if (TRACE)
+ fprintf(outputFile, "Gist is FALSE\n");
+ if (computeGist) {
+ nMemories = 0;
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nGEQs = 0;
+ nEQs = 0;
+ int neweq = newEQ();
+ assert(neweq == 0);
+ eqnnzero(&EQs[neweq], nVars);
+ EQs[neweq].color = EQ_RED;
+ EQs[neweq].coef[0] = 1;
+ }
+ mayBeRed--;
+ return redFalse;
+ }
+
+ if (!computeGist && nMemories)
+ return redConstraints;
+ if (normalize() == normalize_false) {
+ if (TRACE)
+ fprintf(outputFile, "Gist is FALSE\n");
+ if (computeGist) {
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nMemories = 0;
+ nGEQs = 0;
+ nEQs = 0;
+ int neweq = newEQ();
+ assert(neweq == 0);
+ eqnnzero(&EQs[neweq], nVars);
+ EQs[neweq].color = EQ_RED;
+ EQs[neweq].coef[0] = 1;
+ }
+ mayBeRed--;
+ return redFalse;
+ }
+
+ result = 0;
+ for (e = nGEQs - 1; e >= 0; e--) if (GEQs[e].color == EQ_RED) result = 1;
+ for (e = nEQs - 1; e >= 0; e--) if (EQs[e].color == EQ_RED) result = 1;
+ if (nMemories) result = 1;
+ if (!result) {
+ if (computeGist) {
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nGEQs = 0;
+ nMemories = 0;
+ nEQs = 0;
+ }
+ mayBeRed--;
+ return noRed;
+ }
+
+ result = simplifyProblem(effort?1:0,1,0);
+#ifdef GIST_CHECK
+ if (!r1 && TRACE && result)
+ fprintf(outputFile, "Gist is False but not detected\n");
+#endif
+ if (!result) {
+ if (TRACE)
+ fprintf(outputFile, "Gist is FALSE\n");
+ if (computeGist) {
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nGEQs = 0;
+ nEQs = 0;
+ int neweq = newEQ();
+ assert(neweq == 0);
+ nMemories = 0;
+ eqnnzero(&EQs[neweq], nVars);
+ EQs[neweq].color = EQ_RED;
+ EQs[neweq].coef[0] = 1;
+ }
+ mayBeRed--;
+ return redFalse;
+ }
+
+ freeRedEliminations();
+
+ result = 0;
+ for (e = nGEQs - 1; e >= 0; e--) if (GEQs[e].color == EQ_RED) result = 1;
+ for (e = nEQs - 1; e >= 0; e--) if (EQs[e].color == EQ_RED) result = 1;
+ if (nMemories) result = 1;
+ if (!result) {
+ if (computeGist) {
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nGEQs = 0;
+ nMemories = 0;
+ nEQs = 0;
+ }
+ mayBeRed--;
+ return noRed;
+ }
+
+ if (effort && (computeGist || !nMemories)) {
+ if (TRACE)
+ fprintf(outputFile, "*** Doing potentially expensive elimination tests for red equations [\n");
+ quickRedKill(computeGist);
+ checkGistInvariant();
+ result = nMemories;
+ for (e = nGEQs - 1; e >= 0; e--) if (GEQs[e].color == EQ_RED) result++;
+ for (e = nEQs - 1; e >= 0; e--) if (EQs[e].color == EQ_RED) result++;
+ if (result && effort > 1 && (computeGist || !nMemories)) {
+ expensiveRedKill();
+ result = nMemories;
+ for (e = nGEQs-1; e >= 0; e--) if (GEQs[e].color == EQ_RED) result++;
+ for (e = nEQs-1; e >= 0; e--) if (EQs[e].color == EQ_RED) result++;
+ }
+
+ if (!result) {
+ if (TRACE)
+ fprintf(outputFile, "]******************** Redudant Red Equations eliminated!!\n");
+ if (computeGist) {
+ nGEQs = 0;
+ nEQs = 0;
+ resurrectSubs();
+ nGEQs = 0;
+ nMemories = 0;
+ nEQs = 0;
+ }
+ mayBeRed--;
+ return noRed;
+ }
+
+ if (TRACE) fprintf(outputFile, "]******************** Red Equations remain\n");
+ if (DEBUG) printProblem();
+ }
+
+ if (computeGist) {
+ resurrectSubs();
+ cleanoutWildcards();
+
+
+ // Restore saved modular equations into EQs without affecting the problem
+ if (nEQs+nModularEQs > allocEQs) {
+ allocEQs = padEQs(allocEQs, nEQs+nModularEQs);
+ eqn *new_eqs = new eqn[allocEQs];
+ for (int e = 0; e < nEQs; e++)
+ eqnncpy(&(new_eqs[e]), &(EQs[e]), nVars);
+ delete[] EQs;
+ EQs= new_eqs;
+ }
+
+ for (int e = 0; e < nModularEQs; e++) {
+ eqnncpy(&(EQs[nEQs+e]), &(ModularEQs[e]), old_nVars);
+ EQs[nEQs+e].color = EQ_RED;
+ Tuple<coef_t> t(safeVars);
+ for (int i = 1; i <= safeVars; i++)
+ t[i] = ModularEQs[e].coef[var[i]];
+ for (int i = 1; i <= safeVars; i++)
+ EQs[nEQs+e].coef[i] = t[i];
+ }
+
+
+ // Now simplify modular equations using Chinese remainder theorem -- by chun 12/10/2006
+ for (int e = 0; e < nEQs; e++)
+ if (EQs[e].color == EQ_RED) {
+ int wild_pos = -1;
+ for (int i = safeVars+1; i <= nVars; i++)
+ if (EQs[e].coef[i] != 0) {
+ wild_pos = i;
+ break;
+ }
+
+ if (wild_pos == -1)
+ continue;
+
+ for (int e2 = e+1; e2 < nEQs+nModularEQs; e2++)
+ if (EQs[e2].color == EQ_RED) {
+ int wild_pos2 = -1;
+ for (int i = safeVars+1; i <= ((e2<nEQs)?nVars:old_nVars); i++)
+ if (EQs[e2].coef[i] != 0) {
+ wild_pos2 = i;
+ break;
+ }
+
+ if (wild_pos2 == -1)
+ continue;
+
+ coef_t g = gcd(abs(EQs[e].coef[wild_pos]), abs(EQs[e2].coef[wild_pos2]));
+ coef_t g2 = 1;
+ coef_t g3;
+ EQs[e].color = EQs[e2].color = EQ_BLACK;
+ while ((g3 = factor(g)) != 1) {
+ coef_t gg = g2 * g3;
+ g = g/g3;
+
+ bool match = true;
+ coef_t c = EQs[e].coef[0];
+ coef_t c2 = EQs[e2].coef[0];
+ bool change_sign = false;
+ for (int i = 1; i <= safeVars; i++) {
+ coef_t coef = int_mod_hat(EQs[e].coef[i], gg);
+ coef_t coef2 = int_mod_hat(EQs[e2].coef[i], gg);
+
+ if (coef == 0 && coef2 == 0)
+ continue;
+
+ if (change_sign && coef == -coef2)
+ continue;
+
+ if (!change_sign) {
+ if (coef == coef2)
+ continue;
+ else if (coef == -coef2) {
+ change_sign = true;
+ continue;
+ }
+ }
+
+ if (coef != 0) {
+ coef_t t = query_variable_mod(i, gg/gcd(abs(coef), gg), EQ_RED, nModularEQs, old_nVars);
+ if (t == posInfinity) {
+ match = false;
+ break;
+ }
+
+ c += coef * t;
+ }
+ if (coef2 != 0) {
+ coef_t t = query_variable_mod(i, gg/gcd(abs(coef2), gg), EQ_RED, nModularEQs, old_nVars);
+ if (t == posInfinity) {
+ match = false;
+ break;
+ }
+
+ c2 += coef2 * t;
+ }
+ }
+ if ((change_sign && int_mod_hat(c, gg) != -int_mod_hat(c2, gg)) ||
+ (!change_sign && int_mod_hat(c, gg) != int_mod_hat(c2, gg)))
+ match = false;
+
+ if (match)
+ g2 = gg;
+ }
+ EQs[e].color = EQs[e2].color = EQ_RED;
+
+ if (g2 == 1)
+ continue;
+
+ if (g2 == abs(EQs[e].coef[wild_pos])) {
+ EQs[e].color = EQ_BLACK;
+ break;
+ }
+ else if (e2 < nEQs && g2 == abs(EQs[e2].coef[wild_pos2]))
+ EQs[e2].color = EQ_BLACK;
+ else {
+ coef_t g4 = abs(EQs[e].coef[wild_pos])/g2;
+ while (lcm(g2, g4) != abs(EQs[e].coef[wild_pos])) {
+ assert(lcm(g2, g4) < abs(EQs[e].coef[wild_pos]));
+ g4 *= abs(EQs[e].coef[wild_pos])/lcm(g2, g4);
+ }
+
+ for (int i = 0; i <= safeVars; i++)
+ EQs[e].coef[i] = int_mod_hat(EQs[e].coef[i], g4);
+ EQs[e].coef[wild_pos] = (EQs[e].coef[wild_pos]>0?1:-1)*g4;
+ }
+ }
+ }
+
+ deleteBlack();
+ }
+
+ setExternals();
+ mayBeRed--;
+ assert(nMemories == 0);
+ return redConstraints;
+}
+
+
+void Problem::convertEQstoGEQs(bool excludeStrides) {
+ int i;
+ int e;
+ if (DBUG)
+ fprintf(outputFile, "Converting all EQs to GEQs\n");
+ simplifyProblem(0,0,0);
+ for(e=0;e<nEQs;e++) {
+ bool isStride = 0;
+ int e2 = newGEQ();
+ if (excludeStrides)
+ for(i = safeVars+1; i <= nVars; i++)
+ isStride = isStride || (EQs[e].coef[i] != 0);
+ if (isStride) continue;
+ eqnncpy(&GEQs[e2], &EQs[e], nVars);
+ GEQs[e2].touched = 1;
+ e2 = newGEQ();
+ eqnncpy(&GEQs[e2], &EQs[e], nVars);
+ GEQs[e2].touched = 1;
+ for (i = 0; i <= nVars; i++)
+ GEQs[e2].coef[i] = -GEQs[e2].coef[i];
+ }
+ // If we have eliminated all EQs, can set nEQs to 0
+ // If some strides are left, we don't know the position of them in the EQs
+ // array, so decreasing nEQs might remove wrong EQs -- we just leave them
+ // all in. (could sort the EQs to move strides to the front, but too hard.)
+ if (!excludeStrides) nEQs=0;
+ if (DBUG)
+ printProblem();
+}
+
+
+void Problem::convertEQtoGEQs(int eq) {
+ int i;
+ if (DBUG)
+ fprintf(outputFile, "Converting EQ %d to GEQs\n",eq);
+ int e2 = newGEQ();
+ eqnncpy(&GEQs[e2], &EQs[eq], nVars);
+ GEQs[e2].touched = 1;
+ e2 = newGEQ();
+ eqnncpy(&GEQs[e2], &EQs[eq], nVars);
+ GEQs[e2].touched = 1;
+ for (i = 0; i <= nVars; i++)
+ GEQs[e2].coef[i] = -GEQs[e2].coef[i];
+ if (DBUG)
+ printProblem();
+}
+
+
+/*
+ * Calculate value of variable modulo integer from problem's equation
+ * set plus additional saved modular equations embedded in the same
+ * EQs array (hinted by nModularEQs) if available. If there is no
+ * solution, return posInfinity.
+ */
+coef_t Problem::query_variable_mod(int v, coef_t factor, int color, int nModularEQs, int nModularVars) const {
+ if (safeVars < v)
+ return posInfinity;
+
+ Tuple<bool> working_on(safeVars);
+ for (int i = 1; i <= safeVars; i++)
+ working_on[i] = false;
+
+ return query_variable_mod(v, factor, color, nModularEQs, nModularVars, working_on);
+}
+
+coef_t Problem::query_variable_mod(int v, coef_t factor, int color, int nModularEQs, int nModularVars, Tuple<bool> &working_on) const {
+ working_on[v] = true;
+
+ for (int e = 0; e < nEQs+nModularEQs; e++)
+ if (EQs[e].color == color) {
+ coef_t coef = int_mod_hat(EQs[e].coef[v], factor);
+ if (abs(coef) != 1)
+ continue;
+
+ bool wild_factored = true;
+ for (int i = safeVars+1; i <= ((e<nEQs)?nVars:nModularVars); i++)
+ if (int_mod_hat(EQs[e].coef[i], factor) != 0) {
+ wild_factored = false;
+ break;
+ }
+ if (!wild_factored)
+ continue;
+
+ coef_t result = 0;
+ for (int i = 1; i <= safeVars; i++)
+ if (i != v) {
+ coef_t p = int_mod_hat(EQs[e].coef[i], factor);
+
+ if (p == 0)
+ continue;
+
+ if (working_on[i] == true) {
+ result = posInfinity;
+ break;
+ }
+
+ coef_t q = query_variable_mod(i, factor, color, nModularEQs, nModularVars, working_on);
+ if (q == posInfinity) {
+ result = posInfinity;
+ break;
+ }
+ result += p*q;
+ }
+
+ if (result != posInfinity) {
+ result += EQs[e].coef[0];
+ if (coef == 1)
+ result = -result;
+ working_on[v] = false;
+
+ return int_mod_hat(result, factor);
+ }
+ }
+
+ working_on[v] = false;
+ return posInfinity;
+}
+
+
+
+#ifdef GIST_CHECK
+enum compareAnswer {apparentlyEqual, mightNotBeEqual, NotEqual};
+
+static compareAnswer checkEquiv(Problem *p1, Problem *p2) {
+ int r1,r2;
+
+ p1->check();
+ p2->check();
+ p1->resurrectSubs();
+ p2->resurrectSubs();
+ p1->check();
+ p2->check();
+ p1->putVariablesInStandardOrder();
+ p2->putVariablesInStandardOrder();
+ p1->check();
+ p2->check();
+ p1->ordered_elimination(0);
+ p2->ordered_elimination(0);
+ p1->check();
+ p2->check();
+ r1 = p1->simplifyProblem(1,1,0);
+ r2 = p2->simplifyProblem(1,1,0);
+ p1->check();
+ p2->check();
+
+ if (!r1 || !r2) {
+ if (r1 == r2) return apparentlyEqual;
+ return NotEqual;
+ }
+ if (p1->nVars != p2->nVars
+ || p1->nGEQs != p2->nGEQs
+ || p1->nSUBs != p2->nSUBs
+ || p1->checkSum() != p2->checkSum()) {
+ r1 = p1->simplifyProblem(0,1,1);
+ r2 = p2->simplifyProblem(0,1,1);
+ assert(r1 && r2);
+ p1->check();
+ p2->check();
+ if (p1->nVars != p2->nVars
+ || p1->nGEQs != p2->nGEQs
+ || p1->nSUBs != p2->nSUBs
+ || p1->checkSum() != p2->checkSum()) {
+ r1 = p1->simplifyProblem(0,1,2);
+ r2 = p2->simplifyProblem(0,1,2);
+ p1->check();
+ p2->check();
+ assert(r1 && r2);
+ if (p1->nVars != p2->nVars
+ || p1->nGEQs != p2->nGEQs
+ || p1->nSUBs != p2->nSUBs
+ || p1->checkSum() != p2->checkSum()) {
+ p1->check();
+ p2->check();
+ p1->resurrectSubs();
+ p2->resurrectSubs();
+ p1->check();
+ p2->check();
+ p1->putVariablesInStandardOrder();
+ p2->putVariablesInStandardOrder();
+ p1->check();
+ p2->check();
+ p1->ordered_elimination(0);
+ p2->ordered_elimination(0);
+ p1->check();
+ p2->check();
+ r1 = p1->simplifyProblem(1,1,0);
+ r2 = p2->simplifyProblem(1,1,0);
+ p1->check();
+ p2->check();
+ }
+ }
+ }
+
+ if (p1->nVars != p2->nVars
+ || p1->nSUBs != p2->nSUBs
+ || p1->nGEQs != p2->nGEQs
+ || p1->nSUBs != p2->nSUBs) return NotEqual;
+ if (p1->checkSum() != p2->checkSum()) return mightNotBeEqual;
+ return apparentlyEqual;
+}
+#endif
+
+void Problem::checkGistInvariant() const {
+#ifdef GIST_CHECK
+ Problem new_answer;
+ int r;
+
+ check();
+ fullAnswer.check();
+ context.check();
+
+ if (safeVars < nVars) {
+ if (DBUG) {
+ fprintf(outputFile,"Can't check gist invariant due to wildcards\n");
+ printProblem();
+ }
+ return;
+ }
+ if (DBUG) {
+ fprintf(outputFile,"Checking gist invariant on: [\n");
+ printProblem();
+ }
+
+ new_answer = *this;
+ new_answer->resurrectSubs();
+ new_answer->cleanoutWildcards();
+ if (DEBUG) {
+ fprintf(outputFile,"which is: \n");
+ printProblem();
+ }
+ deleteBlack(&new_answer);
+ turnRedBlack(&new_answer);
+ if (DEBUG) {
+ fprintf(outputFile,"Black version of answer: \n");
+ printProblem(&new_answer);
+ }
+ problem_merge(&new_answer,&context);
+
+ r = checkEquiv(&full_answer,&new_answer);
+ if (r != apparentlyEqual) {
+ fprintf(outputFile,"GIST INVARIANT REQUIRES MANUAL CHECK:[\n");
+ fprintf(outputFile,"Original problem:\n");
+ printProblem(&redProblem);
+
+ fprintf(outputFile,"Context:\n");
+ printProblem(&context);
+
+ fprintf(outputFile,"Computed gist:\n");
+ printProblem();
+
+ fprintf(outputFile,"Combined answer:\n");
+ printProblem(&full_answer);
+
+ fprintf(outputFile,"Context && red constraints:\n");
+ printProblem(&new_answer);
+ fprintf(outputFile,"]\n");
+ }
+
+ if (DBUG) {
+ fprintf(outputFile,"] Done checking gist invariant on\n");
+ }
+#endif
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_eq.cc b/omegalib/omega/src/omega_core/oc_eq.cc
new file mode 100644
index 0000000..dc595ea
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_eq.cc
@@ -0,0 +1,653 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+ Solve equalities.
+
+ Notes:
+
+ History:
+ *****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+
+namespace omega {
+
+void Problem::simplifyStrideConstraints() {
+ int e, e2, i;
+ if (DBUG)
+ fprintf(outputFile, "Checking for stride constraints\n");
+ for (i = safeVars + 1; i <= nVars; i++) {
+ if (DBUG)
+ fprintf(outputFile, "checking %s\n", variable(i));
+ for (e = 0; e < nGEQs; e++)
+ if (GEQs[e].coef[i])
+ break;
+ if (e >= nGEQs) {
+ if (DBUG)
+ fprintf(outputFile, "%s passed GEQ test\n", variable(i));
+ e2 = -1;
+ for (e = 0; e < nEQs; e++)
+ if (EQs[e].coef[i]) {
+ if (e2 == -1)
+ e2 = e;
+ else {
+ e2 = -1;
+ break;
+ }
+ }
+ if (e2 >= 0) {
+ if (DBUG) {
+ fprintf(outputFile, "Found stride constraint: ");
+ printEQ(&EQs[e2]);
+ fprintf(outputFile, "\n");
+ }
+ /* Is a stride constraint */
+ coef_t g = abs(EQs[e2].coef[i]);
+ assert(g>0);
+ int j;
+ for (j = 0; j <= nVars; j++)
+ if (i != j)
+ EQs[e2].coef[j] = int_mod_hat(EQs[e2].coef[j], g);
+ }
+ }
+ }
+}
+
+void Problem::doMod(coef_t factor, int e, int j) {
+ /* Solve e = factor alpha for x_j and substitute */
+ int k;
+ eqn eq;
+ coef_t nFactor;
+
+ int alpha;
+
+ // if (j > safeVars) alpha = j;
+ // else
+ if (EQs[e].color) {
+ rememberRedConstraint(&EQs[e], redEQ, 0);
+ EQs[e].color = EQ_BLACK;
+ }
+ alpha = addNewUnprotectedWildcard();
+ eqnncpy(&eq, &EQs[e], nVars);
+ newVar = alpha;
+
+ if (DEBUG) {
+ fprintf(outputFile, "doing moding: ");
+ fprintf(outputFile, "Solve ");
+ printTerm(&eq, 1);
+ fprintf(outputFile, " = " coef_fmt " %s for %s and substitute\n",
+ factor, variable(alpha), variable(j));
+ }
+ for (k = nVars; k >= 0; k--)
+ eq.coef[k] = int_mod_hat(eq.coef[k], factor);
+ nFactor = eq.coef[j];
+ assert(nFactor == 1 || nFactor == -1);
+ eq.coef[alpha] = factor / nFactor;
+ if (DEBUG) {
+ fprintf(outputFile, "adjusted: ");
+ fprintf(outputFile, "Solve ");
+ printTerm(&eq, 1);
+ fprintf(outputFile, " = 0 for %s and substitute\n", variable(j));
+ }
+
+ eq.coef[j] = 0;
+ substitute(&eq, j, nFactor);
+ newVar = -1;
+ deleteVariable(j);
+ for (k = nVars; k >= 0; k--) {
+ assert(EQs[e].coef[k] % factor == 0);
+ EQs[e].coef[k] = EQs[e].coef[k] / factor;
+ }
+ if (DEBUG) {
+ fprintf(outputFile, "Mod-ing and normalizing produces:\n");
+ printProblem();
+ }
+}
+
+void Problem::substitute(eqn *sub, int i, coef_t c) {
+ int e, j;
+ coef_t k;
+ int recordSubstitution = (i <= safeVars && var[i] >= 0);
+
+ redType clr;
+ if (sub->color)
+ clr = redEQ;
+ else
+ clr = notRed;
+
+ assert(c == 1 || c == -1);
+
+ if (DBUG || doTrace) {
+ if (sub->color)
+ fprintf(outputFile, "RED SUBSTITUTION\n");
+ fprintf(outputFile, "substituting using %s := ", variable(i));
+ printTerm(sub, -c);
+ fprintf(outputFile, "\n");
+ printVars();
+ }
+#ifndef NDEBUG
+ if (i > safeVars && clr) {
+ bool unsafeSub = false;
+ for (e = nEQs - 1; e >= 0; e--)
+ if (!(EQs[e].color || !EQs[e].coef[i]))
+ unsafeSub = true;
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (!(GEQs[e].color || !GEQs[e].coef[i]))
+ unsafeSub = true;
+ for (e = nSUBs - 1; e >= 0; e--)
+ if (SUBs[e].coef[i])
+ unsafeSub = true;
+ if (unsafeSub) {
+ fprintf(outputFile, "UNSAFE RED SUBSTITUTION\n");
+ fprintf(outputFile, "substituting using %s := ", variable(i));
+ printTerm(sub, -c);
+ fprintf(outputFile, "\n");
+ printProblem();
+ assert(0 && "UNSAFE RED SUBSTITUTION");
+ }
+ }
+#endif
+
+ for (e = nEQs - 1; e >= 0; e--) {
+ eqn *eq;
+ eq = &(EQs[e]);
+ k = eq->coef[i];
+ if (k != 0) {
+ k = check_mul(k, c); // Should be k = k/c, but same effect since abs(c) == 1
+ eq->coef[i] = 0;
+ for (j = nVars; j >= 0; j--) {
+ eq->coef[j] -= check_mul(sub->coef[j], k);
+ }
+ }
+ if (DEBUG) {
+ printEQ(eq);
+ fprintf(outputFile, "\n");
+ }
+ }
+ for (e = nGEQs - 1; e >= 0; e--) {
+ int zero;
+ eqn *eq;
+ eq = &(GEQs[e]);
+ k = eq->coef[i];
+ if (k != 0) {
+ k = check_mul(k, c); // Should be k = k/c, but same effect since abs(c) == 1
+ eq->touched = true;
+ eq->coef[i] = 0;
+ zero = 1;
+ for (j = nVars; j >= 0; j--) {
+ eq->coef[j] -= check_mul(sub->coef[j], k);
+ if (j > 0 && eq->coef[j])
+ zero = 0;
+ }
+ if (zero && clr != notRed && !eq->color) {
+ coef_t z = int_div(eq->coef[0], abs(k));
+ if (DBUG || doTrace) {
+ fprintf(outputFile,
+ "Black inequality matches red substitution\n");
+ if (z < 0)
+ fprintf(outputFile, "System is infeasible\n");
+ else if (z > 0)
+ fprintf(outputFile, "Black inequality is redundant\n");
+ else {
+ fprintf(outputFile,
+ "Black constraint partially implies red equality\n");
+ if (k < 0) {
+ fprintf(outputFile, "Black constraints tell us ");
+ assert(sub->coef[i] == 0);
+ sub->coef[i] = c;
+ printTerm(sub, 1);
+ sub->coef[i] = 0;
+ fprintf(outputFile, "<= 0\n");
+ } else {
+ fprintf(outputFile, "Black constraints tell us ");
+ assert(sub->coef[i] == 0);
+ sub->coef[i] = c;
+ printTerm(sub, 1);
+ sub->coef[i] = 0;
+ fprintf(outputFile, " >= 0\n");
+ }
+ }
+ }
+ if (z == 0) {
+ if (k < 0) {
+ if (clr == redEQ)
+ clr = redGEQ;
+ else if (clr == redLEQ)
+ clr = notRed;
+ } else {
+ if (clr == redEQ)
+ clr = redLEQ;
+ else if (clr == redGEQ)
+ clr = notRed;
+ }
+ }
+
+ }
+ }
+ if (DEBUG) {
+ printGEQ(eq);
+ fprintf(outputFile, "\n");
+ }
+ }
+ if (i <= safeVars && clr) {
+ assert(sub->coef[i] == 0);
+ sub->coef[i] = c;
+ rememberRedConstraint(sub, clr, 0);
+ sub->coef[i] = 0;
+ }
+
+ if (recordSubstitution) {
+ int s = nSUBs++;
+ int kk;
+ eqn *eq = &(SUBs[s]);
+ for (kk = nVars; kk >= 0; kk--)
+ eq->coef[kk] = check_mul(-c, (sub->coef[kk]));
+ eq->key = var[i];
+ if (DEBUG) {
+ fprintf(outputFile, "Recording substition as: ");
+ printSubstitution(s);
+ fprintf(outputFile, "\n");
+ }
+ }
+ if (DEBUG) {
+ fprintf(outputFile, "Ready to update subs\n");
+ if (sub->color)
+ fprintf(outputFile, "RED SUBSTITUTION\n");
+ fprintf(outputFile, "substituting using %s := ", variable(i));
+ printTerm(sub, -c);
+ fprintf(outputFile, "\n");
+ printVars();
+ }
+
+ for (e = nSUBs - 1; e >= 0; e--) {
+ eqn *eq = &(SUBs[e]);
+ k = eq->coef[i];
+ if (k != 0) {
+ k = check_mul(k, c); // Should be k = k/c, but same effect since abs(c) == 1
+ eq->coef[i] = 0;
+ for (j = nVars; j >= 0; j--) {
+ eq->coef[j] -= check_mul(sub->coef[j], k);
+ }
+ }
+ if (DEBUG) {
+ fprintf(outputFile, "updated sub (" coef_fmt "): ", c);
+ printSubstitution(e);
+ fprintf(outputFile, "\n");
+ }
+ }
+
+ if (DEBUG) {
+ fprintf(outputFile, "---\n\n");
+ printProblem();
+ fprintf(outputFile, "===\n\n");
+ }
+}
+
+
+void Problem::doElimination(int e, int i) {
+ if (DBUG || doTrace)
+ fprintf(outputFile, "eliminating variable %s\n", variable(i));
+
+ eqn sub;
+ eqnncpy(&sub, &EQs[e], nVars);
+ coef_t c = sub.coef[i];
+ sub.coef[i] = 0;
+
+ if (c == 1 || c == -1) {
+ substitute(&sub, i, c);
+ } else {
+ coef_t a = abs(c);
+ if (TRACE)
+ fprintf(outputFile,
+ "performing non-exact elimination, c = " coef_fmt "\n", c);
+ if (DBUG)
+ printProblem();
+ assert(inApproximateMode);
+
+ for (int e2 = nEQs - 1; e2 >= 0; e2--) {
+ eqn *eq = &(EQs[e2]);
+ coef_t k = eq->coef[i];
+ if (k != 0) {
+ coef_t l = lcm(abs(k), a);
+ coef_t scale1 = l / abs(k);
+ for (int j = nVars; j >= 0; j--)
+ eq->coef[j] = check_mul(eq->coef[j], scale1);
+ eq->coef[i] = 0;
+ coef_t scale2 = l / c;
+ if (k < 0)
+ scale2 = -scale2;
+ for (int j = nVars; j >= 0; j--)
+ eq->coef[j] -= check_mul(sub.coef[j], scale2);
+ eq->color |= sub.color;
+ }
+ }
+ for (int e2 = nGEQs - 1; e2 >= 0; e2--) {
+ eqn *eq = &(GEQs[e2]);
+ coef_t k = eq->coef[i];
+ if (k != 0) {
+ coef_t l = lcm(abs(k), a);
+ coef_t scale1 = l / abs(k);
+ for (int j = nVars; j >= 0; j--)
+ eq->coef[j] = check_mul(eq->coef[j], scale1);
+ eq->coef[i] = 0;
+ coef_t scale2 = l / c;
+ if (k < 0)
+ scale2 = -scale2;
+ for (int j = nVars; j >= 0; j--)
+ eq->coef[j] -= check_mul(sub.coef[j], scale2);
+ eq->color |= sub.color;
+ eq->touched = 1;
+ }
+ }
+ for (int e2 = nSUBs - 1; e2 >= 0; e2--)
+ if (SUBs[e2].coef[i]) {
+ eqn *eq = &(EQs[e2]);
+ assert(0);
+ // We can't handle this since we can't multiply
+ // the coefficient of the left-hand side
+ assert(!sub.color);
+ for (int j = nVars; j >= 0; j--)
+ eq->coef[j] = check_mul(eq->coef[j], a);
+ coef_t k = eq->coef[i];
+ eq->coef[i] = 0;
+ for (int j = nVars; j >= 0; j--)
+ eq->coef[j] -= check_mul(sub.coef[j], k / c);
+ }
+ }
+ deleteVariable(i);
+}
+
+int Problem::solveEQ() {
+ check();
+
+ // Reorder equations according to complexity.
+ {
+ int delay[nEQs];
+
+ for (int e = 0; e < nEQs; e++) {
+ delay[e] = 0;
+ if (EQs[e].color)
+ delay[e] += 8;
+ int nonunitWildCards = 0;
+ int unitWildCards = 0;
+ for (int i = nVars; i > safeVars; i--)
+ if (EQs[e].coef[i]) {
+ if (EQs[e].coef[i] == 1 || EQs[e].coef[i] == -1)
+ unitWildCards++;
+ else
+ nonunitWildCards++;
+ }
+ int unit = 0;
+ int nonUnit = 0;
+ for (int i = safeVars; i > 0; i--)
+ if (EQs[e].coef[i]) {
+ if (EQs[e].coef[i] == 1 || EQs[e].coef[i] == -1)
+ unit++;
+ else
+ nonUnit++;
+ }
+ if (unitWildCards == 1 && nonunitWildCards == 0)
+ delay[e] += 0;
+ else if (unitWildCards >= 1 && nonunitWildCards == 0)
+ delay[e] += 1;
+ else if (inApproximateMode && nonunitWildCards > 0)
+ delay[e] += 2;
+ else if (unit == 1 && nonUnit == 0 && nonunitWildCards == 0)
+ delay[e] += 3;
+ else if (unit > 1 && nonUnit == 0 && nonunitWildCards == 0)
+ delay[e] += 4;
+ else if (unit >= 1 && nonunitWildCards <= 1)
+ delay[e] += 5;
+ else
+ delay[e] += 6;
+ }
+
+ for (int e = 0; e < nEQs; e++) {
+ int e2, slowest;
+ slowest = e;
+ for (e2 = e + 1; e2 < nEQs; e2++)
+ if (delay[e2] > delay[slowest])
+ slowest = e2;
+ if (slowest != e) {
+ int tmp = delay[slowest];
+ delay[slowest] = delay[e];
+ delay[e] = tmp;
+ eqn eq;
+ eqnncpy(&eq, &EQs[slowest], nVars);
+ eqnncpy(&EQs[slowest], &EQs[e], nVars);
+ eqnncpy(&EQs[e], &eq, nVars);
+ }
+ }
+ }
+
+ // Eliminate all equations.
+ while (nEQs != 0) {
+ int e = nEQs - 1;
+ eqn *eq = &(EQs[e]);
+ coef_t g, g2;
+
+ assert(mayBeRed || !eq->color);
+
+ check();
+
+ // get gcd of coefficients of all unprotected variables
+ g2 = 0;
+ for (int i = nVars; i > safeVars; i--)
+ if (eq->coef[i] != 0) {
+ g2 = gcd(abs(eq->coef[i]), g2);
+ if (g2 == 1)
+ break;
+ }
+
+ // get gcd of coefficients of all variables
+ g = g2;
+ if (g != 1)
+ for (int i = safeVars; i >= 1; i--)
+ if (eq->coef[i] != 0) {
+ g = gcd(abs(eq->coef[i]), g);
+ if (g == 1)
+ break;
+ }
+
+ // approximate mode bypass integer modular test; in Farkas(),
+ // existential variable lambda's are rational numbers.
+ if (inApproximateMode && g2 != 0)
+ g = gcd(abs(eq->coef[0]), g);
+
+ // simple test to see if the equation is satisfiable
+ if (g == 0) {
+ if (eq->coef[0] != 0) {
+ return (false);
+ } else {
+ nEQs--;
+ continue;
+ }
+ } else if (abs(eq->coef[0]) % g != 0) {
+ return (false);
+ }
+
+ // set gcd of all coefficients to 1
+ if (g != 1) {
+ for (int i = nVars; i >= 0; i--)
+ eq->coef[i] /= g;
+ g2 = g2 / g;
+ }
+
+ // exact elimination of unit coefficient variable
+ if (g2 != 0) { // for constraint with unprotected variable
+ int i;
+ for (i = nVars; i > safeVars; i--)
+ if (abs(eq->coef[i]) == 1)
+ break;
+ if (i > safeVars) {
+ nEQs--;
+ doElimination(e, i);
+ continue;
+ }
+ } else { // for constraint without unprotected variable
+
+ // pick the unit coefficient variable with complex inequalites
+ // to eliminate, this will make inequalities tighter. e.g.
+ // {[t4,t6,t10]:exists (alpha: 0<=t6<=3 && t10=4alpha+t6 &&
+ // 64t4<=t10<=64t4+15)}
+ int unit_var;
+ int cost = -1;
+
+ for (int i = safeVars; i > 0; i--)
+
+ if (abs(eq->coef[i]) == 1) {
+ int cur_cost = 0;
+ for (int j = 0; j < nGEQs; j++)
+ if (GEQs[j].coef[i] != 0) {
+ for (int k = safeVars; k > 0; k--)
+ if (GEQs[j].coef[k] != 0) {
+ if (abs(GEQs[j].coef[k]) != 1){
+
+ cur_cost += 3;
+
+ }
+ else
+ cur_cost += 1;
+ }
+ }
+
+ if (cur_cost > cost) {
+ cost = cur_cost;
+ unit_var = i;
+ }
+
+ }
+
+ if (cost != -1) {
+ nEQs--;
+ doElimination(e, unit_var);
+ continue;
+ }
+
+
+ }
+
+ // check if there is an unprotected variable as wildcard
+ if (g2 > 0) {
+ int pos = 0;
+ coef_t g3;
+ for (int k = nVars; k > safeVars; k--)
+ if (eq->coef[k] != 0) {
+ int e2;
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (EQs[e2].coef[k])
+ break;
+ if (e2 >= 0)
+ continue;
+ for (e2 = nGEQs - 1; e2 >= 0; e2--)
+ if (GEQs[e2].coef[k])
+ break;
+ if (e2 >= 0)
+ continue;
+ for (e2 = nSUBs - 1; e2 >= 0; e2--)
+ if (SUBs[e2].coef[k])
+ break;
+ if (e2 >= 0)
+ continue;
+
+ if (pos == 0) {
+ g3 = abs(eq->coef[k]);
+ pos = k;
+ } else {
+ if (abs(eq->coef[k]) < g3) {
+ g3 = abs(eq->coef[k]);
+ pos = k;
+ }
+ }
+ }
+
+ if (pos != 0) {
+ bool change = false;
+ for (int k2 = nVars; k2 >= 0; k2--)
+ if (k2 != pos && eq->coef[k2] != 0) {
+ coef_t t = int_mod_hat(eq->coef[k2], g3);
+ if (t != eq->coef[k2]) {
+ eq->coef[k2] = t;
+ change = true;
+ }
+ }
+
+ // strength reduced, try this equation again
+ if (change) {
+ // nameWildcard(pos);
+ continue;
+ }
+ }
+ }
+
+ // insert new stride constraint
+ if (g2 > 1 && !(inApproximateMode && !inStridesAllowedMode)) {
+ int newvar = addNewProtectedWildcard();
+ int neweq = newEQ();
+ assert(neweq == e+1);
+ // we were working on highest-numbered EQ
+ eqnnzero(&EQs[neweq], nVars);
+ eqnncpy(&EQs[neweq], eq, safeVars);
+
+ for (int k = nVars; k >= 0; k--) {
+ EQs[neweq].coef[k] = int_mod_hat(EQs[neweq].coef[k], g2);
+ }
+ if (EQs[e].color)
+ rememberRedConstraint(&EQs[neweq], redStride, g2);
+ EQs[neweq].coef[newvar] = g2;
+ EQs[neweq].color = EQ_BLACK;
+ continue;
+ }
+
+ // inexact elimination of unprotected variable
+ if (g2 > 0 && inApproximateMode) {
+ int pos = 0;
+ for (int k = nVars; k > safeVars; k--)
+ if (eq->coef[k] != 0) {
+ pos = k;
+ break;
+ }
+ assert(pos > safeVars);
+
+ // special handling for wildcard used in breaking down
+ // diophantine equation
+ if (abs(eq->coef[pos]) > 1) {
+ int e2;
+ for (e2 = nSUBs - 1; e2 >= 0; e2--)
+ if (SUBs[e2].coef[pos])
+ break;
+ if (e2 >= 0) {
+ protectWildcard(pos);
+ continue;
+ }
+ }
+
+ nEQs--;
+ doElimination(e, pos);
+ continue;
+ }
+
+ // now solve linear diophantine equation using least remainder
+ // algorithm
+ {
+ coef_t factor = (posInfinity); // was MAXINT
+ int pos = 0;
+ for (int k = nVars; k > (g2 > 0 ? safeVars : 0); k--)
+ if (eq->coef[k] != 0 && factor > abs(eq->coef[k]) + 1) {
+ factor = abs(eq->coef[k]) + 1;
+ pos = k;
+ }
+ assert(pos > (g2>0?safeVars:0));
+ doMod(factor, e, pos);
+ continue;
+ }
+ }
+
+ assert(nEQs == 0);
+ return (OC_SOLVE_UNKNOWN);
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_exp_kill.cc b/omegalib/omega/src/omega_core/oc_exp_kill.cc
new file mode 100644
index 0000000..fdb2718
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_exp_kill.cc
@@ -0,0 +1,297 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+ Expensive inequality elimination.
+
+ Notes:
+
+ History:
+ 03/31/09 Use BoolSet, Chun Chen
+*****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+#include <basic/BoolSet.h>
+#include <vector>
+
+namespace omega {
+
+int Problem::expensiveKill() {
+ int e;
+ if (TRACE) fprintf(outputFile,"Performing expensive kill tests: [\n");
+ if (DBUG) printProblem();
+ Problem tmpProblem;
+ int oldTrace = trace;
+ int constraintsRemoved = 0;
+
+ trace = 0;
+ conservative++;
+
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (!GEQs[e].essential) {
+ if (DBUG) {
+ fprintf(outputFile, "checking equation %d to see if it is redundant: ", e);
+ printGEQ(&(GEQs[e]));
+ fprintf(outputFile, "\n");
+ }
+ tmpProblem = *this;
+ tmpProblem.negateGEQ(e);
+ tmpProblem.varsOfInterest = 0;
+ tmpProblem.nSUBs = 0;
+ tmpProblem.nMemories = 0;
+ tmpProblem.safeVars = 0;
+ tmpProblem.variablesFreed = 0;
+ tmpProblem.isTemporary = true;
+
+ if (!tmpProblem.solve(false)) {
+ if (DBUG)
+ fprintf(outputFile, "redundant!\n");
+ constraintsRemoved++;
+ deleteGEQ(e);
+ }
+ }
+
+ if (constraintsRemoved) {
+ if (TRACE) fprintf(outputFile,"%d Constraints removed!!\n",constraintsRemoved);
+ }
+
+ trace = oldTrace;
+ conservative--;
+ if (TRACE) fprintf(outputFile,"] expensive kill tests done\n");
+ return 1;
+}
+
+int Problem::expensiveRedKill() {
+ int e;
+ if (TRACE) fprintf(outputFile,"Performing expensive red kill tests: [\n");
+ Problem tmpProblem;
+ int oldTrace = trace;
+ int constraintsRemoved = 0;
+
+ trace = 0;
+ conservative++;
+
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (!GEQs[e].essential && GEQs[e].color) {
+ if (DEBUG) {
+ fprintf(outputFile, "checking equation %d to see if it is redundant: ", e);
+ printGEQ(&(GEQs[e]));
+ fprintf(outputFile, "\n");
+ }
+ tmpProblem = *this;
+ tmpProblem.negateGEQ(e);
+ tmpProblem.varsOfInterest = 0;
+ tmpProblem.nSUBs = 0;
+ tmpProblem.nMemories = 0;
+ tmpProblem.safeVars = 0;
+ tmpProblem.variablesFreed = 0;
+ tmpProblem.isTemporary = true;
+ tmpProblem.turnRedBlack();
+ if (!tmpProblem.solve(false)) {
+ constraintsRemoved++;
+ deleteGEQ(e);
+ }
+ }
+
+ if (constraintsRemoved) {
+ if (TRACE) fprintf(outputFile,"%d Constraints removed!!\n",constraintsRemoved);
+ }
+
+ trace = oldTrace;
+ conservative--;
+ if (TRACE) fprintf(outputFile,"] expensive red kill tests done\n");
+ return 1;
+}
+
+
+int Problem::expensiveEqualityCheck() {
+ int e;
+ return 1;
+ if (TRACE) fprintf(outputFile,"Performing expensive equality tests: [\n");
+ Problem tmpProblem;
+ int oldTrace = trace;
+ int equalitiesFound = 0;
+
+ trace = 0;
+ conservative++;
+
+ for (e = nGEQs - 1; e >= 0; e--) {
+ if (DEBUG) {
+ fprintf(outputFile, "checking equation %d to see if it is an equality: ", e);
+ printGEQ(&(GEQs[e]));
+ fprintf(outputFile, "\n");
+ }
+ tmpProblem = *this;
+ tmpProblem.GEQs[e].coef[0]--;
+ tmpProblem.varsOfInterest = 0;
+ tmpProblem.nSUBs = 0;
+ tmpProblem.nMemories = 0;
+ tmpProblem.safeVars = 0;
+ tmpProblem.variablesFreed = 0;
+ tmpProblem.isTemporary = true;
+ if (!tmpProblem.solve(false)) {
+ int neweq = newEQ();
+ eqnncpy(&EQs[neweq], &GEQs[e], nVars);
+ equalitiesFound++;
+ addingEqualityConstraint(neweq);
+ }
+ }
+ if (equalitiesFound) {
+ if (TRACE) fprintf(outputFile,"%d Equalities found!!\n",equalitiesFound);
+ }
+
+ trace = oldTrace;
+ conservative--;
+ if (equalitiesFound) {
+ if (!solveEQ()) return 0;
+ if (!normalize()) return 0;
+ }
+ if (TRACE) fprintf(outputFile,"] expensive equality tests done\n");
+ return 1;
+}
+
+
+void Problem::quickRedKill(int computeGist) {
+ if (DBUG) {
+ fprintf(outputFile, "in quickRedKill: [\n");
+ printProblem();
+ }
+
+ noteEssential(0);
+ int moreToDo = chainKill(1,0);
+
+#ifdef NDEBUG
+ if (!moreToDo) {
+ if (DBUG) fprintf(outputFile, "] quickRedKill\n");
+ return;
+ }
+#endif
+
+ int isDead[nGEQs];
+ int deadCount = 0;
+ std::vector<BoolSet<> > P(nGEQs, BoolSet<>(nVars)), Z(nGEQs, BoolSet<>(nVars)), N(nGEQs, BoolSet<>(nVars));
+ BoolSet<> PP, PZ, PN; /* possible Positives, possible zeros & possible negatives */
+ BoolSet<> MZ; /* must zeros */
+
+ int equationsToKill = 0;
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ isDead[e] = 0;
+ if (GEQs[e].color && !GEQs[e].essential) equationsToKill++;
+ if (GEQs[e].color && GEQs[e].essential && !computeGist)
+ if (!moreToDo) {
+ if (DBUG) fprintf(outputFile, "] quickRedKill\n");
+ return;
+ }
+ for (int i = nVars; i >= 1; i--) {
+ if (GEQs[e].coef[i] > 0)
+ P[e].set(i-1);
+ else if (GEQs[e].coef[i] < 0)
+ N[e].set(i-1);
+ else
+ Z[e].set(i-1);
+ }
+ }
+
+ if (!equationsToKill)
+ if (!moreToDo) {
+ if (DBUG) fprintf(outputFile, "] quickRedKill\n");
+ return;
+ }
+
+ for (int e = nGEQs - 1; e > 0; e--)
+ if (!isDead[e])
+ for (int e2 = e - 1; e2 >= 0; e2--)
+ if (!isDead[e2]) {
+ coef_t a = 0;
+ int i, j;
+ for (i = nVars; i > 1; i--) {
+ for (j = i - 1; j > 0; j--) {
+ a = (GEQs[e].coef[i] * GEQs[e2].coef[j] - GEQs[e2].coef[i] * GEQs[e].coef[j]);
+ if (a != 0)
+ goto foundPair;
+ }
+ }
+ continue;
+
+ foundPair:
+ if (DEBUG) {
+ fprintf(outputFile, "found two equations to combine, i = %s, ", variable(i));
+ fprintf(outputFile, "j = %s, alpha = " coef_fmt "\n", variable(j), a);
+ printGEQ(&(GEQs[e]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n");
+ }
+
+ MZ = (Z[e] & Z[e2]);
+ PZ = MZ | (P[e] & N[e2]) | (N[e] & P[e2]);
+ PP = P[e] | P[e2];
+ PN = N[e] | N[e2];
+
+ for (int e3 = nGEQs - 1; e3 >= 0; e3--)
+ if (e3 != e && e3 != e2 && GEQs[e3].color && !GEQs[e3].essential) {
+ coef_t alpha1, alpha2, alpha3;
+
+ if (!PZ.imply(Z[e3]) || MZ.imply(~Z[e3])) continue;
+ if (!PP.imply(P[e3]) || !PN.imply(N[e3])) continue;
+
+ if (a > 0) {
+ alpha1 = GEQs[e2].coef[j] * GEQs[e3].coef[i] - GEQs[e2].coef[i] * GEQs[e3].coef[j];
+ alpha2 = -(GEQs[e].coef[j] * GEQs[e3].coef[i] - GEQs[e].coef[i] * GEQs[e3].coef[j]);
+ alpha3 = a;
+ }
+ else {
+ alpha1 = -(GEQs[e2].coef[j] * GEQs[e3].coef[i] - GEQs[e2].coef[i] * GEQs[e3].coef[j]);
+ alpha2 = -(-(GEQs[e].coef[j] * GEQs[e3].coef[i] - GEQs[e].coef[i] * GEQs[e3].coef[j]));
+ alpha3 = -a;
+ }
+
+ if (alpha1 > 0 && alpha2 > 0) {
+ if (DEBUG) {
+ fprintf(outputFile, "alpha1 = " coef_fmt ", alpha2 = " coef_fmt "; comparing against: ", alpha1, alpha2);
+ printGEQ(&(GEQs[e3]));
+ fprintf(outputFile, "\n");
+ }
+ coef_t c;
+ int k;
+ for (k = nVars; k >= 0; k--) {
+ c = alpha1 * GEQs[e].coef[k] + alpha2 * GEQs[e2].coef[k];
+ if (DEBUG) {
+ if (k>0)
+ fprintf(outputFile, " %s: " coef_fmt ", " coef_fmt "\n", variable(k), c, alpha3 * GEQs[e3].coef[k]);
+ else fprintf(outputFile, " constant: " coef_fmt ", " coef_fmt "\n", c, alpha3 * GEQs[e3].coef[k]);
+ }
+ if (c != alpha3 * GEQs[e3].coef[k])
+ break;
+ }
+ if (k < 0 || (k == 0 && c < alpha3 * (GEQs[e3].coef[k]+1))) {
+ if (DEBUG) {
+ deadCount++;
+ fprintf(outputFile, "red equation#%d is dead (%d dead so far, %d remain)\n", e3, deadCount, nGEQs - deadCount);
+ printGEQ(&(GEQs[e]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e3]));
+ fprintf(outputFile, "\n");
+ assert(moreToDo);
+ }
+ isDead[e3] = 1;
+ }
+ }
+ }
+ }
+
+ for (int e = nGEQs - 1; e >= 0; e--)
+ if (isDead[e])
+ deleteGEQ(e);
+
+ if (DBUG) {
+ fprintf(outputFile,"] quickRedKill\n");
+ printProblem();
+ }
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_global.cc b/omegalib/omega/src/omega_core/oc_global.cc
new file mode 100644
index 0000000..17d8a0c
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_global.cc
@@ -0,0 +1,45 @@
+#include <omega/omega_core/oc_i.h>
+
+namespace omega {
+
+const int Problem::min_alloc = 10;
+const int Problem::first_alloc_pad = 5;
+
+int omega_core_debug = 0; // 3: full debugging info
+
+int maxEQs = 100; // original 35, increased by chun
+int maxGEQs = 200; // original 70, increased by chun
+
+int newVar = -1;
+int findingImplicitEqualities = 0;
+int firstCheckForRedundantEquations = 0;
+int doItAgain;
+int conservative = 0;
+FILE *outputFile = stderr; /* printProblem writes its output to this file */
+char wildName[200][20];
+int nextWildcard = 0;
+int trace = 1;
+int depth = 0;
+int headerLevel;
+int inApproximateMode = 0;
+int inStridesAllowedMode = 0;
+int addingOuterEqualities = 0;
+int outerColor = 0;
+int reduceWithSubs = 1;
+int pleaseNoEqualitiesInSimplifiedProblems = 0;
+Problem *originalProblem = noProblem;
+int omegaInitialized = 0;
+int mayBeRed = 0;
+
+
+// Hash table is used to hash all inequalties for all problems. It
+// persists across problems for quick problem merging in case. When
+// the table is filled to 1/3 full, it is flushed and the filling
+// process starts all over again.
+int packing[maxVars];
+int hashVersion = 0;
+eqn hashMaster[hashTableSize];
+int fastLookup[maxKeys*2];
+int nextKey;
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_print.cc b/omegalib/omega/src/omega_core/oc_print.cc
new file mode 100644
index 0000000..7934713
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_print.cc
@@ -0,0 +1,686 @@
+#include <omega/omega_core/oc_i.h>
+
+namespace omega {
+
+int print_in_code_gen_style = 0;
+
+void Problem::initializeVariables() const {
+ Problem *p = (Problem *)this;
+ int i;
+ assert(!p->variablesInitialized);
+ for (i = p->nVars; i >= 0; i--)
+ p->forwardingAddress[i] = p->var[i] = i;
+ p->variablesInitialized = 1;
+}
+
+
+std::string Problem::print_term_to_string(const eqn *e, int c) const {
+ std::string s="";
+ int i;
+ int first;
+ int n = nVars;
+ int wentFirst = -1;
+ first = 1;
+ for (i = 1; i <= n; i++)
+ if (c * e->coef[i] > 0) {
+ first = 0;
+ wentFirst = i;
+
+ if (c * e->coef[i] == 1)
+ s+= variable(i);
+ else {
+ s += to_string(c * e->coef[i]);
+ if (print_in_code_gen_style) s += "*";
+ s += variable(i);
+ }
+ break;
+ }
+ for (i = 1; i <= n; i++)
+ if (i != wentFirst && c * e->coef[i] != 0) {
+ if (!first && c * e->coef[i] > 0)
+ s += "+";
+
+ first = 0;
+
+ if (c * e->coef[i] == 1)
+ s += variable(i);
+ else if (c * e->coef[i] == -1) {
+ s += "-"; s += variable(i);
+ }
+ else {
+ s += to_string(c * e->coef[i]);
+ if (print_in_code_gen_style) s += "*";
+ s += variable(i);
+ }
+ }
+ if (!first && c * e->coef[0] > 0)
+ s += "+";
+ if (first || c * e->coef[0] != 0)
+ s += to_string(c * e->coef[0]);
+ return s;
+}
+
+
+void Problem::printTerm(const eqn * e, int c) const {
+ std::string s = print_term_to_string(e, c);
+ fprintf(outputFile, "%s", s.c_str());
+}
+
+
+void Problem::printSub(int v) const {
+ std::string s = print_sub_to_string(v);
+ fprintf(outputFile, "%s", s.c_str());
+}
+
+
+std::string Problem::print_sub_to_string(int v) const {
+ std::string s;
+
+ if (v > 0)
+ s = variable(v);
+ else
+ s = print_term_to_string(&SUBs[-v-1], 1);
+ return s;
+}
+
+
+void Problem::clearSubs() {
+ nSUBs = 0;
+ nMemories = 0;
+}
+
+
+void Problem::printEqn(const eqn *e, int test, int extra) const {
+ char buf[maxVars * 12 + 180]; // original buf[maxVars * 12 + 80]
+
+ sprintEqn(buf, e, test, extra);
+ fprintf(outputFile, "%s", buf);
+}
+
+
+std::string Problem::printEqnToString(const eqn *e, int test, int extra) const {
+ char buf[maxVars * 12 + 180]; // original buf[maxVars * 12 + 80]
+ sprintEqn(buf, e, test, extra);
+ return std::string(buf);
+}
+
+
+void Problem::sprintEqn(char *str, const eqn *e, int test, int extra) const {
+ int i;
+ int first;
+ int n = nVars + extra;
+ int isLT;
+
+ isLT = test && e->coef[0] == -1;
+ if (isLT)
+ isLT = 1;
+#if 0
+ if (test) {
+ if (DEBUG && e->touched) {
+ sprintf(str, "!");
+ while (*str)
+ str++;
+ }
+ else if (DBUG && !e->touched && e->key != 0) {
+ sprintf(str, "%d: ", e->key);
+ while (*str)
+ str++;
+ }
+ }
+#endif
+ if (e->color) {
+ sprintf(str, "[");
+ while (*str)
+ str++;
+ }
+ first = 1;
+ for (i = isLT; i <= n; i++)
+ if (e->coef[i] < 0) {
+ if (!first) {
+ sprintf(str, "+");
+ while (*str)
+ str++;
+ }
+ else
+ first = 0;
+ if (i == 0) {
+ sprintf(str, coef_fmt, -e->coef[i]);
+ while (*str)
+ str++;
+ }
+ else if (e->coef[i] == -1) {
+ sprintf(str, "%s", variable(i));
+ while (*str)
+ str++;
+ }
+ else {
+ if (print_in_code_gen_style)
+ sprintf(str, coef_fmt "*%s", -e->coef[i], variable(i));
+ else sprintf(str, coef_fmt "%s", -e->coef[i], variable(i));
+ while (*str)
+ str++;
+ }
+ }
+ if (first) {
+ if (isLT) {
+ sprintf(str, "1");
+ isLT = 0;
+ }
+ else
+ sprintf(str, "0");
+ while (*str)
+ str++;
+ }
+ if (test == 0) {
+ if (print_in_code_gen_style) sprintf(str, " == ");
+ else sprintf(str, " = ");
+ while (*str)
+ str++;
+ }
+ else {
+ if (isLT)
+ sprintf(str, " < ");
+ else
+ sprintf(str, " <= ");
+ while (*str)
+ str++;
+ }
+
+ first = 1;
+ for (i = 0; i <= n; i++)
+ if (e->coef[i] > 0) {
+ if (!first) {
+ sprintf(str, "+");
+ while (*str)
+ str++;
+ }
+ else
+ first = 0;
+ if (i == 0) {
+ sprintf(str, coef_fmt , e->coef[i]);
+ while (*str)
+ str++;
+ }
+ else if (e->coef[i] == 1) {
+ sprintf(str, "%s", variable(i));
+ while (*str)
+ str++;
+ }
+ else {
+ if (print_in_code_gen_style)
+ sprintf(str, coef_fmt "*%s", e->coef[i], variable(i));
+ else
+ sprintf(str, coef_fmt "%s", e->coef[i], variable(i));
+ while (*str)
+ str++;
+ }
+ }
+ if (first) {
+ sprintf(str, "0");
+ while (*str)
+ str++;
+ }
+ if (e->color) {
+ sprintf(str, "]");
+ while (*str)
+ str++;
+ }
+}
+
+
+void Problem::printSubstitution(int s) const {
+ const eqn * eq = &(SUBs[s]);
+ assert(eq->key > 0);
+ fprintf(outputFile, "%s := ", orgVariable(eq->key));
+ printTerm(eq, 1);
+}
+
+
+void Problem::printVars(int /*debug*/) const {
+ int i;
+ fprintf(outputFile, "variables = ");
+ if (safeVars > 0)
+ fprintf(outputFile, "(");
+ for (i = 1; i <= nVars; i++) {
+ fprintf(outputFile, "%s", variable(i));
+ if (i == safeVars)
+ fprintf(outputFile, ")");
+ if (i < nVars)
+ fprintf(outputFile, ", ");
+ }
+ fprintf(outputFile, "\n");
+ /*
+ fprintf(outputFile, "forward addresses = ");
+ if (safeVars > 0)
+ fprintf(outputFile, "(");
+ for (i = 1; i <= nVars; i++)
+ {
+ int v = forwardingAddress[i];
+ if (v > 0) fprintf(outputFile, "%s", variable(i));
+ else fprintf(outputFile, "*");
+ if (i == safeVars)
+ fprintf(outputFile, ")");
+ if (i < nVars)
+ fprintf(outputFile, ", ");
+ };
+ fprintf(outputFile, "\n");
+ */
+}
+
+
+void printHeader() {
+ int i;
+ for(i=0; i<headerLevel; i++) {
+ fprintf(outputFile, ". ");
+ }
+}
+
+
+void Problem::printProblem(int debug) const {
+ int e;
+
+ if (!variablesInitialized)
+ initializeVariables();
+ if (debug) {
+ printHeader();
+ fprintf(outputFile, "#vars = %d, #EQ's = %d, #GEQ's = %d, # SUB's = %d, ofInterest = %d\n",
+ nVars,nEQs,nGEQs,nSUBs,varsOfInterest);
+ printHeader();
+ printVars(debug);
+ }
+ for (e = 0; e < nEQs; e++) {
+ printHeader();
+ printEQ(&EQs[e]);
+ fprintf(outputFile, "\n");
+ }
+ for (e = 0; e < nGEQs; e++) {
+ printHeader();
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile, "\n");
+ }
+ for (e = 0; e < nSUBs; e++) {
+ printHeader();
+ printSubstitution(e);
+ fprintf(outputFile, "\n");
+ }
+
+ for (e = 0; e < nMemories; e++) {
+ int i;
+ printHeader();
+ switch(redMemory[e].kind) {
+ case notRed:
+ fprintf(outputFile,"notRed: ");
+ break;
+ case redGEQ:
+ fprintf(outputFile,"Red: 0 <= ");
+ break;
+ case redLEQ:
+ fprintf(outputFile,"Red ??: 0 >= ");
+ break;
+ case redEQ:
+ fprintf(outputFile,"Red: 0 == ");
+ break;
+ case redStride:
+ fprintf(outputFile,"Red stride " coef_fmt ": ", redMemory[e].stride);
+ break;
+ }
+ fprintf(outputFile," " coef_fmt, redMemory[e].constantTerm);
+ for(i=0;i< redMemory[e].length; i++)
+ if(redMemory[e].coef[i] >= 0)
+ fprintf(outputFile,"+" coef_fmt "%s", redMemory[e].coef[i], orgVariable(redMemory[e].var[i]));
+ else
+ fprintf(outputFile,"-" coef_fmt "%s", -redMemory[e].coef[i], orgVariable(redMemory[e].var[i]));
+ fprintf(outputFile, "\n");
+ }
+ fflush(outputFile);
+
+ CHECK_FOR_DUPLICATE_VARIABLE_NAMES;
+}
+
+
+void Problem::printRedEquations() const {
+ int e;
+
+ if (!variablesInitialized)
+ initializeVariables();
+ printVars(1);
+ for (e = 0; e < nEQs; e++) {
+ if (EQs[e].color == EQ_RED) {
+ printEQ(&EQs[e]);
+ fprintf(outputFile, "\n");
+ }
+ }
+ for (e = 0; e < nGEQs; e++) {
+ if (GEQs[e].color == EQ_RED) {
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile, "\n");
+ }
+ }
+ for (e = 0; e < nSUBs; e++) {
+ if (SUBs[e].color) {
+ printSubstitution(e);
+ fprintf(outputFile, "\n");
+ }
+ }
+ fflush(outputFile);
+}
+
+
+int Problem::prettyPrintProblem() const {
+ std::string s = prettyPrintProblemToString();
+ fprintf(outputFile, "%s", s.c_str());
+ fflush(outputFile);
+ return 0;
+}
+
+
+std::string Problem::prettyPrintProblemToString() const {
+ std::string s="";
+ int e;
+ int v;
+ int live[maxmaxGEQs];
+ int v1, v2, v3;
+ int t, change;
+ int stuffPrinted = 0;
+ const char *connector = " && ";
+
+ typedef enum {
+ none, le, lt
+ } partialOrderType;
+
+ partialOrderType po[maxVars][maxVars];
+ int poE[maxVars][maxVars];
+ int lastLinks[maxVars];
+ int firstLinks[maxVars];
+ int chainLength[maxVars];
+ int chain[maxVars];
+ int varCount[maxVars];
+ int i, m, multiprint;
+
+
+ if (!variablesInitialized)
+ initializeVariables();
+
+ if (nVars > 0) {
+ for (e = 0; e < nEQs; e++) {
+ if (stuffPrinted)
+ s += connector;
+ stuffPrinted = 1;
+ s += print_EQ_to_string(&EQs[e]);
+ }
+
+ for (e = 0; e < nGEQs; e++) {
+ live[e] = true;
+ varCount[e] = 0;
+ for (v = 1; v <= nVars; v++)
+ if (GEQs[e].coef[v]) varCount[e]++;
+ }
+
+ if (!print_in_code_gen_style)
+ while (1) {
+ for (v = 1; v <= nVars; v++) {
+ lastLinks[v] = firstLinks[v] = 0;
+ chainLength[v] = 0;
+ for (v2 = 1; v2 <= nVars; v2++)
+ po[v][v2] = none;
+ }
+
+ for (e = 0; e < nGEQs; e++)
+ if (live[e] && varCount[e] <= 2) {
+ for (v = 1; v <= nVars; v++) {
+ if (GEQs[e].coef[v] == 1)
+ firstLinks[v]++;
+ else if (GEQs[e].coef[v] == -1)
+ lastLinks[v]++;
+ }
+
+ v1 = nVars;
+ while (v1 > 0 && GEQs[e].coef[v1] == 0)
+ v1--;
+ v2 = v1 - 1;
+ while (v2 > 0 && GEQs[e].coef[v2] == 0)
+ v2--;
+ v3 = v2 - 1;
+ while (v3 > 0 && GEQs[e].coef[v3] == 0)
+ v3--;
+
+ if (GEQs[e].coef[0] > 0 || GEQs[e].coef[0] < -1
+ || v2 <= 0 || v3 > 0
+ || GEQs[e].coef[v1] * GEQs[e].coef[v2] != -1) {
+ /* Not a partial order relation */
+
+ }
+ else {
+ if (GEQs[e].coef[v1] == 1) {
+ v3 = v2;
+ v2 = v1;
+ v1 = v3;
+ }
+ /* relation is v1 <= v2 or v1 < v2 */
+ po[v1][v2] = ((GEQs[e].coef[0] == 0) ? le : lt);
+ poE[v1][v2] = e;
+ }
+ }
+ for (v = 1; v <= nVars; v++)
+ chainLength[v] = lastLinks[v];
+
+ /*
+ * printf("\n\nPartial order:\n"); printf(" "); for (v1 = 1; v1 <= nVars; v1++)
+ * printf("%7s",variable(v1)); printf("\n"); for (v1 = 1; v1 <= nVars; v1++) { printf("%6s:
+ * ",variable(v1)); for (v2 = 1; v2 <= nVars; v2++) switch (po[v1][v2]) { case none: printf(" ");
+ * break; case le: printf(" <= "); break; case lt: printf(" < "); break; } printf("\n"); }
+ */
+
+
+ /* Just in case nVars <= 0 */
+ change = false;
+ for (t = 0; t < nVars; t++) {
+ change = false;
+ for (v1 = 1; v1 <= nVars; v1++)
+ for (v2 = 1; v2 <= nVars; v2++)
+ if (po[v1][v2] != none &&
+ chainLength[v1] <= chainLength[v2]) {
+ chainLength[v1] = chainLength[v2] + 1;
+ change = true;
+ }
+ }
+
+ if (change) {
+ /* caught in cycle */
+
+#if 0
+ printf("\n\nPartial order:\n"); printf(" ");
+ for (v1 = 1; v1 <= nVars; v1++) printf("%7s",variable(v1)); printf("\n");
+ for (v1 = 1; v1 <= nVars; v1++) {
+ printf("%6s: ",variable(v1));
+ for (v2 = 1; v2 <= nVars; v2++) switch (po[v1][v2]) {
+ case none: printf(" "); break;
+ case le: printf(" <= "); break;
+ case lt: printf(" < "); break;
+ }
+ printf("\n");
+ }
+
+ printProblem(1);
+#endif
+ break;
+ }
+
+ for (v1 = 1; v1 <= nVars; v1++)
+ if (chainLength[v1] == 0)
+ firstLinks[v1] = 0;
+
+ v = 1;
+ for (v1 = 2; v1 <= nVars; v1++)
+ if (chainLength[v1] + firstLinks[v1] > chainLength[v] + firstLinks[v])
+ v = v1;
+
+ if (chainLength[v] + firstLinks[v] == 0)
+ break;
+
+ if (stuffPrinted)
+ s += connector;
+ stuffPrinted = 1;
+ /* chain starts at v */
+ /* print firstLinks */
+ {
+ coef_t tmp;
+ int first;
+ first = 1;
+ for (e = 0; e < nGEQs; e++)
+ if (live[e] && GEQs[e].coef[v] == 1 && varCount[e] <= 2) {
+ if (!first)
+ s += ", ";
+ tmp = GEQs[e].coef[v];
+ ((Problem *)this)->
+ GEQs[e].coef[v] = 0;
+ s += print_term_to_string(&GEQs[e], -1);
+ ((Problem *)this)->
+ GEQs[e].coef[v] = tmp;
+ live[e] = false;
+ first = 0;
+ }
+ if (!first)
+ s += " <= ";
+ }
+
+
+ /* find chain */
+ chain[0] = v;
+ m = 1;
+ while (1) {
+ /* print chain */
+ for (v2 = 1; v2 <= nVars; v2++)
+ if (po[v][v2] && chainLength[v] == 1 + chainLength[v2])
+ break;
+ if (v2 > nVars)
+ break;
+ chain[m++] = v2;
+ v = v2;
+ }
+
+ s += variable(chain[0]);
+
+ multiprint = 0;
+ for (i = 1; i < m; i++) {
+ v = chain[i - 1];
+ v2 = chain[i];
+ if (po[v][v2] == le)
+ s += " <= ";
+ else
+ s += " < ";
+ s += variable(v2);
+ live[poE[v][v2]] = false;
+ if (!multiprint && i < m - 1)
+ for (v3 = 1; v3 <= nVars; v3++) {
+ if (v == v3 || v2 == v3)
+ continue;
+ if (po[v][v2] != po[v][v3])
+ continue;
+ if (po[v2][chain[i + 1]] != po[v3][chain[i + 1]])
+ continue;
+ s += ","; s += variable(v3);
+ live[poE[v][v3]] = false;
+ live[poE[v3][chain[i + 1]]] = false;
+ multiprint = 1;
+ }
+ else
+ multiprint = 0;
+ }
+
+ v = chain[m - 1];
+ /* print lastLinks */
+ {
+ coef_t tmp;
+ int first;
+ first = 1;
+ for (e = 0; e < nGEQs; e++)
+ if (live[e] && GEQs[e].coef[v] == -1 && varCount[e] <= 2) {
+ if (!first)
+ s += ", ";
+ else
+ s += " <= ";
+ tmp = GEQs[e].coef[v];
+ ((Problem *)this)->
+ GEQs[e].coef[v] = 0;
+ s += print_term_to_string(&GEQs[e], 1);
+ ((Problem *)this)->
+ GEQs[e].coef[v] = tmp;
+ live[e] = false;
+ first = 0;
+ }
+ }
+ }
+
+
+ for (e = 0; e < nGEQs; e++)
+ if (live[e]) {
+ if (stuffPrinted)
+ s += connector;
+ stuffPrinted = 1;
+ s += print_GEQ_to_string(&GEQs[e]);
+ }
+
+ for (e = 0; e < nSUBs; e++) {
+ const eqn * eq = &SUBs[e];
+ if (stuffPrinted)
+ s += connector;
+ stuffPrinted = 1;
+ if (eq->key > 0) {
+ s += orgVariable(eq->key); s += " := ";
+ }
+ else {
+ s += "#"; s += to_string(eq->key); s += " := ";
+ }
+ s += print_term_to_string(eq, 1);
+ }
+ }
+ return s;
+}
+
+
+int Problem::prettyPrintRedEquations() const {
+ int e, stuffPrinted = 0;
+ const char *connector = " && ";
+
+ if (!variablesInitialized)
+ initializeVariables();
+
+ for (e = 0; e < nEQs; e++) {
+ if (EQs[e].color == EQ_RED) {
+ if (stuffPrinted)
+ fprintf(outputFile, "%s", connector);
+ stuffPrinted = 1;
+ ((Problem *)this)->
+ EQs[e].color = EQ_BLACK;
+ printEQ(&EQs[e]);
+ ((Problem *)this)->
+ EQs[e].color = EQ_RED;
+ }
+ }
+ for (e = 0; e < nGEQs; e++) {
+ if (GEQs[e].color == EQ_RED) {
+ if (stuffPrinted)
+ fprintf(outputFile, "%s", connector);
+ stuffPrinted = 1;
+ ((Problem *)this)->
+ GEQs[e].color = EQ_BLACK;
+ printGEQ(&GEQs[e]);
+ ((Problem *)this)->
+ GEQs[e].color = EQ_RED;
+ }
+ }
+ for (e = 0; e < nSUBs; e++) {
+ if (SUBs[e].color) {
+ if (stuffPrinted)
+ fprintf(outputFile, "%s", connector);
+ stuffPrinted = 1;
+ printSubstitution(e);
+ }
+ }
+ fflush(outputFile);
+
+ return 0;
+}
+
+}
diff --git a/omegalib/omega/src/omega_core/oc_problems.cc b/omegalib/omega/src/omega_core/oc_problems.cc
new file mode 100644
index 0000000..8b6e04c
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_problems.cc
@@ -0,0 +1,198 @@
+#include <omega/omega_core/oc_i.h>
+#include <basic/omega_error.h>
+
+namespace omega {
+
+Problem::~Problem() {
+ delete[] EQs;
+ delete[] GEQs;
+}
+
+
+void check_number_EQs(int n) {
+ if (n < 0) {
+ fprintf(stderr,"ERROR: nEQs < 0??\n");
+ exit(1);
+ }
+
+ if (n > maxmaxEQs) {
+ // clear global variables
+ inApproximateMode = 0;
+ outerColor = 0;
+
+ throw presburger_error("\nERROR:\n"
+ "An attempt was made to set the number of available equality constraints to " + to_string(n) + ".\n"
+ "The maximum number of equality constraints in a conjunction is " + to_string(maxmaxEQs) + ".\n"
+ "This limit can be changed by redefining maxmaxEQs in oc.h and recompiling.\n\n");
+
+ // fprintf(stderr, "\nERROR:\n");
+ // fprintf(stderr, "An attempt was made to set the number of available equality constraints to %d.\n", n);
+ // fprintf(stderr, "The maximum number of equality constraints in a conjunction is %d.\n", maxmaxEQs);
+ // fprintf(stderr, "This limit can be changed by redefining maxmaxEQs in oc.h and recompiling.\n\n");
+ // exit(2);
+ }
+}
+
+void check_number_GEQs(int n) {
+ if (n < 0) {
+ fprintf(stderr,"ERROR: nGEQs < 0??\n");
+ exit(1);
+ }
+
+ if (n > maxmaxGEQs) {
+ // clear global variables
+ inApproximateMode = 0;
+ outerColor = 0;
+
+ throw presburger_error("\nERROR:\n"
+ "An attempt was made to set the number of available inequality constraints to " + to_string(n) +".\n"
+ "The maximum number of inequality constraints in a conjunction is " + to_string(maxmaxGEQs) +".\n"
+ "This limit can be changed by redefining maxmaxGEQs in oc.h and recompiling.\n\n");
+
+ // fprintf(stderr, "\nERROR:\n");
+ // fprintf(stderr, "An attempt was made to set the number of available inequality constraints to %d.\n", n);
+ // fprintf(stderr, "The maximum number of inequality constraints in a conjunction is %d.\n", maxmaxGEQs);
+ // fprintf(stderr, "This limit can be changed by redefining maxmaxGEQs in oc.h and recompiling.\n\n");
+ // exit(2);
+ }
+}
+
+
+void check_number_EQs_GEQs(int e, int g) {
+ check_number_EQs(e);
+ check_number_GEQs(g);
+}
+
+
+Problem::Problem(int in_eqs, int in_geqs) {
+ check_number_EQs_GEQs(in_eqs, in_geqs);
+ allocEQs = padEQs(in_eqs);
+ allocGEQs = padGEQs(in_geqs);
+ assert(allocEQs > 0 && allocGEQs > 0);
+ EQs = new eqn[allocEQs];
+ GEQs = new eqn[allocGEQs];
+ nVars = 0;
+ hashVersion = omega::hashVersion;
+ variablesInitialized = 0;
+ variablesFreed = 0;
+ varsOfInterest = 0;
+ safeVars = 0;
+ nEQs = 0;
+ nGEQs = 0;
+ nSUBs = 0;
+ nMemories = 0;
+ isTemporary = false;
+}
+
+Problem::Problem(const Problem & p2) {
+ allocEQs = padEQs(p2.nEQs); // Don't over-allocate; p2 might have too many!
+ allocGEQs = padGEQs(p2.nGEQs);
+ assert(allocEQs > 0 && allocGEQs > 0);
+ EQs = new eqn[allocEQs];
+ GEQs = new eqn[allocGEQs];
+ int e, i;
+ nVars = p2.nVars;
+ hashVersion = p2.hashVersion;
+ variablesInitialized = p2.variablesInitialized;
+ variablesFreed = p2.variablesFreed;
+ varsOfInterest = p2.varsOfInterest;
+ safeVars = p2.safeVars;
+ nEQs = p2.nEQs;
+ isTemporary = p2.isTemporary;
+ //nSUBs = 0;
+ for (e = p2.nEQs - 1; e >= 0; e--)
+ eqnncpy(&(EQs[e]), &(p2.EQs[e]), p2.nVars);
+ nGEQs = p2.nGEQs;
+ for (e = p2.nGEQs - 1; e >= 0; e--)
+ eqnncpy(&(GEQs[e]), &(p2.GEQs[e]), p2.nVars);
+ for (i = 0; i <= p2.nVars; i++)
+ var[i] = p2.var[i];
+ for (i = 0; i <= maxVars; i++)
+ forwardingAddress[i] = p2.forwardingAddress[i];
+ //nMemories = 0;
+ get_var_name = p2.get_var_name;
+ getVarNameArgs = p2.getVarNameArgs;
+}
+
+Problem & Problem::operator=(const Problem & p2) {
+ if (this != &p2) {
+ if(allocEQs < p2.nEQs) {
+ delete[] EQs;
+ allocEQs = padEQs(p2.nEQs);
+ EQs = new eqn[allocEQs];
+ }
+ if(allocGEQs < p2.nGEQs) {
+ delete[] GEQs;
+ allocGEQs = padGEQs(p2.nGEQs);
+ GEQs = new eqn[allocGEQs];
+ }
+ int e, i;
+ nVars = p2.nVars;
+ hashVersion = p2.hashVersion;
+ variablesInitialized = p2.variablesInitialized;
+ variablesFreed = p2.variablesFreed;
+ varsOfInterest = p2.varsOfInterest;
+ safeVars = p2.safeVars;
+ nEQs = p2.nEQs;
+ isTemporary = p2.isTemporary;
+ //nSUBs = 0;
+ for (e = p2.nEQs - 1; e >= 0; e--)
+ eqnncpy(&(EQs[e]), &(p2.EQs[e]), p2.nVars);
+ nGEQs = p2.nGEQs;
+ for (e = p2.nGEQs - 1; e >= 0; e--)
+ eqnncpy(&(GEQs[e]), &(p2.GEQs[e]), p2.nVars);
+ for (i = 0; i <= p2.nVars; i++)
+ var[i] = p2.var[i];
+ for (i = 0; i <= maxVars; i++)
+ forwardingAddress[i] = p2.forwardingAddress[i];
+ //nMemories = 0;
+ get_var_name = p2.get_var_name;
+ getVarNameArgs = p2.getVarNameArgs;
+ }
+ return *this;
+}
+
+
+void Problem::zeroVariable(int i) {
+ int e;
+ for (e = nGEQs - 1; e >= 0; e--) GEQs[e].coef[i] = 0;
+ for (e = nEQs - 1; e >= 0; e--) EQs[e].coef[i] = 0;
+ for (e = nSUBs - 1; e >= 0; e--) SUBs[e].coef[i] = 0;
+}
+
+/* Functions for allocating EQ's and GEQ's */
+
+int Problem::newGEQ() {
+ if (++nGEQs > allocGEQs) {
+ check_number_GEQs(nGEQs);
+ allocGEQs = padGEQs(allocGEQs, nGEQs);
+ assert(allocGEQs >= nGEQs);
+ eqn *new_geqs = new eqn[allocGEQs];
+ for (int e = nGEQs - 2; e >= 0; e--)
+ eqnncpy(&(new_geqs[e]), &(GEQs[e]), nVars);
+ delete[] GEQs;
+ GEQs = new_geqs;
+ }
+// problem->GEQs[nGEQs-1].color = black;
+// eqnnzero(&problem->GEQs[nGEQs-1],problem->nVars);
+ return nGEQs-1;
+}
+
+int Problem::newEQ() {
+ if (++nEQs > allocEQs) {
+ check_number_EQs(nEQs);
+ allocEQs = padEQs(allocEQs, nEQs);
+ assert(allocEQs >= nEQs);
+ eqn *new_eqs = new eqn[allocEQs];
+ for (int e = nEQs - 2; e >= 0; e--)
+ eqnncpy(&(new_eqs[e]), &(EQs[e]), nVars);
+ delete[] EQs;
+ EQs = new_eqs;
+ }
+// Could do this here, but some calls to newEQ do a copy instead of a zero;
+// problem->EQs[nEQs-1].color = black;
+// eqnnzero(&problem->EQs[nEQs-1],problem->nVars);
+ return nEQs-1;
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_query.cc b/omegalib/omega/src/omega_core/oc_query.cc
new file mode 100644
index 0000000..528b297
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_query.cc
@@ -0,0 +1,478 @@
+#include <omega/omega_core/oc_i.h>
+
+namespace omega {
+
+void Problem::unprotectVariable( int v) {
+ int e, j, i;
+ coef_t t;
+ i = forwardingAddress[v];
+ if (i < 0) {
+ i = -1 - i;
+ nSUBs--;
+ if (i < nSUBs) {
+ eqnncpy(&SUBs[i], &SUBs[nSUBs], nVars);
+ forwardingAddress[SUBs[i].key] = -i - 1;
+ }
+ }
+ else {
+ int bringToLife[maxVars];
+ int comingBack = 0;
+ int e2;
+ for (e = nSUBs - 1; e >= 0; e--)
+ if ((bringToLife[e] = (SUBs[e].coef[i] != 0)))
+ comingBack++;
+
+ for (e2 = nSUBs - 1; e2 >= 0; e2--)
+ if (bringToLife[e2]) {
+
+ nVars++;
+ safeVars++;
+ if (safeVars < nVars) {
+ for (e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].coef[nVars] = GEQs[e].coef[safeVars];
+ GEQs[e].coef[safeVars] = 0;
+ }
+ for (e = nEQs - 1; e >= 0; e--) {
+ EQs[e].coef[nVars] = EQs[e].coef[safeVars];
+ EQs[e].coef[safeVars] = 0;
+ }
+ for (e = nSUBs - 1; e >= 0; e--) {
+ SUBs[e].coef[nVars] = SUBs[e].coef[safeVars];
+ SUBs[e].coef[safeVars] = 0;
+ }
+ var[nVars] = var[safeVars];
+ forwardingAddress[var[nVars]] = nVars;
+ }
+ else {
+ for (e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].coef[safeVars] = 0;
+ }
+ for (e = nEQs - 1; e >= 0; e--) {
+ EQs[e].coef[safeVars] = 0;
+ }
+ for (e = nSUBs - 1; e >= 0; e--) {
+ SUBs[e].coef[safeVars] = 0;
+ }
+ }
+
+ var[safeVars] = SUBs[e2].key;
+ forwardingAddress[SUBs[e2].key] = safeVars;
+
+ int neweq = newEQ();
+ eqnncpy(&(EQs[neweq]), &(SUBs[e2]), nVars);
+ EQs[neweq].coef[safeVars] = -1;
+ if (e2 < nSUBs - 1)
+ eqnncpy(&(SUBs[e2]), &(SUBs[nSUBs - 1]), nVars);
+ nSUBs--;
+ }
+
+ if (i < safeVars) {
+ j = safeVars;
+ for (e = nSUBs - 1; e >= 0; e--) {
+ t = SUBs[e].coef[j];
+ SUBs[e].coef[j] = SUBs[e].coef[i];
+ SUBs[e].coef[i] = t;
+ }
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[j] != GEQs[e].coef[i]) {
+ GEQs[e].touched = true;
+ t = GEQs[e].coef[j];
+ GEQs[e].coef[j] = GEQs[e].coef[i];
+ GEQs[e].coef[i] = t;
+ }
+ for (e = nEQs - 1; e >= 0; e--) {
+ t = EQs[e].coef[j];
+ EQs[e].coef[j] = EQs[e].coef[i];
+ EQs[e].coef[i] = t;
+ }
+ {
+ short t;
+ t = var[j];
+ var[j] = var[i];
+ var[i] = t;
+ }
+ forwardingAddress[var[i]] = i;
+ forwardingAddress[var[j]] = j;
+ }
+ safeVars--;
+ }
+ chainUnprotect();
+}
+
+void Problem::constrainVariableSign( int color, int i, int sign) {
+ int nV = nVars;
+ int e, k, j;
+
+ k = forwardingAddress[i];
+ if (k < 0) {
+ k = -1 - k;
+
+ if (sign != 0) {
+ e = newGEQ();
+ eqnncpy(&GEQs[e], &SUBs[k], nVars);
+ for (j = 0; j <= nV; j++)
+ GEQs[e].coef[j] *= sign;
+ GEQs[e].coef[0]--;
+ GEQs[e].touched = 1;
+ GEQs[e].color = color;
+ }
+ else {
+ e = newEQ();
+ eqnncpy(&EQs[e], &SUBs[k], nVars);
+ EQs[e].color = color;
+ }
+ }
+ else if (sign != 0) {
+ e = newGEQ();
+ eqnnzero(&GEQs[e], nVars);
+ GEQs[e].coef[k] = sign;
+ GEQs[e].coef[0] = -1;
+ GEQs[e].touched = 1;
+ GEQs[e].color = color;
+ }
+ else {
+ e = newEQ();
+ eqnnzero(&EQs[e], nVars);
+ EQs[e].coef[k] = 1;
+ EQs[e].color = color;
+ }
+ /*
+ unprotectVariable(i);
+ return (simplifyProblem(0,1,0));
+ */
+}
+
+void Problem::constrainVariableValue( int color, int i, int value) {
+ int e, k;
+
+ k = forwardingAddress[i];
+ if (k < 0) {
+ k = -1 - k;
+
+ e = newEQ();
+ eqnncpy(&EQs[e], &SUBs[k], nVars);
+ EQs[e].coef[0] -= value;
+
+ }
+ else {
+ e = newEQ();
+ eqnnzero(&EQs[e], nVars);
+ EQs[e].coef[k] = 1;
+ EQs[e].coef[0] = -value;
+ }
+ EQs[e].color = color;
+}
+
+// Analyze v1-v2
+void Problem:: query_difference(int v1, int v2, coef_t &lowerBound, coef_t &upperBound, bool &guaranteed) {
+ int nV = nVars;
+ int e,i,e2;
+
+ coef_t lb1,ub1;
+ coef_t lb2,ub2;
+ assert(nSUBs == 0);
+ lowerBound = negInfinity;
+ lb1 = lb2 = negInfinity;
+ upperBound = posInfinity;
+ ub1 = ub2 = posInfinity;
+ guaranteed = true;
+ for (e = nEQs - 1; e >= 0; e--) {
+ if (EQs[e].coef[v1] == 0 && EQs[e].coef[v2] == 0)
+ continue;
+ for(i=nV;i>0;i--)
+ if (EQs[e].coef[i] && i!=v1 && i != v2) {
+ break;
+ }
+ if (i != 0) {
+ if (i > safeVars) {
+ // check to see if this variable appears anywhere else
+ for(e2 = nEQs-1; e2>=0;e2--) if (e != e2 && EQs[e2].coef[i]) break;
+ if (e2 < 0)
+ for(e2 = nGEQs-1; e2>=0;e2--) if (e != e2 && GEQs[e2].coef[i]) break;
+ if (e2 < 0)
+ for(e2 = nSUBs-1; e2>=0;e2--) if (e != e2 && SUBs[e2].coef[i]) break;
+ if (e2 >= 0) guaranteed = false;
+ }
+ else guaranteed = false;
+ continue;
+ }
+ if (EQs[e].coef[v1]*EQs[e].coef[v2] == -1) {
+ // found exact difference
+ coef_t d = - EQs[e].coef[v1] * EQs[e].coef[0];
+ set_max(lowerBound, d);
+ set_min(upperBound, d);
+ guaranteed =true;
+ return;
+ }
+ else if (EQs[e].coef[v1] == 0)
+ lb2 = ub2 = -EQs[e].coef[0]/ EQs[e].coef[v2];
+ else if (EQs[e].coef[v2] == 0)
+ lb1 = ub1 = -EQs[e].coef[0]/ EQs[e].coef[v1];
+ else guaranteed = false;
+ }
+
+ bool isDead[maxmaxGEQs];
+
+ for (e = nGEQs - 1; e >= 0; e--) isDead[e] = false;
+ int tryAgain = 1;
+ while (tryAgain) {
+ tryAgain = 0;
+ for (i = nVars; i > 0;i--)
+ if (i!= v1 && i != v2) {
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (!isDead[e] && GEQs[e].coef[i])
+ break;
+ if (e < 0)
+ e2 = e;
+ else if (GEQs[e].coef[i] > 0) {
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (!isDead[e2] && GEQs[e2].coef[i] < 0)
+ break;
+ }
+ else {
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (!isDead[e2] && GEQs[e2].coef[i] > 0)
+ break;
+ }
+ if (e2 < 0) {
+ int e3;
+ for (e3 = nSUBs - 1; e3 >= 0; e3--)
+ if (SUBs[e3].coef[i])
+ break;
+ if (e3 >= 0)
+ continue;
+ for (e3 = nEQs - 1; e3 >= 0; e3--)
+ if (EQs[e3].coef[i])
+ break;
+ if (e3 >= 0)
+ continue;
+ if (e >= 0) {
+ isDead[e] = true;
+ for (e--; e >= 0; e--)
+ if (GEQs[e].coef[i]) isDead[e] = true;
+ }
+ }
+ }
+ }
+
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (!isDead[e]) {
+ if (GEQs[e].coef[v1] == 0 && GEQs[e].coef[v2] == 0)
+ continue;
+ for(i=nV;i>0;i--)
+ if (GEQs[e].coef[i] && i!=v1 && i != v2)
+ break;
+ if (i != 0) {
+ guaranteed = false;
+ continue;
+ }
+ if (GEQs[e].coef[v1]*GEQs[e].coef[v2] == -1) {
+ // found relative difference
+ if (GEQs[e].coef[v1] == 1) {
+ // v1 - v2 + c >= 0
+ set_max(lowerBound, - GEQs[e].coef[0]);
+ }
+ else {
+ // v2 - v1 + c >= 0
+ // c >= v1-v2
+ set_min(upperBound, GEQs[e].coef[0]);
+ }
+ }
+ else if (GEQs[e].coef[v1] == 0 && GEQs[e].coef[v2] > 0)
+ lb2 = -GEQs[e].coef[0]/ GEQs[e].coef[v2];
+ else if (GEQs[e].coef[v1] == 0 && GEQs[e].coef[v2] < 0)
+ ub2 = -GEQs[e].coef[0]/ GEQs[e].coef[v2];
+ else if (GEQs[e].coef[v2] == 0 && GEQs[e].coef[v1] > 0)
+ lb1 = -GEQs[e].coef[0]/ GEQs[e].coef[v1];
+ else if (GEQs[e].coef[v2] == 0 && GEQs[e].coef[v1] < 0)
+ ub1 = -GEQs[e].coef[0]/ GEQs[e].coef[v1];
+ else guaranteed = false;
+ }
+
+ // ub1-lb2 >= v1-v2 >= lb1-ub2
+
+ if (negInfinity < lb2 && ub1 < posInfinity) set_min(upperBound, ub1-lb2);
+ if (negInfinity < lb1 && ub2 < posInfinity) set_max(lowerBound, lb1-ub2);
+ if (lowerBound >= upperBound) guaranteed = 1;
+}
+
+
+int Problem::queryVariable(int i, coef_t *lowerBound, coef_t *upperBound) {
+ int nV = nVars;
+ int e, j;
+ int isSimple;
+ int coupled = false;
+ for(j=1;j<=safeVars;j++)
+ if (var[j] > 0)
+ assert(forwardingAddress[var[j]] == j);
+
+ assert(i > 0);
+ i = forwardingAddress[i];
+ assert(i != 0);
+
+ (*lowerBound) = negInfinity;
+ (*upperBound) = posInfinity;
+
+ if (i < 0) {
+ int easy = true;
+ i = -i - 1;
+ for (j = 1; j <= nV; j++)
+ if (SUBs[i].coef[j] != 0)
+ easy = false;
+ if (easy) {
+ *upperBound = *lowerBound = SUBs[i].coef[0];
+ return (false);
+ }
+ return (true);
+ }
+
+ for (e = nSUBs - 1; e >= 0; e--)
+ if (SUBs[e].coef[i] != 0)
+ coupled = true;
+
+ for (e = nEQs - 1; e >= 0; e--)
+ if (EQs[e].coef[i] != 0) {
+ isSimple = true;
+ for (j = 1; j <= nV; j++)
+ if (i != j && EQs[e].coef[j] != 0) {
+ isSimple = false;
+ coupled = true;
+ break;
+ }
+ if (!isSimple)
+ continue;
+ else {
+ *lowerBound = *upperBound = -EQs[e].coef[i] * EQs[e].coef[0];
+ return (false);
+ }
+ }
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[i] != 0) {
+ if (GEQs[e].key == i) {
+ set_max(*lowerBound, -GEQs[e].coef[0]);
+ }
+ else if (GEQs[e].key == -i) {
+ set_min(*upperBound, GEQs[e].coef[0]);
+ }
+ else
+ coupled = true;
+ }
+ return (coupled);
+}
+
+int Problem::query_variable_bounds(int i, coef_t *l, coef_t *u) {
+ int coupled;
+ *l = negInfinity;
+ *u = posInfinity;
+ coupled = queryVariable(i, l, u);
+ if (!coupled || (nVars == 1 && forwardingAddress[i] == 1))
+ return 0;
+ if (abs(forwardingAddress[i]) == 1 && nVars + nSUBs == 2 && nEQs + nSUBs == 1) {
+ int couldBeZero;
+ queryCoupledVariable(i, l, u, &couldBeZero, negInfinity, posInfinity);
+ return 0;
+ }
+ return 1;
+}
+
+void Problem::queryCoupledVariable(int i, coef_t *l, coef_t *u, int *couldBeZero, coef_t lowerBound, coef_t upperBound) {
+ int e;
+ coef_t b1, b2;
+ const eqn *eqn;
+ coef_t sign;
+ int v;
+
+ if (abs(forwardingAddress[i]) != 1 || nVars + nSUBs != 2 || nEQs + nSUBs != 1) {
+ fprintf(outputFile, "queryCoupledVariablecalled with bad parameters\n");
+ printProblem();
+ exit(2);
+ }
+
+ if (forwardingAddress[i] == -1) {
+ eqn = &SUBs[0];
+ sign = 1;
+ v = 1;
+ }
+ else {
+ eqn = &EQs[0];
+ sign = -eqn->coef[1];
+ v = 2;
+ }
+
+ /* Variable i is defined in terms of variable v */
+
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[v] != 0) {
+ if (GEQs[e].coef[v] == 1) {
+ set_max(lowerBound, -GEQs[e].coef[0]);
+ }
+ else {
+ set_min(upperBound, GEQs[e].coef[0]);
+ }
+ }
+ /* lowerBound and upperBound are bounds on the value of v */
+
+ if (lowerBound > upperBound) {
+ *l = posInfinity;
+ *u = negInfinity;
+ *couldBeZero = 0;
+ return;
+ }
+ if (lowerBound == negInfinity) {
+ if (eqn->coef[v] > 0)
+ b1 = sign * negInfinity;
+ else
+ b1 = -sign * negInfinity;
+ }
+ else
+ b1 = sign * (eqn->coef[0] + eqn->coef[v] * lowerBound);
+ if (upperBound == posInfinity) {
+ if (eqn->coef[v] > 0)
+ b2 = sign * posInfinity;
+ else
+ b2 = -sign * posInfinity;
+ }
+ else
+ b2 = sign * (eqn->coef[0] + eqn->coef[v] * upperBound);
+
+ /* b1 and b2 are bounds on the value of i (don't know which is upper bound) */
+ if (b1 <= b2) {
+ set_max(*l, b1);
+ set_min(*u, b2);
+ }
+ else {
+ set_max(*l, b2);
+ set_min(*u, b1);
+ }
+ *couldBeZero = *l <= 0 && 0 <= *u && int_mod(eqn->coef[0], abs(eqn->coef[v])) == 0;
+}
+
+
+int Problem::queryVariableSigns(int i, int dd_lt, int dd_eq, int dd_gt, coef_t lowerBound, coef_t upperBound, bool *distKnown, coef_t *dist) {
+ int result;
+ coef_t l, u;
+ int couldBeZero;
+
+ l = negInfinity;
+ u = posInfinity;
+
+ queryVariable(i, &l, &u);
+ queryCoupledVariable(i, &l, &u, &couldBeZero, lowerBound, upperBound);
+ result = 0;
+ if (l < 0)
+ result |= dd_gt;
+ if (u > 0)
+ result |= dd_lt;
+ if (couldBeZero)
+ result |= dd_eq;
+ if (l == u) {
+ *distKnown = 1;
+ *dist = l;
+ }
+ else {
+ *distKnown = 0;
+ }
+ return (result);
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_quick_kill.cc b/omegalib/omega/src/omega_core/oc_quick_kill.cc
new file mode 100644
index 0000000..1b988d4
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_quick_kill.cc
@@ -0,0 +1,775 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+ Quick inequality elimination.
+
+ Notes:
+
+ History:
+ 03/31/09 Use BoolSet, Chun Chen
+*****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+#include <vector>
+#include <algorithm>
+#include <basic/BoolSet.h>
+
+namespace omega {
+
+int Problem::combineToTighten() {
+ int effort = min(12+5*(nVars-safeVars),23);
+
+ if (DBUG) {
+ fprintf(outputFile, "\nin combineToTighten (%d,%d):\n",effort,nGEQs);
+ printProblem();
+ fprintf(outputFile, "\n");
+ }
+ if (nGEQs > effort) {
+ if (TRACE) {
+ fprintf(outputFile, "too complicated to tighten\n");
+ }
+ return 1;
+ }
+
+ for(int e = 1; e < nGEQs; e++) {
+ for(int e2 = 0; e2 < e; e2++) {
+ coef_t g = 0;
+
+ bool has_wildcard = false;
+ bool has_wildcard2 = false;
+ for (int i = nVars; i > safeVars; i--) {
+ coef_t a = GEQs[e].coef[i];
+ coef_t b = GEQs[e2].coef[i];
+ g = gcd(g, abs(a+b));
+ if (a != 0)
+ has_wildcard = true;
+ if (b != 0)
+ has_wildcard2 = true;
+ }
+
+ coef_t c, c2;
+ if ((has_wildcard && !has_wildcard2) || (!has_wildcard && has_wildcard2))
+ c = 0;
+ else
+ c = -1;
+ for (int i = safeVars; i >= 1; i--) {
+ coef_t a = GEQs[e].coef[i];
+ coef_t b = GEQs[e2].coef[i];
+ if (a != 0 || b != 0) {
+ g = gcd(g, abs(a+b));
+
+ if (c < 0) {
+ if (g == 1)
+ break;
+ }
+ else if ((a>0 && b<0) || (a<0 && b>0)) {
+ if (c == 0) {
+ try {
+ coef_t prod = lcm(abs(a), abs(b));
+ c = prod/abs(a);
+ c2 = prod/abs(b);
+ }
+ catch (std::overflow_error) {
+ c = -1;
+ }
+ }
+ else {
+ if (c*a+c2*b != 0)
+ c = -1;
+ }
+ }
+ else {
+ c = -1;
+ }
+ }
+ }
+
+ bool done_unit_combine = false;
+ if (g > 1 && (GEQs[e].coef[0] + GEQs[e2].coef[0]) % g != 0) {
+ int e3 = newGEQ();
+ for(int i = nVars; i >= 1; i--) {
+ GEQs[e3].coef[i] = (GEQs[e].coef[i] + GEQs[e2].coef[i])/g;
+ }
+ GEQs[e3].coef[0] = int_div(GEQs[e].coef[0] + GEQs[e2].coef[0], g);
+ GEQs[e3].color = GEQs[e].color || GEQs[e2].color;
+ GEQs[e3].touched = 1;
+ if (DBUG) {
+ fprintf(outputFile, "Combined ");
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile,"\n and ");
+ printGEQ(&GEQs[e2]);
+ fprintf(outputFile,"\n to get #%d: ",e3);
+ printGEQ(&GEQs[e3]);
+ fprintf(outputFile,"\n\n");
+ }
+
+ done_unit_combine = true;
+ if (nGEQs > effort+5 || nGEQs > maxmaxGEQs-10) goto doneCombining;
+ }
+
+ if (c > 0 && !(c == 1 && c2 == 1 && done_unit_combine)) {
+ bool still_has_wildcard = false;
+ coef_t p[nVars-safeVars];
+ for (int i = nVars; i > safeVars; i--) {
+ p[i-safeVars-1] = c * GEQs[e].coef[i] + c2 * GEQs[e2].coef[i];
+ if (p[i-safeVars-1] != 0)
+ still_has_wildcard = true;
+ }
+ if (still_has_wildcard) {
+ int e3 = newGEQ();
+ for(int i = nVars; i > safeVars; i--)
+ GEQs[e3].coef[i] = p[i-safeVars-1];
+ for (int i = safeVars; i > 0; i--)
+ GEQs[e3].coef[i] = 0;
+ GEQs[e3].coef[0] = c * GEQs[e].coef[0] + c2 * GEQs[e2].coef[0];
+ GEQs[e3].color = GEQs[e].color || GEQs[e2].color;
+ GEQs[e3].touched = 1;
+ if (DBUG) {
+ fprintf(outputFile, "Additionally combined ");
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile,"\n and ");
+ printGEQ(&GEQs[e2]);
+ fprintf(outputFile,"\n to get #%d: ",e3);
+ printGEQ(&GEQs[e3]);
+ fprintf(outputFile,"\n\n");
+ }
+
+ if (nGEQs > effort+5 || nGEQs > maxmaxGEQs-10) goto doneCombining;
+ }
+ }
+ }
+ }
+
+doneCombining:
+ if (normalize() == normalize_false) return 0;
+ while (nEQs) {
+ if (!solveEQ()) return 0;
+ if (normalize() == normalize_false) return 0;
+ }
+ return 1;
+}
+
+
+void Problem::noteEssential(int onlyWildcards) {
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].essential = 0;
+ GEQs[e].varCount = 0;
+ }
+ if (onlyWildcards) {
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].essential = 1;
+ for (int i = nVars; i > safeVars; i--)
+ if (GEQs[e].coef[i] < -1 || GEQs[e].coef[i] > 1) {
+ GEQs[e].essential = 0;
+ break;
+ }
+ }
+ }
+ for (int i = nVars; i >= 1; i--) {
+ int onlyLB = -1;
+ int onlyUB = -1;
+ for (int e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[i] > 0) {
+ GEQs[e].varCount ++;
+ if (onlyLB == -1) onlyLB = e;
+ else onlyLB = -2;
+ }
+ else if (GEQs[e].coef[i] < 0) {
+ GEQs[e].varCount ++;
+ if (onlyUB == -1) onlyUB = e;
+ else onlyUB = -2;
+ }
+ if (onlyUB >= 0) {
+ if (DBUG) {
+ fprintf(outputFile,"only UB: ");
+ printGEQ(&GEQs[onlyUB]);
+ fprintf(outputFile,"\n");
+ }
+ GEQs[onlyUB].essential = 1;
+ }
+ if (onlyLB >= 0) {
+ if (DBUG) {
+ fprintf(outputFile,"only LB: ");
+ printGEQ(&GEQs[onlyLB]);
+ fprintf(outputFile,"\n");
+ }
+ GEQs[onlyLB].essential = 1;
+ }
+ }
+ for (int e = nGEQs - 1; e >= 0; e--)
+ if (!GEQs[e].essential && GEQs[e].varCount > 1) {
+ int i1,i2,i3;
+ for (i1 = nVars; i1 >= 1; i1--) if (GEQs[e].coef[i1]) break;
+ for (i2 = i1-1; i2 >= 1; i2--) if (GEQs[e].coef[i2]) break;
+ for (i3 = i2-1; i3 >= 1; i3--) if (GEQs[e].coef[i3]) break;
+ assert(i2 >= 1);
+ int e2;
+ for (e2 = nGEQs - 1; e2 >= 0; e2--)
+ if (e!=e2) {
+ coef_t crossProduct;
+ crossProduct = GEQs[e].coef[i1]*GEQs[e2].coef[i1];
+ crossProduct += GEQs[e].coef[i2]*GEQs[e2].coef[i2];
+ for (int i = i3; i >= 1; i--)
+ if (GEQs[e2].coef[i])
+ crossProduct += GEQs[e].coef[i]*GEQs[e2].coef[i];
+ if (crossProduct > 0) {
+ if (DBUG) fprintf(outputFile,"Cross product of %d and %d is " coef_fmt "\n", e, e2, crossProduct);
+ break;
+ }
+ }
+ if (e2 < 0) GEQs[e].essential = 1;
+ }
+ if (DBUG) {
+ fprintf(outputFile,"Computed essential equations\n");
+ fprintf(outputFile,"essential equations:\n");
+ for (int e = 0; e < nGEQs; e++)
+ if (GEQs[e].essential) {
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile,"\n");
+ }
+ fprintf(outputFile,"potentially redundant equations:\n");
+ for (int e = 0; e < nGEQs; e++)
+ if (!GEQs[e].essential) {
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile,"\n");
+ }
+ }
+}
+
+
+int Problem::findDifference(int e, int &v1, int &v2) {
+ // if 1 returned, eqn E is of form v1 -coef >= v2
+ for(v1=1;v1<=nVars;v1++)
+ if (GEQs[e].coef[v1]) break;
+ for(v2=v1+1;v2<=nVars;v2++)
+ if (GEQs[e].coef[v2]) break;
+ if (v2 > nVars) {
+ if (GEQs[e].coef[v1] == -1) {
+ v2 = v1;
+ v1 = 0;
+ return 1;
+ }
+ if (GEQs[e].coef[v1] == 1) {
+ v2 = 0;
+ return 1;
+ }
+ return 0;
+ }
+ if (GEQs[e].coef[v1] * GEQs[e].coef[v2] != -1) return 0;
+ if (GEQs[e].coef[v1] < 0) std::swap(v1,v2);
+ return 1;
+}
+
+
+namespace {
+ struct succListStruct {
+ int num;
+ int notEssential;
+ int var[maxVars];
+ coef_t diff[maxVars];
+ int eqn[maxVars];
+ };
+}
+
+
+int Problem::chainKill(int color, int onlyWildcards) {
+ int v1,v2,e;
+ int essentialPred[maxVars];
+ int redundant[maxmaxGEQs];
+ int inChain[maxVars];
+ int goodStartingPoint[maxVars];
+ int tryToEliminate[maxmaxGEQs];
+ int triedDoubleKill = 0;
+
+ succListStruct succ[maxVars];
+
+restart:
+
+ int anyToKill = 0;
+ int anyKilled = 0;
+ int canHandle = 0;
+
+ for(v1=0;v1<=nVars;v1++) {
+ succ[v1].num = 0;
+ succ[v1].notEssential = 0;
+ goodStartingPoint[v1] = 0;
+ inChain[v1] = -1;
+ essentialPred[v1] = 0;
+ }
+
+ int essentialEquations = 0;
+ for (e = 0; e < nGEQs; e++) {
+ redundant[e] = 0;
+ tryToEliminate[e] = !GEQs[e].essential;
+ if (GEQs[e].essential) essentialEquations++;
+ if (color && !GEQs[e].color) tryToEliminate[e] = 0;
+ }
+ if (essentialEquations == nGEQs) return 0;
+ if (2*essentialEquations < nVars) return 1;
+
+ for (e = 0; e < nGEQs; e++)
+ if (tryToEliminate[e] && GEQs[e].varCount <= 2 && findDifference(e,v1,v2)) {
+ assert(v1 == 0 || GEQs[e].coef[v1] == 1);
+ assert(v2 == 0 || GEQs[e].coef[v2] == -1);
+ succ[v2].notEssential++;
+ int s = succ[v2].num++;
+ succ[v2].eqn[s] = e;
+ succ[v2].var[s] = v1;
+ succ[v2].diff[s] = -GEQs[e].coef[0];
+ goodStartingPoint[v2] = 1;
+ anyToKill++;
+ canHandle++;
+ }
+ if (!anyToKill) {
+ return canHandle < nGEQs;
+ }
+ for (e = 0; e < nGEQs; e++)
+ if (!tryToEliminate[e] && GEQs[e].varCount <= 2 && findDifference(e,v1,v2)) {
+ assert(v1 == 0 || GEQs[e].coef[v1] == 1);
+ assert(v2 == 0 || GEQs[e].coef[v2] == -1);
+ int s = succ[v2].num++;
+ essentialPred[v1]++;
+ succ[v2].eqn[s] = e;
+ succ[v2].var[s] = v1;
+ succ[v2].diff[s] = -GEQs[e].coef[0];
+ canHandle++;
+ }
+
+
+ if (DBUG) {
+ int s;
+ fprintf(outputFile,"In chainkill: [\n");
+ for(v1 = 0;v1<=nVars;v1++) {
+ fprintf(outputFile,"#%d <= %s: ",essentialPred[v1],variable(v1));
+ for(s=0;s<succ[v1].notEssential;s++)
+ fprintf(outputFile," %s(" coef_fmt ") ",variable(succ[v1].var[s]), succ[v1].diff[s]);
+ for(;s<succ[v1].num;s++)
+ fprintf(outputFile," %s[" coef_fmt "] ",variable(succ[v1].var[s]), succ[v1].diff[s]);
+ fprintf(outputFile,"\n");
+ }
+ }
+
+ for(;v1<=nVars;v1++)
+ if (succ[v1].num == 1 && succ[v1].notEssential == 1) {
+ succ[v1].notEssential--;
+ essentialPred[succ[v1].var[succ[v1].notEssential]]++;
+ }
+
+ if (DBUG) fprintf(outputFile,"Trying quick double kill:\n");
+ int s1a,s1b,s2;
+ int v3;
+ for(v1 = 0;v1<=nVars;v1++)
+ for(s1a=0;s1a<succ[v1].notEssential;s1a++) {
+ v3 = succ[v1].var[s1a];
+ for(s1b=0;s1b<succ[v1].num;s1b++)
+ if (s1a != s1b) {
+ v2 = succ[v1].var[s1b];
+ for(s2=0;s2<succ[v2].num;s2++)
+ if (succ[v2].var[s2] == v3 && succ[v1].diff[s1b] + succ[v2].diff[s2] >= succ[v1].diff[s1a]) {
+ if (DBUG) {
+ fprintf(outputFile,"quick double kill: ");
+ printGEQ(&GEQs[succ[v1].eqn[s1a]]);
+ fprintf(outputFile,"\n");
+ }
+ redundant[succ[v1].eqn[s1a]] = 1;
+ anyKilled++;
+ anyToKill--;
+ goto nextVictim;
+ }
+ }
+ nextVictim: v1 = v1;
+ }
+ if (anyKilled) {
+ for (e = nGEQs-1; e >= 0;e--)
+ if (redundant[e]) {
+ if (DBUG) {
+ fprintf(outputFile,"Deleting ");
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile,"\n");
+ }
+ deleteGEQ(e);
+ }
+
+ if (!anyToKill) return canHandle < nGEQs;
+ noteEssential(onlyWildcards);
+ triedDoubleKill = 1;
+ goto restart;
+ }
+
+ for(v1 = 0;v1<=nVars;v1++)
+ if (succ[v1].num == succ[v1].notEssential && succ[v1].notEssential > 0) {
+ succ[v1].notEssential--;
+ essentialPred[succ[v1].var[succ[v1].notEssential]]++;
+ }
+
+ while (1) {
+ int chainLength;
+ int chain[maxVars];
+ coef_t distance[maxVars];
+ // pick a place to start
+ for(v1 = 0;v1<=nVars;v1++)
+ if (essentialPred[v1] == 0 && succ[v1].num > succ[v1].notEssential)
+ break;
+ if (v1 > nVars)
+ for(v1 = 0;v1<=nVars;v1++)
+ if (goodStartingPoint[v1] && succ[v1].num > succ[v1].notEssential)
+ break;
+ if (v1 > nVars) break;
+
+ chainLength = 1;
+ chain[0] = v1;
+ distance[0] = 0;
+ inChain[v1] = 0;
+ int s;
+
+ while (succ[v1].num > succ[v1].notEssential) {
+ s = succ[v1].num-1;
+ if (inChain[succ[v1].var[s]] >= 0) {
+ // Found cycle, don't do anything with them yet
+ break;
+ }
+ succ[v1].num = s;
+
+ distance[chainLength]= distance[chainLength-1] + succ[v1].diff[s];
+ v1 = chain[chainLength] = succ[v1].var[s];
+ essentialPred[v1]--;
+ assert(essentialPred[v1] >= 0);
+ inChain[v1] = chainLength;
+ chainLength++;
+ }
+
+
+ int c;
+ if (DBUG) {
+ fprintf(outputFile,"Found chain: \n");
+ for (c = 0; c < chainLength; c++)
+ fprintf(outputFile,"%s:" coef_fmt " ",variable(chain[c]), distance[c]);
+ fprintf(outputFile,"\n");
+ }
+
+
+ for (c = 0; c < chainLength; c++) {
+ v1 = chain[c];
+ for(s=0;s<succ[v1].notEssential;s++) {
+ if (DBUG)
+ fprintf(outputFile,"checking for %s + " coef_fmt " <= %s \n", variable(v1), succ[v1].diff[s], variable(succ[v1].var[s]));
+ if (inChain[succ[v1].var[s]] > c+1) {
+ if (DBUG)
+ fprintf(outputFile,"%s + " coef_fmt " <= %s is in chain\n", variable(v1), distance[inChain[succ[v1].var[s]]]- distance[c], variable(succ[v1].var[s]));
+ if ( distance[inChain[succ[v1].var[s]]]- distance[c] >= succ[v1].diff[s]) {
+ if (DBUG)
+ fprintf(outputFile,"%s + " coef_fmt " <= %s is redundant\n", variable(v1),succ[v1].diff[s], variable(succ[v1].var[s]));
+ redundant[succ[v1].eqn[s]] = 1;
+ }
+ }
+ }
+ }
+ for (c = 0; c < chainLength; c++)
+ inChain[chain[c]] = -1;
+ }
+
+ for (e = nGEQs-1; e >= 0;e--)
+ if (redundant[e]) {
+ if (DBUG) {
+ fprintf(outputFile,"Deleting ");
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile,"\n");
+ }
+ deleteGEQ(e);
+ anyKilled = 1;
+ }
+
+ if (anyKilled) noteEssential(onlyWildcards);
+
+ if (anyKilled && DBUG) {
+ fprintf(outputFile,"\nResult:\n");
+ printProblem();
+ }
+ if (DBUG) {
+ fprintf(outputFile,"] end chainkill\n");
+ printProblem();
+ }
+ return canHandle < nGEQs;
+}
+
+
+namespace {
+ struct varCountStruct {
+ int e;
+ int safeVarCount;
+ int wildVarCount;
+ varCountStruct(int e_, int count1_, int count2_) {
+ e = e_;
+ safeVarCount = count1_;
+ wildVarCount = count2_; }
+ };
+ bool operator<(const varCountStruct &a, const varCountStruct &b) {
+ if (a.wildVarCount < b.wildVarCount)
+ return true;
+ else if (a.wildVarCount > b.wildVarCount)
+ return false;
+ else
+ return a.safeVarCount < b.safeVarCount;
+ }
+}
+
+
+//
+// Deduct redundant inequalities by combination of any two inequalities.
+// Return value: 0 (no solution),
+// 1 (nothing killed),
+// 2 (some inequality killed).
+//
+int Problem::quickKill(int onlyWildcards, bool desperate) {
+ if (!onlyWildcards && !combineToTighten())
+ return 0;
+ noteEssential(onlyWildcards);
+ int moreToDo = chainKill(0, onlyWildcards);
+
+#ifdef NDEBUG
+ if (!moreToDo) return 1;
+#endif
+
+
+ if (!desperate && nGEQs > 256) { // original 60, increased by chun
+ if (TRACE) {
+ fprintf(outputFile, "%d inequalities are too complicated to quick kill\n", nGEQs);
+ }
+ return 1;
+ }
+
+ if (DBUG) {
+ fprintf(outputFile, "in eliminate Redudant:\n");
+ printProblem();
+ }
+
+ int isDead[nGEQs];
+ std::vector<varCountStruct> killOrder;
+ std::vector<BoolSet<> > P(nGEQs, BoolSet<>(nVars)), Z(nGEQs, BoolSet<>(nVars)), N(nGEQs, BoolSet<>(nVars));
+ BoolSet<> PP, PZ, PN; // possible Positives, possible zeros & possible negatives
+
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ isDead[e] = 0;
+ int safeVarCount = 0;
+ int wildVarCount = 0;
+ for (int i = nVars; i >= 1; i--) {
+ if (GEQs[e].coef[i] == 0)
+ Z[e].set(i-1);
+ else {
+ if (i > safeVars)
+ wildVarCount++;
+ else
+ safeVarCount++;
+ if (GEQs[e].coef[i] < 0)
+ N[e].set(i-1);
+ else
+ P[e].set(i-1);
+ }
+ }
+
+ if (!GEQs[e].essential || wildVarCount > 0)
+ killOrder.push_back(varCountStruct(e, safeVarCount, wildVarCount));
+ }
+
+ sort(killOrder.begin(), killOrder.end());
+
+ if (DEBUG) {
+ fprintf(outputFile,"Prefered kill order:\n");
+ for (int e3I = killOrder.size()-1; e3I >= 0; e3I--) {
+ fprintf(outputFile,"%2d: ",nGEQs-1-e3I);
+ printGEQ(&GEQs[killOrder[e3I].e]);
+ fprintf(outputFile,"\n");
+ }
+ }
+
+ int e3U = killOrder.size()-1;
+ while (e3U >= 0) {
+ // each round of elimination is for inequalities of same complexity and rounds are at descending complexity order
+ int e3L = e3U-1;
+ for(; e3L >= 0; e3L--)
+ if (killOrder[e3L].safeVarCount+killOrder[e3L].wildVarCount != killOrder[e3U].safeVarCount + killOrder[e3U].wildVarCount)
+ break;
+
+ // check if e3 can be eliminated from combination of e1 and e2
+ for (int e1 = 0; e1 < nGEQs; e1++)
+ if (!isDead[e1])
+ for (int e2 = e1+1; e2 < nGEQs; e2++)
+ if (!isDead[e2]) {
+ coef_t alpha = 0;
+ int p, q;
+ for (p = nVars; p > 1; p--)
+ for (q = p - 1; q > 0; q--) {
+ try {
+ alpha = check_mul(GEQs[e1].coef[p], GEQs[e2].coef[q]) - check_mul(GEQs[e2].coef[p], GEQs[e1].coef[q]);
+ }
+ catch (std::overflow_error) {
+ continue;
+ }
+ if (alpha != 0)
+ goto foundPQ;
+ }
+ continue;
+
+ foundPQ:
+ PZ = (Z[e1] & Z[e2]) | (P[e1] & N[e2]) | (N[e1] & P[e2]);
+ PP = P[e1] | P[e2];
+ PN = N[e1] | N[e2];
+ if (DEBUG) {
+ fprintf(outputFile,"Considering combination of ");
+ printGEQ(&(GEQs[e1]));
+ fprintf(outputFile," and ");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile,"\n");
+ }
+
+ for (int e3I = e3U; e3I > e3L; e3I--) {
+ int e3 = killOrder[e3I].e;
+ if (!isDead[e3] && e3 != e1 && e3 != e2)
+ try {
+ coef_t alpha1, alpha2, alpha3;
+
+ if (!PZ.imply(Z[e3]))
+ goto nextE3;
+
+ alpha1 = check_mul(GEQs[e2].coef[q], GEQs[e3].coef[p]) - check_mul(GEQs[e2].coef[p], GEQs[e3].coef[q]);
+ alpha2 = -(check_mul(GEQs[e1].coef[q], GEQs[e3].coef[p]) - check_mul(GEQs[e1].coef[p], GEQs[e3].coef[q]));
+ alpha3 = alpha;
+
+ if (alpha1 < 0) {
+ alpha1 = -alpha1;
+ alpha2 = -alpha2;
+ alpha3 = -alpha3;
+ }
+ if (alpha1 == 0 || alpha2 <= 0)
+ goto nextE3;
+
+ {
+ coef_t g = gcd(gcd(alpha1, alpha2), abs(alpha3));
+ alpha1 /= g;
+ alpha2 /= g;
+ alpha3 /= g;
+ }
+
+ if (DEBUG) {
+ fprintf(outputFile, coef_fmt "e1 + " coef_fmt "e2 = " coef_fmt "e3: ",alpha1,alpha2,alpha3);
+ printGEQ(&(GEQs[e3]));
+ fprintf(outputFile,"\n");
+ }
+
+ if (alpha3 > 0) { // trying to prove e3 is redundant
+ if (!GEQs[e3].color && (GEQs[e1].color || GEQs[e2].color)) {
+ goto nextE3;
+ }
+ if (!PP.imply(P[e3]) | !PN.imply(N[e3]))
+ goto nextE3;
+
+ // verify alpha1*v1+alpha2*v2 = alpha3*v3
+ for (int k = nVars; k >= 1; k--)
+ if (check_mul(alpha3, GEQs[e3].coef[k]) != check_mul(alpha1, GEQs[e1].coef[k]) + check_mul(alpha2, GEQs[e2].coef[k]))
+ goto nextE3;
+
+ coef_t c = check_mul(alpha1, GEQs[e1].coef[0]) + check_mul(alpha2, GEQs[e2].coef[0]);
+ if (c < check_mul(alpha3, (GEQs[e3].coef[0] + 1))) {
+ if (DBUG) {
+ fprintf(outputFile, "found redundant inequality\n");
+ fprintf(outputFile, "alpha1, alpha2, alpha3 = " coef_fmt "," coef_fmt "," coef_fmt "\n", alpha1, alpha2, alpha3);
+ printGEQ(&(GEQs[e1]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n=> ");
+ printGEQ(&(GEQs[e3]));
+ fprintf(outputFile, "\n\n");
+ assert(moreToDo);
+ }
+
+ isDead[e3] = 1;
+ }
+ }
+ else { // trying to prove e3 <= 0 or e3 = 0
+ if (!PN.imply(P[e3]) | !PP.imply(N[e3]))
+ goto nextE3;
+
+ // verify alpha1*v1+alpha2*v2 = alpha3*v3
+ for (int k = nVars; k >= 1; k--)
+ if (check_mul(alpha3, GEQs[e3].coef[k]) != check_mul(alpha1, GEQs[e1].coef[k]) + check_mul(alpha2, GEQs[e2].coef[k]))
+ goto nextE3;
+
+ if (DEBUG) {
+ fprintf(outputFile,"All but constant term checked\n");
+ }
+ coef_t c = check_mul(alpha1, GEQs[e1].coef[0]) + check_mul(alpha2, GEQs[e2].coef[0]);
+ if (DEBUG) {
+ fprintf(outputFile,"All but constant term checked\n");
+ fprintf(outputFile,"Constant term is " coef_fmt " vs " coef_fmt "\n",
+ alpha3*GEQs[e3].coef[0],
+ alpha3*(GEQs[e3].coef[0]-1));
+ }
+ if (c < check_mul(alpha3, (GEQs[e3].coef[0]))) {
+ // we just proved e3 < 0, so no solutions exist
+ if (DBUG) {
+ fprintf(outputFile, "found implied over tight inequality\n");
+ fprintf(outputFile, "alpha1, alpha2, alpha3 = " coef_fmt "," coef_fmt "," coef_fmt "\n", alpha1, alpha2, -alpha3);
+ printGEQ(&(GEQs[e1]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n=> not ");
+ printGEQ(&(GEQs[e3]));
+ fprintf(outputFile, "\n\n");
+ }
+ return 0;
+ }
+ else if (!GEQs[e3].color && (GEQs[e1].color || GEQs[e2].color)) {
+ goto nextE3;
+ }
+ else if (c < check_mul(alpha3, (GEQs[e3].coef[0] - 1))) {
+ // we just proved e3 <= 0, so e3 = 0
+ if (DBUG) {
+ fprintf(outputFile, "found implied tight inequality\n");
+ fprintf(outputFile, "alpha1, alpha2, alpha3 = " coef_fmt "," coef_fmt "," coef_fmt "\n", alpha1, alpha2, -alpha3);
+ printGEQ(&(GEQs[e1]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n=> inverse ");
+ printGEQ(&(GEQs[e3]));
+ fprintf(outputFile, "\n\n");
+ }
+ int neweq = newEQ();
+ eqnncpy(&EQs[neweq], &GEQs[e3], nVars);
+ addingEqualityConstraint(neweq);
+ isDead[e3] = 1;
+ }
+ }
+ nextE3:;
+ }
+ catch (std::overflow_error) {
+ continue;
+ }
+ }
+ }
+
+ e3U = e3L;
+ }
+
+ bool anything_killed = false;
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ if (isDead[e]) {
+ anything_killed = true;
+ deleteGEQ(e);
+ }
+ }
+
+ if (DBUG) {
+ fprintf(outputFile,"\nResult:\n");
+ printProblem();
+ }
+
+ if (anything_killed)
+ return 2;
+ else
+ return 1;
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_simple.cc b/omegalib/omega/src/omega_core/oc_simple.cc
new file mode 100644
index 0000000..0e492db
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_simple.cc
@@ -0,0 +1,1373 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+ Support functions for solving a problem.
+
+ Notes:
+
+ History:
+ 10/13/08 Complete back substitution process, Chun Chen.
+ 05/28/09 Extend normalize process to handle redundancy involving
+ wilddcards, Chun Chen
+*****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+#include <basic/BoolSet.h>
+#include <algorithm>
+#include <vector>
+
+namespace omega {
+
+int checkIfSingleVar(eqn* e, int i) {
+ for (; i > 0; i--)
+ if (e->coef[i]) {
+ i--;
+ break;
+ }
+ for (; i > 0; i--)
+ if (e->coef[i])
+ break;
+ return (i == 0);
+}
+
+
+int singleVarGEQ(eqn* e) {
+ return !e->touched && e->key != 0 && -maxVars <= e->key && e->key <= maxVars;
+}
+
+
+void checkVars(int nVars) {
+ if (nVars > maxVars) {
+ fprintf(stderr, "\nERROR:\n");
+ fprintf(stderr, "An attempt was made to create a conjunction with %d variables.\n", nVars);
+ fprintf(stderr, "The current limit on variables in a single conjunction is %d.\n", maxVars);
+ fprintf(stderr, "This limit can be changed by changing the #define of maxVars in oc.h.\n\n");
+ exit(2);
+ }
+}
+
+
+void Problem::difficulty(int &numberNZs, coef_t &maxMinAbsCoef, coef_t &sumMinAbsCoef) const {
+ numberNZs=0;
+ maxMinAbsCoef=0;
+ sumMinAbsCoef=0;
+ for (int e = 0; e < nGEQs; e++) {
+ coef_t maxCoef = 0;
+ for(int i = 1;i <= nVars;i++)
+ if (GEQs[e].coef[i]!=0) {
+ coef_t a = abs(GEQs[e].coef[i]);
+ maxCoef = max(maxCoef,a);
+ numberNZs++;
+ }
+ coef_t nextCoef = 0;
+ for(int i = 1;i <= nVars;i++)
+ if (GEQs[e].coef[i]!=0) {
+ coef_t a = abs(GEQs[e].coef[i]);
+ if (a < maxCoef) nextCoef = max(nextCoef,a);
+ else if (a == maxCoef) maxCoef = 0x7fffffff;
+ }
+ maxMinAbsCoef = max(maxMinAbsCoef,nextCoef);
+ sumMinAbsCoef += nextCoef;
+ }
+
+ for (int e = 0; e < nEQs; e++) {
+ coef_t maxCoef = 0;
+ for(int i = 1;i <= nVars;i++)
+ if (EQs[e].coef[i]!=0) {
+ coef_t a = abs(EQs[e].coef[i]);
+ maxCoef = max(maxCoef,a);
+ numberNZs++;
+ }
+ coef_t nextCoef = 0;
+ for(int i = 1;i <= nVars;i++)
+ if (EQs[e].coef[i]!=0) {
+ coef_t a = abs(EQs[e].coef[i]);
+ if (a < maxCoef) nextCoef = max(nextCoef,a);
+ else if (a == maxCoef) maxCoef = 0x7fffffff;
+ }
+ maxMinAbsCoef = max(maxMinAbsCoef,nextCoef);
+ sumMinAbsCoef += nextCoef;
+ }
+}
+
+int Problem::countRedGEQs() const {
+ int result = 0;
+ for (int e = 0; e < nGEQs; e++)
+ if (GEQs[e].color == EQ_RED) result++;
+ return result;
+}
+
+int Problem::countRedEQs() const {
+ int result = 0;
+ for (int e = 0; e < nEQs; e++)
+ if (EQs[e].color == EQ_RED) result++;
+ return result;
+}
+
+int Problem::countRedEquations() const {
+ int result = 0;
+ for (int e = 0; e < nEQs; e++)
+ if (EQs[e].color == EQ_RED) {
+ int i;
+ for (i = nVars; i > 0; i--) if (EQs[e].coef[i]) break;
+ if (i == 0 && EQs[e].coef[0] != 0) return 0;
+ else result+=2;
+ }
+ for (int e = 0; e < nGEQs; e++)
+ if (GEQs[e].color == EQ_RED) result+=1;
+ for (int e = 0; e < nMemories; e++)
+ switch(redMemory[e].kind ) {
+ case redEQ:
+ case redStride:
+ e++;
+ case redLEQ:
+ case redGEQ:
+ e++;
+ case notRed:
+ ; /* avoid warning about notRed not handled */
+ }
+ return result;
+}
+
+void Problem::deleteBlack() {
+ int RedVar[maxVars];
+ for(int i = safeVars+1;i <= nVars;i++) RedVar[i] = 0;
+
+ assert(nSUBs == 0);
+
+ for (int e = nEQs-1; e >= 0; e--)
+ if (EQs[e].color != EQ_RED) {
+ eqnncpy(&EQs[e],&EQs[nEQs-1], nVars);
+ nEQs--;
+ }
+ else
+ for(int i = safeVars+1;i <= nVars;i++)
+ if (EQs[e].coef[i]) RedVar[i] = 1;
+
+ for (int e = nGEQs-1; e >= 0; e--)
+ if (GEQs[e].color != EQ_RED) {
+ eqnncpy(&GEQs[e],&GEQs[nGEQs-1], nVars);
+ nGEQs--;
+ }
+ else
+ for(int i = safeVars+1;i <= nVars;i++)
+ if (GEQs[e].coef[i]) RedVar[i] = 1;
+
+ assert(nSUBs == 0);
+
+ for(int i = nVars; i > safeVars;i--) {
+ if (!RedVar[i]) deleteVariable(i);
+ }
+}
+
+
+void Problem::deleteRed() {
+ int BlackVar[maxVars];
+ for(int i = safeVars+1;i <= nVars;i++) BlackVar[i] = 0;
+
+ assert(nSUBs == 0);
+ for (int e = nEQs-1; e >=0; e--)
+ if (EQs[e].color) {
+ eqnncpy(&EQs[e],&EQs[nEQs-1], nVars);
+ nEQs--;
+ }
+ else
+ for(int i = safeVars+1;i <= nVars;i++)
+ if (EQs[e].coef[i]) BlackVar[i] = 1;
+
+ for (int e = nGEQs-1; e >=0; e--)
+ if (GEQs[e].color) {
+ eqnncpy(&GEQs[e],&GEQs[nGEQs-1], nVars);
+ nGEQs--;
+ }
+ else
+ for(int i = safeVars+1;i <= nVars;i++)
+ if (GEQs[e].coef[i]) BlackVar[i] = 1;
+
+ assert(nSUBs == 0);
+
+ for(int i = nVars; i> safeVars;i--) {
+ if (!BlackVar[i]) deleteVariable(i);
+ }
+}
+
+
+void Problem::turnRedBlack() {
+ for (int e = nEQs-1; e >= 0; e--) EQs[e].color = 0;
+ for (int e = nGEQs-1; e >= 0; e--) GEQs[e].color = 0;
+}
+
+
+void Problem::useWildNames() {
+ for(int i = safeVars+1; i <= nVars; i++) nameWildcard(i);
+}
+
+
+void negateCoefficients(eqn* eqn, int nVars) {
+ for (int i = nVars; i >= 0; i--)
+ eqn-> coef[i] = -eqn->coef[i];
+ eqn->touched = true;
+}
+
+
+void Problem::negateGEQ(int e) {
+ negateCoefficients(&GEQs[e],nVars);
+ GEQs[e].coef[0]--;
+}
+
+
+void Problem:: deleteVariable(int i) {
+ if (i < safeVars) {
+ int j = safeVars;
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].touched = true;
+ GEQs[e].coef[i] = GEQs[e].coef[j];
+ GEQs[e].coef[j] = GEQs[e].coef[nVars];
+ }
+ for (int e = nEQs - 1; e >= 0; e--) {
+ EQs[e].coef[i] = EQs[e].coef[j];
+ EQs[e].coef[j] = EQs[e].coef[nVars];
+ }
+ for (int e = nSUBs - 1; e >= 0; e--) {
+ SUBs[e].coef[i] = SUBs[e].coef[j];
+ SUBs[e].coef[j] = SUBs[e].coef[nVars];
+ }
+ var[i] = var[j];
+ var[j] = var[nVars];
+ }
+ else if (i < nVars) {
+ for (int e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[nVars]) {
+ GEQs[e].coef[i] = GEQs[e].coef[nVars];
+ GEQs[e].touched = true;
+ }
+ for (int e = nEQs - 1; e >= 0; e--)
+ EQs[e].coef[i] = EQs[e].coef[nVars];
+ for (int e = nSUBs - 1; e >= 0; e--)
+ SUBs[e].coef[i] = SUBs[e].coef[nVars];
+ var[i] = var[nVars];
+ }
+ if (i <= safeVars)
+ safeVars--;
+ nVars--;
+}
+
+
+void Problem::setInternals() {
+ if (!variablesInitialized) {
+ initializeVariables();
+ }
+
+ var[0] = 0;
+ nextWildcard = 0;
+ for(int i = 1;i <= nVars;i++)
+ if (var[i] < 0)
+ var[i] = --nextWildcard;
+
+ assert(nextWildcard >= -maxWildcards);
+
+ CHECK_FOR_DUPLICATE_VARIABLE_NAMES;
+
+ int v = nSUBs;
+ for(int i = 1;i <= safeVars;i++) if (var[i] > 0) v++;
+ varsOfInterest = v;
+
+ if (nextKey * 3 > maxKeys) {
+ omega::hashVersion++;
+ nextKey = maxVars + 1;
+ for (int e = nGEQs - 1; e >= 0; e--)
+ GEQs[e].touched = true;
+ for (int i = 0; i < hashTableSize; i++)
+ hashMaster[i].touched = -1;
+ hashVersion = omega::hashVersion;
+ }
+ else if (hashVersion != omega::hashVersion) {
+ for (int e = nGEQs - 1; e >= 0; e--)
+ GEQs[e].touched = true;
+ hashVersion = omega::hashVersion;
+ }
+}
+
+
+void Problem::setExternals() {
+ for (int i = 1; i <= safeVars; i++)
+ forwardingAddress[var[i]] = i;
+ for (int i = 0; i < nSUBs; i++)
+ forwardingAddress[SUBs[i].key] = -i - 1;
+}
+
+
+void setOutputFile(FILE * file) {
+ /* sets the file to which printProblem should send its output to "file" */
+
+ outputFile = file;
+}
+
+
+void setPrintLevel(int level) {
+ /* Sets the nber of points printed before constraints in printProblem */
+ headerLevel = level;
+}
+
+
+void Problem::putVariablesInStandardOrder() {
+ for(int i = 1;i <= safeVars;i++) {
+ int b = i;
+ for(int j=i+1;j<=safeVars;j++) {
+ if (var[b] < var[j]) b = j;
+ }
+ if (b != i) swapVars(i,b);
+ }
+}
+
+
+void Problem::nameWildcard(int i) {
+ int j;
+ do {
+ --nextWildcard;
+ if (nextWildcard < -maxWildcards)
+ nextWildcard = -1;
+ var[i] = nextWildcard;
+ for(j = nVars; j > 0;j--) if (i!=j && var[j] == nextWildcard) break;
+ } while (j != 0);
+}
+
+
+int Problem::protectWildcard(int i) {
+ assert (i > safeVars);
+ if (i != safeVars+1) swapVars(i,safeVars+1);
+ safeVars++;
+ nameWildcard(safeVars);
+ return safeVars;
+}
+
+
+int Problem::addNewProtectedWildcard() {
+ int i = ++safeVars;
+ nVars++;
+ if (nVars != i) {
+ for (int e = nGEQs - 1; e >= 0; e--) {
+ if (GEQs[e].coef[i] != 0)
+ GEQs[e].touched = true;
+ GEQs[e].coef[nVars] = GEQs[e].coef[i];
+ }
+ for (int e = nEQs - 1; e >= 0; e--) {
+ EQs[e].coef[nVars] = EQs[e].coef[i];
+ }
+ for (int e = nSUBs - 1; e >= 0; e--) {
+ SUBs[e].coef[nVars] = SUBs[e].coef[i];
+ }
+ var[nVars] = var[i];
+ }
+ for (int e = nGEQs - 1; e >= 0; e--)
+ GEQs[e].coef[i] = 0;
+ for (int e = nEQs - 1; e >= 0; e--)
+ EQs[e].coef[i] = 0;
+ for (int e = nSUBs - 1; e >= 0; e--)
+ SUBs[e].coef[i] = 0;
+ nameWildcard(i);
+ return (i);
+}
+
+
+int Problem::addNewUnprotectedWildcard() {
+ int i = ++nVars;
+ for (int e = nGEQs - 1; e >= 0; e--) GEQs[e].coef[i] = 0;
+ for (int e = nEQs - 1; e >= 0; e--) EQs[e].coef[i] = 0;
+ for (int e = nSUBs - 1; e >= 0; e--) SUBs[e].coef[i] = 0;
+ nameWildcard(i);
+ return i;
+}
+
+
+void Problem::cleanoutWildcards() {
+ bool renormalize = false;
+
+ // substituting wildcard equality
+ for (int e = nEQs-1; e >= 0; e--) {
+ for (int i = nVars; i >= safeVars+1; i--)
+ if (EQs[e].coef[i] != 0) {
+ coef_t c = EQs[e].coef[i];
+ coef_t a = abs(c);
+
+ bool preserveThisConstraint = true;
+ for (int e2 = nEQs-1; e2 >= 0; e2--)
+ if (e2 != e && EQs[e2].coef[i] != 0 && EQs[e2].color >= EQs[e].color) {
+ preserveThisConstraint = preserveThisConstraint && (gcd(a,abs(EQs[e2].coef[i])) != 1);
+ coef_t k = lcm(a, abs(EQs[e2].coef[i]));
+ coef_t coef1 = (EQs[e2].coef[i]>0?1:-1) * k / c;
+ coef_t coef2 = k / abs(EQs[e2].coef[i]);
+ for (int j = nVars; j >= 0; j--)
+ EQs[e2].coef[j] = EQs[e2].coef[j] * coef2 - EQs[e].coef[j] * coef1;
+
+ coef_t g = 0;
+ for (int j = nVars; j >= 0; j--) {
+ g = gcd(abs(EQs[e2].coef[j]), g);
+ if (g == 1)
+ break;
+ }
+ if (g != 0 && g != 1)
+ for (int j = nVars; j >= 0; j--)
+ EQs[e2].coef[j] /= g;
+ }
+
+ for (int e2 = nGEQs-1; e2 >= 0; e2--)
+ if (GEQs[e2].coef[i] != 0 && GEQs[e2].color >= EQs[e].color) {
+ coef_t k = lcm(a, abs(GEQs[e2].coef[i]));
+ coef_t coef1 = (GEQs[e2].coef[i]>0?1:-1) * k / c;
+ coef_t coef2 = k / abs(GEQs[e2].coef[i]);
+ for (int j = nVars; j >= 0; j--)
+ GEQs[e2].coef[j] = GEQs[e2].coef[j] * coef2 - EQs[e].coef[j] * coef1;
+
+ GEQs[e2].touched = 1;
+ renormalize = true;
+ }
+
+ for (int e2 = nSUBs-1; e2 >= 0; e2--)
+ if (SUBs[e2].coef[i] != 0 && SUBs[e2].color >= EQs[e].color) {
+ coef_t k = lcm(a, abs(SUBs[e2].coef[i]));
+ coef_t coef1 = (SUBs[e2].coef[i]>0?1:-1) * k / c;
+ coef_t coef2 = k / abs(SUBs[e2].coef[i]);
+ for (int j = nVars; j >= 0; j--)
+ SUBs[e2].coef[j] = SUBs[e2].coef[j] * coef2 - EQs[e].coef[j] * coef1;
+
+ coef_t g = 0;
+ for (int j = nVars; j >= 0; j--) {
+ g = gcd(abs(SUBs[e2].coef[j]), g);
+ if (g == 1)
+ break;
+ }
+ if (g != 0 && g != 1)
+ for (int j = nVars; j >= 0; j--)
+ SUBs[e2].coef[j] /= g;
+ }
+
+ // remove redundent wildcard equality
+ if (!preserveThisConstraint) {
+ if (e < nEQs-1)
+ eqnncpy (&EQs[e], &EQs[nEQs-1], nVars);
+ nEQs--;
+ deleteVariable(i);
+ }
+
+ break;
+ }
+ }
+
+ // remove multi-wildcard equality in approximation mode
+ if (inApproximateMode)
+ for (int e = nEQs-1; e >= 0; e--)
+ for (int i = nVars; i >= safeVars+1; i--)
+ if (EQs[e].coef[i] != 0) {
+ int j = i-1;
+ for (; j >= safeVars+1; j--)
+ if (EQs[e].coef[j] != 0)
+ break;
+
+ if (j != safeVars) {
+ if (e < nEQs-1)
+ eqnncpy (&EQs[e], &EQs[nEQs-1], nVars);
+ nEQs--;
+ }
+
+ break;
+ }
+
+ if (renormalize)
+ normalize();
+}
+
+
+void Problem:: check() const {
+#ifndef NDEBUG
+ int v = nSUBs;
+ checkVars(nVars+1);
+ for(int i = 1; i <= safeVars; i++) if (var[i] > 0) v++;
+ assert(v == varsOfInterest);
+ for(int e = 0; e < nGEQs; e++) assert(GEQs[e].touched || GEQs[e].key != 0);
+ if(!mayBeRed) {
+ for(int e = 0; e < nEQs; e++) assert(!EQs[e].color);
+ for(int e = 0; e < nGEQs; e++) assert(!GEQs[e].color);
+ }
+ else
+ for(int i = safeVars+1; i <= nVars; i++) {
+ int isBlack = 0;
+ int isRed = 0;
+ for(int e = 0; e < nEQs; e++)
+ if (EQs[e].coef[i]) {
+ if (EQs[e].color) isRed = 1;
+ else isBlack = 1;
+ }
+ for(int e = 0; e < nGEQs; e++)
+ if (GEQs[e].coef[i]) {
+ if (GEQs[e].color) isRed = 1;
+ else isBlack = 1;
+ }
+ if (isBlack && isRed && 0) {
+ fprintf(outputFile,"Mixed Red and Black variable:\n");
+ printProblem();
+ }
+ }
+#endif
+}
+
+
+void Problem::rememberRedConstraint(eqn *e, redType type, coef_t stride) {
+ // Check if this is really a stride constraint
+ if (type == redEQ && newVar == nVars && e->coef[newVar]) {
+ type = redStride;
+ stride = e->coef[newVar];
+ }
+ // else for(int i = safeVars+1; i <= nVars; i++) assert(!e->coef[i]); // outdated -- by chun 10/30/2008
+
+ assert(type != notRed);
+ assert(type == redStride || stride == 0);
+
+ if (TRACE) {
+ fprintf(outputFile,"being asked to remember red constraint:\n");
+ switch(type) {
+ case notRed: fprintf(outputFile,"notRed: ");
+ break;
+ case redGEQ: fprintf(outputFile,"Red: 0 <= ");
+ break;
+ case redLEQ: fprintf(outputFile,"Red: 0 >= ");
+ break;
+ case redEQ: fprintf(outputFile,"Red: 0 == ");
+ break;
+ case redStride: fprintf(outputFile,"Red stride " coef_fmt ": ",stride);
+ break;
+ }
+ printTerm(e,1);
+ fprintf(outputFile,"\n");
+ printProblem();
+ fprintf(outputFile,"----\n");
+ }
+
+ // Convert redLEQ to redGEQ
+ eqn mem;
+ eqnncpy(&mem,e, nVars);
+ e = &mem;
+ if (type == redLEQ) {
+ for(int i = 0; i <= safeVars; i++)
+ e->coef[i] = -e->coef[i];
+ type = redGEQ;
+ }
+
+ // Prepare coefficient array for red constraint
+ bool has_wildcard = false;
+ coef_t coef[varsOfInterest-nextWildcard+1];
+ for (int i = 0; i <= varsOfInterest-nextWildcard; i++)
+ coef[i] = 0;
+ for (int i = 0; i <= safeVars; i++) {
+ if (var[i] < 0) {
+ if (e->coef[i] != 0) {
+ coef[varsOfInterest-var[i]] = e->coef[i];
+ has_wildcard = true;
+ }
+ }
+ else
+ coef[var[i]] = e->coef[i];
+ }
+
+ // Sophisticated back substituion for wildcards, use Gaussian elimination
+ // as a fallback if no simple equations available. -- by chun 10/13/2008
+ if (has_wildcard) {
+ // Find substitutions involving wildcard
+ coef_t *repl_subs[nSUBs];
+ int num_wild_in_repl_subs[nSUBs];
+ int num_repl_subs = 0;
+ for (int i = 0; i < nSUBs; i++) {
+ int t = 0;
+ for (int j = 1; j <= safeVars; j++) {
+ if (var[j] < 0 && SUBs[i].coef[j] != 0)
+ t++;
+ }
+ if (t > 0) {
+ repl_subs[num_repl_subs] = new coef_t[varsOfInterest-nextWildcard+1];
+ for (int j = 0; j <= varsOfInterest-nextWildcard; j++)
+ repl_subs[num_repl_subs][j] = 0;
+
+ for (int k = 0; k <= safeVars; k++)
+ repl_subs[num_repl_subs][(var[k]<0)?varsOfInterest-var[k]:var[k]] = SUBs[i].coef[k];
+ repl_subs[num_repl_subs][SUBs[i].key] = -1;
+ num_wild_in_repl_subs[num_repl_subs] = t;
+ num_repl_subs++;
+ }
+ }
+
+ int wild_solved[-nextWildcard+1];
+ bool has_unsolved = false;
+ for (int i = 1; i <= -nextWildcard; i++) {
+ int minimum_wild = 0;
+ int pos;
+ for (int j = 0; j < num_repl_subs; j++)
+ if (repl_subs[j][varsOfInterest+i] != 0 && (minimum_wild == 0 || num_wild_in_repl_subs[j] < minimum_wild)) {
+ minimum_wild = num_wild_in_repl_subs[j];
+ pos = j;
+ }
+
+ if (minimum_wild == 0) {
+ wild_solved[i] = -1;
+ if (coef[varsOfInterest+i] != 0) {
+ fprintf(outputFile,"No feasible back substitutions available\n");
+ printProblem();
+ exit(1);
+ }
+ }
+ else if (minimum_wild == 1)
+ wild_solved[i] = pos;
+ else {
+ wild_solved[i] = -1;
+ if (coef[varsOfInterest+i] != 0)
+ has_unsolved = true;
+ }
+ }
+
+ // Gaussian elimination
+ while (has_unsolved) {
+ for (int i = 0; i < num_repl_subs; i++)
+ if (num_wild_in_repl_subs[i] > 1) {
+ for (int j = 1; j <= -nextWildcard; j++) {
+ if (repl_subs[i][varsOfInterest+j] != 0 && wild_solved[j] >= 0) {
+ int s = wild_solved[j];
+ coef_t l = lcm(abs(repl_subs[i][varsOfInterest+j]), abs(repl_subs[s][varsOfInterest+j]));
+ coef_t scale_1 = l/abs(repl_subs[i][varsOfInterest+j]);
+ coef_t scale_2 = l/abs(repl_subs[s][varsOfInterest+j]);
+ int sign = ((repl_subs[i][varsOfInterest+j]>0)?1:-1) * ((repl_subs[s][varsOfInterest+j]>0)?1:-1);
+ for (int k = 0; k <= varsOfInterest-nextWildcard; k++)
+ repl_subs[i][k] = scale_1*repl_subs[i][k] - sign*scale_2*repl_subs[s][k];
+ num_wild_in_repl_subs[i]--;
+ }
+ }
+
+ if (num_wild_in_repl_subs[i] == 1) {
+ for (int j = 1; j <= -nextWildcard; j++)
+ if (repl_subs[i][varsOfInterest+j] != 0) {
+ assert(wild_solved[j]==-1);
+ wild_solved[j] = i;
+ break;
+ }
+ }
+ else if (num_wild_in_repl_subs[i] > 1) {
+ int pos = 0;
+ for (int j = 1; j <= -nextWildcard; j++)
+ if (repl_subs[i][varsOfInterest+j] != 0) {
+ pos = j;
+ break;
+ }
+ assert(pos > 0);
+
+ for (int j = i+1; j < num_repl_subs; j++)
+ if (repl_subs[j][varsOfInterest+pos] != 0) {
+ coef_t l = lcm(abs(repl_subs[i][varsOfInterest+pos]), abs(repl_subs[j][varsOfInterest+pos]));
+ coef_t scale_1 = l/abs(repl_subs[i][varsOfInterest+pos]);
+ coef_t scale_2 = l/abs(repl_subs[j][varsOfInterest+pos]);
+ int sign = ((repl_subs[i][varsOfInterest+pos]>0)?1:-1) * ((repl_subs[j][varsOfInterest+pos]>0)?1:-1);
+ for (int k = 0; k <= varsOfInterest-nextWildcard; k++)
+ repl_subs[j][k] = scale_2*repl_subs[j][k] - sign*scale_1*repl_subs[i][k];
+
+ num_wild_in_repl_subs[j] = 0;
+ int first_wild = 0;
+ for (int k = 1; k <= -nextWildcard; k++)
+ if (repl_subs[j][varsOfInterest+k] != 0) {
+ num_wild_in_repl_subs[j]++;
+ first_wild = k;
+ }
+
+ if (num_wild_in_repl_subs[j] == 1) {
+ if (wild_solved[first_wild] < 0)
+ wild_solved[first_wild] = j;
+ }
+ }
+ }
+ }
+
+ has_unsolved = false;
+ for (int i = 1; i <= -nextWildcard; i++)
+ if (coef[varsOfInterest+i] != 0 && wild_solved[i] < 0) {
+ has_unsolved = true;
+ break;
+ }
+ }
+
+ // Substitute all widecards in the red constraint
+ for (int i = 1; i <= -nextWildcard; i++) {
+ if (coef[varsOfInterest+i] != 0) {
+ int s = wild_solved[i];
+ assert(s >= 0);
+
+ coef_t l = lcm(abs(coef[varsOfInterest+i]), abs(repl_subs[s][varsOfInterest+i]));
+ coef_t scale_1 = l/abs(coef[varsOfInterest+i]);
+ coef_t scale_2 = l/abs(repl_subs[s][varsOfInterest+i]);
+ int sign = ((coef[varsOfInterest+i]>0)?1:-1) * ((repl_subs[s][varsOfInterest+i]>0)?1:-1);
+ for (int j = 0; j <= varsOfInterest-nextWildcard; j++)
+ coef[j] = scale_1*coef[j] - sign*scale_2*repl_subs[s][j];
+
+ if (scale_1 != 1)
+ stride *= scale_1;
+ }
+ }
+
+ for (int i = 0; i < num_repl_subs; i++)
+ delete []repl_subs[i];
+ }
+
+ // Ready to insert into redMemory
+ int m = nMemories++;
+ redMemory[m].length = 0;
+ redMemory[m].kind = type;
+ redMemory[m].constantTerm = coef[0];
+ for(int i = 1; i <= varsOfInterest; i++)
+ if (coef[i]) {
+ int j = redMemory[m].length++;
+ redMemory[m].coef[j] = coef[i];
+ redMemory[m].var[j] = i;
+ }
+ if (type == redStride) redMemory[m].stride = stride;
+ if (DBUG) {
+ fprintf(outputFile,"Red constraint remembered\n");
+ printProblem();
+ }
+}
+
+void Problem::recallRedMemories() {
+ if (nMemories) {
+ if (TRACE) {
+ fprintf(outputFile,"Recalling red memories\n");
+ printProblem();
+ }
+
+ eqn* e = 0;
+ for(int m = 0; m < nMemories; m++) {
+ switch(redMemory[m].kind) {
+ case redGEQ:
+ {
+ int temporary_eqn = newGEQ();
+ e = &GEQs[temporary_eqn];
+ eqnnzero(e, nVars);
+ e->touched = 1;
+ break;
+ }
+ case redEQ:
+ {
+ int temporary_eqn = newEQ();
+ e = &EQs[temporary_eqn];
+ eqnnzero(e, nVars);
+ break;
+ }
+ case redStride:
+ {
+ int temporary_eqn = newEQ();
+ e = &EQs[temporary_eqn];
+ eqnnzero(e, nVars);
+ int i = addNewUnprotectedWildcard();
+ e->coef[i] = -redMemory[m].stride;
+ break;
+ }
+ default:
+ assert(0);
+ }
+ e->color = EQ_RED;
+ e->coef[0] = redMemory[m].constantTerm;
+ for(int i = 0; i < redMemory[m].length; i++) {
+ int v = redMemory[m].var[i];
+ assert(var[forwardingAddress[v]] == v);
+ e->coef[forwardingAddress[v]] = redMemory[m].coef[i];
+ }
+ }
+
+ nMemories = 0;
+ if (TRACE) {
+ fprintf(outputFile,"Red memories recalled\n");
+ printProblem();
+ }
+ }
+}
+
+void Problem::swapVars(int i, int j) {
+ if (DEBUG) {
+ use_ugly_names++;
+ fprintf(outputFile, "Swapping %d and %d\n", i, j);
+ printProblem();
+ use_ugly_names--;
+ }
+ std::swap(var[i], var[j]);
+ for (int e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[i] != GEQs[e].coef[j]) {
+ GEQs[e].touched = true;
+ coef_t t = GEQs[e].coef[i];
+ GEQs[e].coef[i] = GEQs[e].coef[j];
+ GEQs[e].coef[j] = t;
+ }
+ for (int e = nEQs - 1; e >= 0; e--)
+ if (EQs[e].coef[i] != EQs[e].coef[j]) {
+ coef_t t = EQs[e].coef[i];
+ EQs[e].coef[i] = EQs[e].coef[j];
+ EQs[e].coef[j] = t;
+ }
+ for (int e = nSUBs - 1; e >= 0; e--)
+ if (SUBs[e].coef[i] != SUBs[e].coef[j]) {
+ coef_t t = SUBs[e].coef[i];
+ SUBs[e].coef[i] = SUBs[e].coef[j];
+ SUBs[e].coef[j] = t;
+ }
+ if (DEBUG) {
+ use_ugly_names++;
+ fprintf(outputFile, "Swapping complete \n");
+ printProblem();
+ fprintf(outputFile, "\n");
+ use_ugly_names--;
+ }
+}
+
+void Problem::addingEqualityConstraint(int e) {
+ if (addingOuterEqualities && originalProblem != noProblem &&
+ originalProblem != this && !conservative) {
+ int e2 = originalProblem->newEQ();
+ if (TRACE)
+ fprintf(outputFile, "adding equality constraint %d to outer problem\n", e2);
+ eqnnzero(&originalProblem->EQs[e2], originalProblem->nVars);
+ for (int i = nVars; i >= 1; i--) {
+ int j;
+ for (j = originalProblem->nVars; j >= 1; j--)
+ if (originalProblem->var[j] == var[i])
+ break;
+ if (j <= 0 || (outerColor && j > originalProblem->safeVars)) {
+ if (DBUG)
+ fprintf(outputFile, "retracting\n");
+ originalProblem->nEQs--;
+ return;
+ }
+ originalProblem->EQs[e2].coef[j] = EQs[e].coef[i];
+ }
+ originalProblem->EQs[e2].coef[0] = EQs[e].coef[0];
+
+ originalProblem->EQs[e2].color = outerColor;
+ if (DBUG)
+ originalProblem->printProblem();
+ }
+}
+
+
+// Initialize hash codes for inequalities, remove obvious redundancy.
+// Case 1:
+// a1*x1+a2*x2+...>=c (1)
+// a1*x2+a2*x2+...>=c' (2)
+// if c>=c' then (2) is redundant, and vice versa.
+//
+// case 2:
+// a1*x1+a2*x2+...>=c (1)
+// a1*x1+a2*x2+...<=c' (2)
+// if c=c' then add equality of (1) or (2),
+// if c>c' then no solution.
+//
+// Finally it calls extended normalize process which handles
+// wildcards in redundacy removal.
+normalizeReturnType Problem::normalize() {
+ int i, j;
+ bool coupledSubscripts = false;
+
+ check();
+
+ for (int e = 0; e < nGEQs; e++) {
+ if (!GEQs[e].touched) {
+ if (!singleVarGEQ(&GEQs[e]))
+ coupledSubscripts = true;
+ }
+ else { // normalize e
+ coef_t g;
+ int topVar;
+ int i0;
+ coef_t hashCode;
+
+ {
+ int *p = &packing[0];
+ for (int k = 1; k <= nVars; k++)
+ if (GEQs[e].coef[k]) {
+ *(p++) = k;
+ }
+ topVar = (p - &packing[0]) - 1;
+ }
+
+ if (topVar == -1) {
+ if (GEQs[e].coef[0] < 0) {
+ // e has no solution
+ return (normalize_false);
+ }
+ deleteGEQ(e);
+ e--;
+ continue;
+ }
+ else if (topVar == 0) {
+ int singleVar = packing[0];
+ g = GEQs[e].coef[singleVar];
+ if (g > 0) {
+ GEQs[e].coef[singleVar] = 1;
+ GEQs[e].key = singleVar;
+ }
+ else {
+ g = -g;
+ GEQs[e].coef[singleVar] = -1;
+ GEQs[e].key = -singleVar;
+ }
+ if (g > 1)
+ GEQs[e].coef[0] = int_div(GEQs[e].coef[0], g);
+ }
+ else {
+ coupledSubscripts = true;
+ i0 = topVar;
+ i = packing[i0--];
+ g = GEQs[e].coef[i];
+ hashCode = g * (i + 3);
+ if (g < 0)
+ g = -g;
+ for (; i0 >= 0; i0--) {
+ coef_t x;
+ i = packing[i0];
+ x = GEQs[e].coef[i];
+ hashCode = hashCode * keyMult * (i + 3) + x;
+ if (x < 0)
+ x = -x;
+ if (x == 1) {
+ g = 1;
+ i0--;
+ break;
+ }
+ else
+ g = gcd(x, g);
+ }
+ for (; i0 >= 0; i0--) {
+ coef_t x;
+ i = packing[i0];
+ x = GEQs[e].coef[i];
+ hashCode = hashCode * keyMult * (i + 3) + x;
+ }
+ if (g > 1) {
+ GEQs[e].coef[0] = int_div(GEQs[e].coef[0], g);
+ i0 = topVar;
+ i = packing[i0--];
+ GEQs[e].coef[i] = GEQs[e].coef[i] / g;
+ hashCode = GEQs[e].coef[i] * (i + 3);
+ for (; i0 >= 0; i0--) {
+ i = packing[i0];
+ GEQs[e].coef[i] = GEQs[e].coef[i] / g;
+ hashCode = hashCode * keyMult * (i + 3) + GEQs[e].coef[i];
+ }
+ }
+
+ {
+ coef_t g2 = abs(hashCode); // get e's hash code
+ j = static_cast<int>(g2 % static_cast<coef_t>(hashTableSize));
+ assert (g2 % (coef_t) hashTableSize == j);
+ while (1) {
+ eqn *proto = &(hashMaster[j]);
+ if (proto->touched == g2) {
+ if (proto->coef[0] == topVar) {
+ if (hashCode >= 0)
+ for (i0 = topVar; i0 >= 0; i0--) {
+ i = packing[i0];
+ if (GEQs[e].coef[i] != proto->coef[i])
+ break;
+ }
+ else
+ for (i0 = topVar; i0 >= 0; i0--) {
+ i = packing[i0];
+ if (GEQs[e].coef[i] != -proto->coef[i])
+ break;
+ }
+
+ if (i0 < 0) {
+ if (hashCode >= 0)
+ GEQs[e].key = proto->key;
+ else
+ GEQs[e].key = -proto->key;
+ break;
+ }
+ }
+ }
+ else if (proto->touched < 0) { //insert e into the empty entry in hash table
+ eqnnzero(proto, nVars);
+ if (hashCode >= 0)
+ for (i0 = topVar; i0 >= 0; i0--) {
+ i = packing[i0];
+ proto->coef[i] = GEQs[e].coef[i];
+ }
+ else
+ for (i0 = topVar; i0 >= 0; i0--) {
+ i = packing[i0];
+ proto->coef[i] = -GEQs[e].coef[i];
+ }
+ proto->coef[0] = topVar;
+ proto->touched = g2;
+ proto->key = nextKey++;
+
+ if (proto->key > maxKeys) {
+ fprintf(outputFile, "too many hash keys generated \n");
+ fflush(outputFile);
+ exit(2);
+ }
+ if (hashCode >= 0)
+ GEQs[e].key = proto->key;
+ else
+ GEQs[e].key = -proto->key;
+ break;
+ }
+ j = (j + 1) % hashTableSize;
+ }
+ }
+ }
+ }
+
+ GEQs[e].touched = false;
+
+ {
+ int eKey = GEQs[e].key;
+ int e2;
+ if (e > 0) {
+ e2 = fastLookup[maxKeys - eKey];
+ if (e2 >= 0 && e2 < e && GEQs[e2].key == -eKey) {
+ // confirm it is indeed a match -- by chun 10/29/2008
+ int k;
+ for (k = nVars; k >= 1; k--)
+ if (GEQs[e2].coef[k] != -GEQs[e].coef[k])
+ break;
+
+ if (k == 0) {
+ if (GEQs[e2].coef[0] < -GEQs[e].coef[0]) {
+ // there is no solution from e and e2
+ return (normalize_false);
+ }
+ else if (GEQs[e2].coef[0] == -GEQs[e].coef[0]) {
+ // reduce e and e2 to an equation
+ int neweq = newEQ();
+ eqnncpy(&EQs[neweq], &GEQs[e], nVars);
+ EQs[neweq].color = GEQs[e].color || GEQs[e2].color;
+ addingEqualityConstraint(neweq);
+ }
+ }
+ }
+
+ e2 = fastLookup[maxKeys + eKey];
+ if (e2 >= 0 && e2 < e && GEQs[e2].key == eKey) {
+ // confirm it is indeed a match -- by chun 10/29/2008
+ int k;
+ for (k = nVars; k >= 1; k--)
+ if (GEQs[e2].coef[k] != GEQs[e].coef[k])
+ break;
+
+ if (k == 0) {
+ if (GEQs[e2].coef[0] > GEQs[e].coef[0] ||
+ (GEQs[e2].coef[0] == GEQs[e].coef[0] && GEQs[e2].color)) {
+ // e2 is redundant
+ GEQs[e2].coef[0] = GEQs[e].coef[0];
+ GEQs[e2].color = GEQs[e].color;
+ deleteGEQ(e);
+ e--;
+ continue;
+ }
+ else {
+ // e is redundant
+ deleteGEQ(e);
+ e--;
+ continue;
+ }
+ }
+ }
+ }
+ fastLookup[maxKeys + eKey] = e;
+ }
+ }
+
+ // bypass entended normalization for temporary problem
+ if (!isTemporary && !inApproximateMode)
+ normalize_ext();
+
+ return coupledSubscripts ? normalize_coupled : normalize_uncoupled;
+}
+
+//
+// Extended normalize process, remove redundancy involving wildcards.
+// e.g.
+// exists alpha, beta:
+// v1+8*alpha<=v2<=15+8*alpha (1)
+// v1+8*beta<=v2<=15+8*beta (2)
+// if there are no other inequalities involving alpha or beta,
+// then either (1) or (2) is redundant. Such case can't be simplified
+// by fourier-motzkin algorithm due to special meanings of existentials.
+//
+void Problem::normalize_ext() {
+ std::vector<BoolSet<> > disjoint_wildcards(nVars-safeVars, BoolSet<>(nVars-safeVars));
+ std::vector<BoolSet<> > wildcards_in_inequality(nVars-safeVars, BoolSet<>(nGEQs));
+ for (int i = 0; i < nVars-safeVars; i++) {
+ disjoint_wildcards[i].set(i);
+ }
+
+ // create disjoint wildcard sets according to inequalities
+ for (int e = 0; e < nGEQs; e++) {
+ std::vector<BoolSet<> >::iterator first_set = disjoint_wildcards.end();
+ for (int i = 0; i < nVars-safeVars; i++)
+ if (GEQs[e].coef[i+safeVars+1] != 0) {
+ wildcards_in_inequality[i].set(e);
+
+ std::vector<BoolSet<> >::iterator cur_set = disjoint_wildcards.end();
+ for (std::vector<BoolSet<> >::iterator j = disjoint_wildcards.begin(); j != disjoint_wildcards.end(); j++)
+ if ((*j).get(i)) {
+ cur_set = j;
+ break;
+ }
+ assert(cur_set!=disjoint_wildcards.end());
+ if (first_set == disjoint_wildcards.end())
+ first_set = cur_set;
+ else if (first_set != cur_set) {
+ *first_set |= *cur_set;
+ disjoint_wildcards.erase(cur_set);
+ }
+ }
+ }
+
+ // do not consider wildcards appearing in equalities
+ for (int e = 0; e < nEQs; e++)
+ for (int i = 0; i < nVars-safeVars; i++)
+ if (EQs[e].coef[i+safeVars+1] != 0) {
+ for (std::vector<BoolSet<> >::iterator j = disjoint_wildcards.begin(); j != disjoint_wildcards.end(); j++)
+ if ((*j).get(i)) {
+ disjoint_wildcards.erase(j);
+ break;
+ }
+ }
+
+ // create disjoint inequality sets
+ std::vector<BoolSet<> > disjoint_inequalities(disjoint_wildcards.size());
+ for (size_t i = 0; i < disjoint_wildcards.size(); i++)
+ for (int j = 0; j < nVars-safeVars; j++)
+ if (disjoint_wildcards[i].get(j))
+ disjoint_inequalities[i] |= wildcards_in_inequality[j];
+
+ // hash the inequality again, this time separate wildcard variables from
+ // regular variables
+ coef_t hash_safe[nGEQs];
+ coef_t hash_wild[nGEQs];
+ for (int e = 0; e < nGEQs; e++) {
+ coef_t hashCode = 0;
+ for (int i = 1; i <= safeVars; i++)
+ if (GEQs[e].coef[i] != 0)
+ hashCode = hashCode * keyMult * (i+3) + GEQs[e].coef[i];
+ hash_safe[e] = hashCode;
+
+ hashCode = 0;
+ for (int i = safeVars+1; i <= nVars; i++)
+ if (GEQs[e].coef[i] != 0)
+ hashCode = hashCode * keyMult + GEQs[e].coef[i];
+ hash_wild[e] = hashCode;
+ }
+
+ // sort hash keys for each disjoint set
+ std::vector<std::vector<std::pair<int, std::pair<coef_t, coef_t> > > > disjoint_hash(disjoint_inequalities.size());
+ for (size_t i = 0; i < disjoint_inequalities.size(); i++)
+ for (int e = 0; e < nGEQs; e++)
+ if (disjoint_inequalities[i].get(e)) {
+ std::vector<std::pair<int, std::pair<coef_t, coef_t> > >::iterator j = disjoint_hash[i].begin();
+ for (; j != disjoint_hash[i].end(); j++)
+ if ((hash_safe[e] > (*j).second.first) ||
+ (hash_safe[e] == (*j).second.first && hash_wild[e] > (*j).second.second))
+ break;
+ disjoint_hash[i].insert(j, std::make_pair(e, std::make_pair(hash_safe[e], hash_wild[e])));
+ }
+
+ // test wildcard equivalance
+ std::vector<bool> is_dead(nGEQs, false);
+ for (size_t i = 0; i < disjoint_wildcards.size(); i++) {
+ if (disjoint_inequalities[i].num_elem() == 0)
+ continue;
+
+ for (size_t j = i+1; j < disjoint_wildcards.size(); j++) {
+ if (disjoint_wildcards[i].num_elem() != disjoint_wildcards[j].num_elem() ||
+ disjoint_hash[i].size() != disjoint_hash[j].size())
+ continue;
+
+ bool match = true;
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ if (disjoint_hash[i][k].second != disjoint_hash[j][k].second) {
+ match = false;
+ break;
+ }
+ }
+ if (!match)
+ continue;
+
+ // confirm same coefficients for regular variables
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ for (int p = 1; p <= safeVars; p++)
+ if (GEQs[disjoint_hash[i][k].first].coef[p] != GEQs[disjoint_hash[j][k].first].coef[p]) {
+ match = false;
+ break;
+ }
+ if (!match)
+ break;
+ }
+ if (!match)
+ continue;
+
+ // now try combinatory wildcard matching
+ std::vector<int> wild_map(nVars-safeVars, -1);
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ int e1 = disjoint_hash[i][k].first;
+ int e2 = disjoint_hash[j][k].first;
+
+ for (int p = 0; p < nVars-safeVars; p++)
+ if (GEQs[e1].coef[p+safeVars+1] != 0) {
+ if (wild_map[p] == -1) {
+ for (int q = 0; q < nVars-safeVars; q++)
+ if (wild_map[q] == -1 &&
+ GEQs[e2].coef[q+safeVars+1] == GEQs[e1].coef[p+safeVars+1]) {
+ wild_map[p] = q;
+ wild_map[q] = p;
+ break;
+ }
+ if (wild_map[p] == -1) {
+ match = false;
+ break;
+ }
+ }
+ else if (GEQs[e2].coef[wild_map[p]+safeVars+1] != GEQs[e1].coef[p+safeVars+1]) {
+ match = false;
+ break;
+ }
+ }
+ if (!match)
+ break;
+
+ for (int p = 0; p < nVars-safeVars; p++)
+ if (GEQs[e2].coef[p+safeVars+1] != 0 &&
+ (wild_map[p] == -1 || GEQs[e2].coef[p+safeVars+1] != GEQs[e1].coef[wild_map[p]+safeVars+1])) {
+ match = false;
+ break;
+ }
+ if (!match)
+ break;
+ }
+ if (!match)
+ continue;
+
+ // check constants
+ int dir = 0;
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ if (GEQs[disjoint_hash[i][k].first].coef[0] > GEQs[disjoint_hash[j][k].first].coef[0]) {
+ if (dir == 0)
+ dir = 1;
+ else if (dir == -1) {
+ match = false;
+ break;
+ }
+ }
+ else if (GEQs[disjoint_hash[i][k].first].coef[0] < GEQs[disjoint_hash[j][k].first].coef[0]) {
+ if (dir == 0)
+ dir = -1;
+ else if (dir == 1) {
+ match = false;
+ break;
+ }
+ }
+ }
+ if (!match)
+ continue;
+
+ // check redness
+ int red_dir = 0;
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ if (GEQs[disjoint_hash[i][k].first].color > GEQs[disjoint_hash[j][k].first].color) {
+ if (red_dir == 0)
+ red_dir = 1;
+ else if (red_dir == -1) {
+ match = false;
+ break;
+ }
+ }
+ else if (GEQs[disjoint_hash[i][k].first].color < GEQs[disjoint_hash[j][k].first].color) {
+ if (red_dir == 0)
+ red_dir = -1;
+ else if (red_dir == 1) {
+ match = false;
+ break;
+ }
+ }
+ }
+ if (!match)
+ continue;
+
+ // remove redundant inequalities
+ if (dir == 1 || (dir == 0 && red_dir == 1)) {
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ GEQs[disjoint_hash[i][k].first].coef[0] = GEQs[disjoint_hash[j][k].first].coef[0];
+ GEQs[disjoint_hash[i][k].first].color = GEQs[disjoint_hash[j][k].first].color;
+ is_dead[disjoint_hash[j][k].first] = true;
+ }
+ }
+ else {
+ for (size_t k = 0; k < disjoint_hash[i].size(); k++) {
+ is_dead[disjoint_hash[j][k].first] = true;
+ }
+ }
+ }
+ }
+
+ // eliminate dead inequalities
+ for (int e = nGEQs-1; e >= 0; e--)
+ if (is_dead[e]) {
+ deleteGEQ(e);
+ }
+}
+
+
+void initializeOmega(void) {
+ if (omegaInitialized)
+ return;
+
+// assert(sizeof(eqn)==sizeof(int)*(headerWords)+sizeof(coef_t)*(1+maxVars));
+ nextWildcard = 0;
+ nextKey = maxVars + 1;
+ for (int i = 0; i < hashTableSize; i++)
+ hashMaster[i].touched = -1;
+
+ sprintf(wildName[1], "__alpha");
+ sprintf(wildName[2], "__beta");
+ sprintf(wildName[3], "__gamma");
+ sprintf(wildName[4], "__delta");
+ sprintf(wildName[5], "__tau");
+ sprintf(wildName[6], "__sigma");
+ sprintf(wildName[7], "__chi");
+ sprintf(wildName[8], "__omega");
+ sprintf(wildName[9], "__pi");
+ sprintf(wildName[10], "__ni");
+ sprintf(wildName[11], "__Alpha");
+ sprintf(wildName[12], "__Beta");
+ sprintf(wildName[13], "__Gamma");
+ sprintf(wildName[14], "__Delta");
+ sprintf(wildName[15], "__Tau");
+ sprintf(wildName[16], "__Sigma");
+ sprintf(wildName[17], "__Chi");
+ sprintf(wildName[18], "__Omega");
+ sprintf(wildName[19], "__Pi");
+
+ omegaInitialized = 1;
+}
+
+//
+// This is experimental (I would say, clinical) fact:
+// If the code below is removed then simplifyProblem cycles.
+//
+class brainDammage {
+public:
+ brainDammage();
+};
+
+brainDammage::brainDammage() {
+ initializeOmega();
+}
+
+static brainDammage Podgorny;
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_solve.cc b/omegalib/omega/src/omega_core/oc_solve.cc
new file mode 100644
index 0000000..c25b6d0
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_solve.cc
@@ -0,0 +1,1378 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+ Solve ineqalities.
+
+ Notes:
+
+ History:
+*****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+
+namespace omega {
+
+static int solveDepth = 0;
+#define maxDead maxmaxGEQs
+
+int Problem::solve(int desiredResult) {
+ assert(omegaInitialized);
+ int result;
+
+ checkVars(nVars+1);
+ assert(nVars >= safeVars);
+ if (desiredResult != OC_SOLVE_SIMPLIFY)
+ safeVars = 0;
+
+ solveDepth++;
+ if (solveDepth > 50) {
+ fprintf(outputFile, "Solve depth = %d, inApprox = %d, aborting\n", solveDepth, inApproximateMode);
+ printProblem();
+ fflush(outputFile);
+
+ if (solveDepth > 60)
+ exit(2);
+ }
+
+ check();
+ do {
+ doItAgain = 0;
+ check();
+ if (solveEQ() == false) {
+ solveDepth--;
+ return (false);
+ }
+ check();
+ if (!nGEQs) {
+ result = true;
+ nVars = safeVars;
+ break;
+ }
+ else
+ result = solveGEQ(desiredResult);
+ check();
+ }
+ while (doItAgain && desiredResult == OC_SOLVE_SIMPLIFY);
+ solveDepth--;
+
+ return (result);
+}
+
+
+// Supporting functions of solveGEQ
+int Problem::smoothWeirdEquations() {
+ int e1, e2, e3, p, q, k;
+ coef_t alpha, alpha1, alpha2, alpha3;
+ coef_t c;
+ int v;
+ int result = 0;
+
+ for (e1 = nGEQs - 1; e1 >= 0; e1--)
+ if (!GEQs[e1].color) {
+ coef_t g = 999999;
+ for (v = nVars; v >= 1; v--)
+ if (GEQs[e1].coef[v] != 0 && abs(GEQs[e1].coef[v]) < g)
+ g = abs(GEQs[e1].coef[v]);
+ if (g > 20) {
+ e3 = newGEQ(); /* Create a scratch GEQ,not part of the prob.*/
+ nGEQs--;
+ for (v = nVars; v >= 1; v--)
+ GEQs[e3].coef[v] = int_div(6 * GEQs[e1].coef[v] + g / 2, g);
+ GEQs[e3].color = EQ_BLACK;
+ GEQs[e3].touched = 1;
+ GEQs[e3].coef[0] = 9997;
+ if (DBUG) {
+ fprintf(outputFile, "Checking to see if we can derive: ");
+ printGEQ(&GEQs[e3]);
+ fprintf(outputFile, "\n from: ");
+ printGEQ(&GEQs[e1]);
+ fprintf(outputFile, "\n");
+ }
+
+
+ for (e2 = nGEQs - 1; e2 >= 0; e2--)
+ if (e1 != e2 && !GEQs[e2].color) {
+ for (p = nVars; p > 1; p--) {
+ for (q = p - 1; q > 0; q--) {
+ alpha = check_mul(GEQs[e1].coef[p], GEQs[e2].coef[q]) - check_mul(GEQs[e2].coef[p], GEQs[e1].coef[q]);
+ if (alpha != 0)
+ goto foundPQ;
+ }
+ }
+ continue;
+
+ foundPQ:
+
+ alpha1 = check_mul(GEQs[e2].coef[q], GEQs[e3].coef[p]) - check_mul(GEQs[e2].coef[p], GEQs[e3].coef[q]);
+ alpha2 = -(check_mul(GEQs[e1].coef[q], GEQs[e3].coef[p]) - check_mul(GEQs[e1].coef[p], GEQs[e3].coef[q]));
+ alpha3 = alpha;
+
+ if (alpha1 * alpha2 <= 0)
+ continue;
+ if (alpha1 < 0) {
+ alpha1 = -alpha1;
+ alpha2 = -alpha2;
+ alpha3 = -alpha3;
+ }
+ if (alpha3 > 0) {
+ /* Trying to prove e3 is redundant */
+
+ /* verify alpha1*v1+alpha2*v2 = alpha3*v3 */
+ for (k = nVars; k >= 1; k--)
+ if (check_mul(alpha3, GEQs[e3].coef[k])
+ != check_mul(alpha1, GEQs[e1].coef[k]) + check_mul(alpha2, GEQs[e2].coef[k]))
+ goto nextE2;
+
+ c = check_mul(alpha1, GEQs[e1].coef[0]) + check_mul(alpha2, GEQs[e2].coef[0]);
+ if (c < check_mul(alpha3, (GEQs[e3].coef[0] + 1)))
+ GEQs[e3].coef[0] = int_div(c, alpha3);
+
+ }
+ nextE2:;
+ }
+ if (GEQs[e3].coef[0] < 9997) {
+ result++;
+#if !defined NDEBUG
+ int e4 =
+#endif
+ newGEQ();
+#if !defined NDEBUG
+ assert(e3 == e4);
+#endif
+ if (DBUG) {
+ fprintf(outputFile, "Smoothing wierd equations; adding:\n");
+ printGEQ(&GEQs[e3]);
+ fprintf(outputFile, "\nto:\n");
+ printProblem();
+ fprintf(outputFile, "\n\n");
+ }
+ }
+ }
+ }
+ return (result);
+}
+
+
+void Problem::analyzeElimination(
+ int &v,
+ int &darkConstraints,
+ int &darkShadowFeasible,
+ int &unit,
+ coef_t &parallelSplinters,
+ coef_t &disjointSplinters,
+ coef_t &lbSplinters,
+ coef_t &ubSplinters,
+ int &parallelLB) {
+
+ parallelSplinters = (posInfinity); // was MAXINT
+ disjointSplinters = 0;
+ lbSplinters = 0;
+ ubSplinters = 0;
+
+ darkConstraints = 0;
+ darkShadowFeasible = 1;
+ coef_t maxUBc = 0;
+ coef_t maxLBc = 0;
+ int e,e2;
+ unit = 0;
+ int exact = 1;
+
+ for (e = nGEQs - 1; e >= 0; e--) {
+ coef_t c = GEQs[e].coef[v];
+
+ if (c < 0) {
+ coef_t Lc, Uc, g, diff, grey;
+
+ set_max(maxUBc, -c);
+ Uc = -c;
+ for (e2 = nGEQs - 1; e2 >= 0; e2--)
+ if (GEQs[e2].coef[v] > 0) {
+ Lc = GEQs[e2].coef[v];
+ g = 0;
+ grey = (Lc - 1) * (Uc - 1);
+
+ for (int j = nVars; j >= 1; j--) {
+ coef_t diff = check_mul(Lc, GEQs[e].coef[j]) + check_mul(Uc, GEQs[e2].coef[j]);
+ if (diff < 0) diff = -diff;
+ g = gcd(g, diff);
+ if (g == 1)
+ break;
+ }
+ diff = check_mul(Lc, GEQs[e].coef[0]) + check_mul(Uc, GEQs[e2].coef[0]);
+ if (g == 0) {
+ if (diff < 0) {
+ /* Real shadow must be true */
+ /* otherwise we would have found it during */
+ /* check for opposing constraints */
+ fprintf(outputFile, "Found conflicting constraints ");
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile," and ");
+ printGEQ(&GEQs[e2]);
+ fprintf(outputFile,"\nin\n");
+ printProblem();
+ assert(diff >= 0);
+ }
+ if (diff < grey) {
+ darkShadowFeasible = 0;
+ if (parallelSplinters > diff+1) {
+ parallelSplinters = diff + 1;
+ parallelLB = e2;
+ }
+ }
+ else {/* dark shadow is true, don't need to worry about this constraint pair */
+ }
+ }
+ else {
+ coef_t splinters= int_div(diff, g) - int_div(diff - grey, g);
+ if (splinters) exact = 0;
+ disjointSplinters += splinters;
+ if (g > 1) unit++;
+ darkConstraints++;
+ }
+ }
+ }
+ else if (c > 0) {
+ set_max(maxLBc, c);
+ } /* else
+ darkConstraints++; */
+ }
+
+ if (darkShadowFeasible) {
+ disjointSplinters++;
+ ubSplinters++;
+ lbSplinters++;
+ }
+ else disjointSplinters = (posInfinity); // was MAXINT
+
+
+ if (!darkShadowFeasible || !exact)
+ for (e = nGEQs - 1; e >= 0; e--) {
+ coef_t c = GEQs[e].coef[v];
+ if (c < -1) {
+ c = -c;
+ ubSplinters += 1+(check_mul(c, maxLBc) - c - maxLBc) / maxLBc;
+ }
+ else if (c > 1) {
+ lbSplinters += 1+ (check_mul(c, maxUBc) - c - maxUBc) / maxUBc;
+ }
+ }
+
+ if (DEBUG) {
+ fprintf(outputFile,"analyzing elimination of %s(%d)\n",variable(v),v);
+ if (darkShadowFeasible)
+ fprintf(outputFile," # dark constraints = %d\n", darkConstraints);
+ else
+ fprintf(outputFile," dark shadow obviously unfeasible\n");
+
+ fprintf(outputFile," " coef_fmt " LB splinters\n", lbSplinters);
+ fprintf(outputFile," " coef_fmt " UB splinters\n", ubSplinters);
+ if (disjointSplinters != (posInfinity))
+ fprintf(outputFile," " coef_fmt " disjoint splinters\n", disjointSplinters);
+ if (parallelSplinters != (posInfinity))
+ fprintf(outputFile," " coef_fmt " parallel splinters\n", parallelSplinters);
+ fprintf(outputFile, "\n");
+ fprintf(outputFile," %3d unit score \n", unit);
+ }
+}
+
+
+void Problem::partialElimination() {
+ if (DBUG) {
+ fprintf(outputFile, "Performing Partial elimination\n");
+ printProblem();
+ }
+ int fv;
+ if (0)
+ fv = 0;
+ else
+ fv = safeVars;
+ bool somethingHappened = false;
+ for (int i = nVars; i > fv; i--) {
+ bool isDead[maxmaxGEQs];
+ int e;
+ for (e = nGEQs-1; e >= 0; e--) isDead[e] = false;
+ int deadEqns[maxDead];
+ int numDead = 0;
+ for (int e1 = nGEQs-1; e1 >= 0; e1--)
+ if (abs(GEQs[e1].coef[i]) == 1) {
+ bool isGood = true;
+ for (int e2 = nGEQs-1; e2 >= 0; e2--)
+ if (check_mul(GEQs[e2].coef[i], GEQs[e1].coef[i]) < 0)
+ if (GEQs[e1].key != -GEQs[e2].key) {
+ coef_t Uc = abs(GEQs[e2].coef[i]);
+ for (int k = nVars; k > fv; k--)
+ if (GEQs[e2].coef[k] + check_mul(GEQs[e1].coef[k], Uc) != 0)
+ isGood = false;
+ }
+ if (isGood) {
+ somethingHappened = true;
+ for (int e2 = nGEQs-1; e2 >= 0; e2--)
+ if (check_mul(GEQs[e2].coef[i], GEQs[e1].coef[i]) < 0) {
+ if (GEQs[e1].key != -GEQs[e2].key) {
+ coef_t Uc = abs(GEQs[e2].coef[i]);
+ int new_eqn;
+ if (numDead == 0) {
+ new_eqn = newGEQ();
+ }
+ else {
+ new_eqn = deadEqns[--numDead];
+ }
+ isDead[new_eqn] = false;
+ if (DBUG) {
+ fprintf(outputFile,"Eliminating constraint on %s\n", variable(i));
+ fprintf(outputFile, "e1 = %d, e2 = %d, gen = %d\n", e1, e2, new_eqn);
+ printGEQ(&(GEQs[e1]));
+ fprintf(outputFile, "\n");
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n");
+ }
+
+ for (int k = nVars; k >= 0; k--)
+ GEQs[new_eqn].coef[k] = GEQs[e2].coef[k] + check_mul(GEQs[e1].coef[k], Uc);
+ GEQs[new_eqn].touched = true;
+ GEQs[new_eqn].color = GEQs[e2].color | GEQs[e1].color;
+ if (DBUG) {
+ fprintf(outputFile, "give ");
+ printGEQ(&(GEQs[new_eqn]));
+ fprintf(outputFile, "\n");
+ }
+ assert(GEQs[new_eqn].coef[i] == 0);
+ }
+ }
+ deadEqns[numDead++] = e1;
+ isDead[e1] = true;
+ if (DEBUG)
+ fprintf(outputFile, "Killed %d\n", e1);
+ }
+ }
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (isDead[e]) {
+ deleteGEQ(e);
+ }
+ }
+ if (somethingHappened && DBUG) {
+ fprintf(outputFile, "Result of Partial elimination\n");
+ printProblem();
+ }
+}
+
+
+int Problem:: solveGEQ(int desiredResult) {
+ int i, j, k, e;
+ int fv;
+ int result;
+ int coupledSubscripts;
+ int eliminateAgain;
+ int smoothed = 0;
+ int triedEliminatingRedundant = 0;
+ j = 0;
+
+ if (desiredResult != OC_SOLVE_SIMPLIFY) {
+ nSUBs = 0;
+ nMemories = 0;
+ safeVars = 0;
+ varsOfInterest = 0;
+ }
+
+solveGEQstart:
+ while (1) {
+ assert(desiredResult == OC_SOLVE_SIMPLIFY || nSUBs == 0);
+ check_number_GEQs(nGEQs);
+
+ if (DEBUG) {
+ fprintf(outputFile, "\nSolveGEQ(%d,%d):\n", desiredResult, pleaseNoEqualitiesInSimplifiedProblems);
+ printProblem();
+ fprintf(outputFile, "\n");
+ }
+
+#ifndef NDEBUG
+ for(e=0;e<nSUBs;e++)
+ for(i=safeVars+1;i<=nVars;i++)
+ assert(!SUBs[e].coef[i]);
+#endif
+
+ check();
+
+ if (nVars == 1) {
+ int uColor = EQ_BLACK;
+ int lColor = EQ_BLACK;
+ coef_t upperBound = posInfinity;
+ coef_t lowerBound = negInfinity;
+ for (e = nGEQs - 1; e >= 0; e--) {
+ coef_t a = GEQs[e].coef[1];
+ coef_t c = GEQs[e].coef[0];
+ /* our equation is ax + c >= 0, or ax >= -c, or c >= -ax */
+ if (a == 0) {
+ if (c < 0) {
+ if (TRACE)
+ fprintf(outputFile, "equations have no solution (G)\n");
+ return (false);
+ }
+ }
+ else if (a > 0) {
+ if (a != 1)
+ c = int_div(c, a);
+ if (lowerBound < -c || (lowerBound == -c && !isRed(&GEQs[e]))) {
+ lowerBound = -c;
+ lColor = GEQs[e].color;
+ }
+ }
+ else {
+ if (a != -1)
+ c = int_div(c, -a);
+ if (upperBound > c || (upperBound == c && !isRed(&GEQs[e]))) {
+ upperBound = c;
+ uColor = GEQs[e].color;
+ }
+ }
+ }
+ if (DEBUG)
+ fprintf(outputFile, "upper bound = " coef_fmt "\n", upperBound);
+ if (DEBUG)
+ fprintf(outputFile, "lower bound = " coef_fmt "\n", lowerBound);
+ if (lowerBound > upperBound) {
+ if (TRACE)
+ fprintf(outputFile, "equations have no solution (H)\n");
+ return (false);
+ }
+ if (desiredResult == OC_SOLVE_SIMPLIFY) {
+ nGEQs = 0;
+ if (safeVars == 1) {
+ if (lowerBound == upperBound && !uColor && !lColor) {
+ int e = newEQ();
+ assert(e == 0);
+ EQs[e].coef[0] = -lowerBound;
+ EQs[e].coef[1] = 1;
+ EQs[e].color = lColor | uColor;
+ return (solve(desiredResult));
+ }
+ else {
+ if (lowerBound > negInfinity) {
+ int e = newGEQ();
+ assert(e == 0);
+ GEQs[e].coef[0] = -lowerBound;
+ GEQs[e].coef[1] = 1;
+ GEQs[e].key = 1;
+ GEQs[e].color = lColor;
+ GEQs[e].touched = 0;
+ }
+ if (upperBound < posInfinity) {
+ int e = newGEQ();
+ GEQs[e].coef[0] = upperBound;
+ GEQs[e].coef[1] = -1;
+ GEQs[e].key = -1;
+ GEQs[e].color = uColor;
+ GEQs[e].touched = 0;
+ }
+ }
+ }
+ else
+ nVars = 0;
+ return (true);
+ }
+ if (originalProblem != noProblem && !lColor && !uColor && !conservative && lowerBound == upperBound) {
+ int e = newEQ();
+ assert(e == 0);
+ EQs[e].coef[0] = -lowerBound;
+ EQs[e].coef[1] = 1;
+ EQs[e].color = EQ_BLACK;
+ addingEqualityConstraint(0);
+ }
+ return (true);
+ }
+
+ if (!variablesFreed) {
+ variablesFreed = 1;
+ if (desiredResult != OC_SOLVE_SIMPLIFY)
+ freeEliminations(0);
+ else
+ freeEliminations(safeVars);
+ if (nVars == 1)
+ continue;
+ }
+
+
+ switch (normalize()) {
+ case normalize_false:
+ return (false);
+ break;
+ case normalize_coupled:
+ coupledSubscripts = true;
+ break;
+ case normalize_uncoupled:
+ coupledSubscripts = false;
+ break;
+ default:
+ coupledSubscripts = false;
+ assert(0 && "impossible case in SolveGEQ");
+ }
+
+
+ if ((doTrace && desiredResult == OC_SOLVE_SIMPLIFY) || DBUG) {
+ fprintf(outputFile, "\nafter normalization:\n");
+ printProblem();
+ fprintf(outputFile, "\n");
+ for(e=0;e<nGEQs;e++) assert(!GEQs[e].touched);
+ fprintf(outputFile, "eliminating variable using fourier-motzkin elimination\n");
+ }
+
+ // eliminating variable using fourier-motzkin elimination
+ do {
+ eliminateAgain = 0;
+
+ if (nEQs > 0)
+ return (solve(desiredResult));
+
+ if (!coupledSubscripts) {
+ if (safeVars == 0)
+ nGEQs = 0;
+ else
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].key > safeVars || -safeVars > GEQs[e].key)
+ deleteGEQ(e);
+ nVars = safeVars;
+ return (true);
+ }
+
+ if (desiredResult != OC_SOLVE_SIMPLIFY)
+ fv = 0;
+ else
+ fv = safeVars;
+
+ if (nVars == 0 || nGEQs == 0) {
+ nGEQs = 0;
+ if (desiredResult == OC_SOLVE_SIMPLIFY)
+ nVars = safeVars;
+ return (true);
+ }
+ if (desiredResult == OC_SOLVE_SIMPLIFY && nVars == safeVars) {
+ return (true);
+ }
+
+
+ if (nGEQs+6 > maxGEQs || nGEQs > 2 * nVars * nVars + 4 * nVars + 10) {
+ if (TRACE)
+ fprintf(outputFile, "TOO MANY EQUATIONS; %d equations, %d variables, ELIMINATING REDUNDANT ONES\n", nGEQs, nVars);
+ if (!quickKill(0,true))
+ return 0;
+ if (nEQs > 0)
+ return (solve(desiredResult));
+ if (TRACE)
+ fprintf(outputFile, "END ELIMINATION OF REDUNDANT EQUATIONS\n");
+ if (DBUG) printProblem();
+ }
+
+
+ {
+ int darkConstraints, darkShadowFeasible, unit, parallelLB;
+ coef_t parallelSplinters, disjointSplinters, lbSplinters, ubSplinters, splinters;
+ coef_t bestScore, score;
+ int bestVar;
+ int exact;
+ int Ue,Le;
+
+ if (desiredResult != OC_SOLVE_SIMPLIFY) fv = 0;
+ else fv = safeVars;
+
+ if (DEBUG) {
+ fprintf(outputFile,"Considering elimination possibilities[ \n");
+ printProblem();
+ }
+
+ analyzeGEQstart:
+ try {
+ bestScore = posInfinity;
+ bestVar = -1;
+ for (i = nVars; i != fv; i--) {
+ analyzeElimination(i, darkConstraints, darkShadowFeasible, unit, parallelSplinters, disjointSplinters, lbSplinters, ubSplinters, parallelLB);
+
+ score = min(min(parallelSplinters,disjointSplinters),
+ min(lbSplinters,ubSplinters));
+ exact = score == 1;
+ score = 10000*(score-1) + darkConstraints;
+ if (score >= posInfinity) // too big the score
+ score = posInfinity - 1;
+ score -= 3*unit;
+
+ if (score < bestScore) {
+ bestScore = score;
+ bestVar = i;
+ if (i > 4 && score < nGEQs) break;
+ }
+ }
+ assert(bestVar>=0);
+ exact = bestScore < 10000;
+ i = bestVar;
+ assert(i<=nVars);
+ analyzeElimination(i, darkConstraints, darkShadowFeasible, unit, parallelSplinters, disjointSplinters, lbSplinters, ubSplinters, parallelLB);
+ if (DEBUG) {
+ fprintf(outputFile,"] Choose to eliminate %s \n",variable(i));
+ }
+ splinters = lbSplinters;
+ if (splinters <= parallelSplinters)
+ parallelSplinters = posInfinity;
+ else splinters = parallelSplinters;
+ if (disjointSplinters == 1) splinters = 1;
+ exact = splinters == 1;
+ if (inApproximateMode) exact = 1;
+ }
+ catch (std::overflow_error) {
+ int result = quickKill(0, true);
+ if (result == 0)
+ return 0;
+ else if (result == 1)
+ return true;
+ else {
+ if (nEQs > 0)
+ return (solve(desiredResult));
+ triedEliminatingRedundant = 1;
+ goto analyzeGEQstart;
+ }
+ }
+
+ if (!triedEliminatingRedundant && darkConstraints > maxGEQs) {
+ if (TRACE)
+ fprintf(outputFile, "Elimination will create TOO MANY EQUATIONS; %d equations, %d variables, %d new constraints, ELIMINATING REDUNDANT ONES\n", nGEQs, nVars,darkConstraints);
+ if (!quickKill(0))
+ return 0;
+ if (nEQs > 0)
+ return (solve(desiredResult));
+ if (TRACE)
+ fprintf(outputFile, "END ELIMINATION OF REDUNDANT EQUATIONS\n");
+ if (DBUG) printProblem();
+
+ triedEliminatingRedundant = 1;
+ eliminateAgain = 1;
+ continue;
+ }
+
+ if (!exact && !triedEliminatingRedundant &&
+ safeVars > 0 && desiredResult == OC_SOLVE_SIMPLIFY) {
+ if (TRACE)
+ fprintf(outputFile, "Trying to produce exact elimination by finding redundant constraints [\n");
+ if (!quickKill(1)) return 0;
+ if (TRACE)
+ fprintf(outputFile, "]\n");
+ triedEliminatingRedundant = 1;
+ eliminateAgain = 1;
+ continue;
+ }
+ triedEliminatingRedundant = 0;
+
+ if (desiredResult == OC_SOLVE_SIMPLIFY && !exact) {
+ partialElimination();
+ switch (normalize()) {
+ case normalize_false:
+ return (false);
+ break;
+ case normalize_coupled:
+ case normalize_uncoupled:
+ break;
+ }
+ if (nEQs) return solveEQ();
+ if (DBUG) fprintf(outputFile,"Stopping short due to non-exact elimination\n");
+ return (true);
+ }
+
+ if ( desiredResult == OC_SOLVE_SIMPLIFY && darkConstraints > maxGEQs) {
+ if (DBUG) fprintf(outputFile,"Stopping short due to overflow of GEQs: %d\n", darkConstraints);
+ return (true);
+ }
+
+ if ((doTrace && desiredResult == OC_SOLVE_SIMPLIFY) || DBUG) {
+ fprintf(outputFile, "going to eliminate %s, (%d)\n", variable(i), i);
+ if (DEBUG)
+ printProblem();
+ fprintf(outputFile, "score = " coef_fmt "/" coef_fmt "\n", bestScore,splinters);
+ }
+
+ if (!exact && desiredResult == OC_SOLVE_SIMPLIFY && parallelSplinters == splinters) {
+ return parallelSplinter(parallelLB, parallelSplinters, desiredResult);
+ }
+
+ // smoothed = 0; // what a bug!!! -- by chun 6/10/2008
+
+ if (i != nVars) {
+ j = nVars;
+ swapVars(i,j);
+
+ i = j;
+ }
+ else if (DEBUG) {
+ printVars();
+ fprintf(outputFile, "No swap needed before eliminating %s(%d/%d)\n",variable(i),i,nVars);
+ for(j=1;j<=i;j++) fprintf(outputFile,"var #%d = %s(%x)\n",j,variable(j),var[j]);
+ printProblem();
+ }
+ nVars--;
+
+ if (exact) {
+ if (nVars == 1) {
+ coef_t upperBound = posInfinity;
+ coef_t lowerBound = negInfinity;
+ int ub_color = 0;
+ int lb_color = 0;
+ coef_t constantTerm, coefficient;
+ int topEqn = nGEQs - 1;
+ coef_t Lc;
+ for (Le = topEqn; Le >= 0; Le--)
+ if ((Lc = GEQs[Le].coef[i]) == 0) {
+ if (GEQs[Le].coef[1] == 1) {
+ constantTerm = -GEQs[Le].coef[0];
+ if (constantTerm > lowerBound || (constantTerm == lowerBound && !isRed(&GEQs[Le]))) {
+ lowerBound = constantTerm;
+ lb_color = GEQs[Le].color;
+ }
+ if (DEBUG) {
+ if (GEQs[Le].color == EQ_BLACK)
+ fprintf(outputFile, " :::=> %s >= " coef_fmt "\n", variable(1), constantTerm);
+ else
+ fprintf(outputFile, " :::=> [%s >= " coef_fmt "]\n", variable(1), constantTerm);
+ }
+ }
+ else {
+ constantTerm = GEQs[Le].coef[0];
+ if (constantTerm < upperBound || (constantTerm == upperBound && !isRed(&GEQs[Le]))) {
+ upperBound = constantTerm;
+ ub_color = GEQs[Le].color;
+ }
+ if (DEBUG) {
+ if (GEQs[Le].color == EQ_BLACK)
+ fprintf(outputFile, " :::=> %s <= " coef_fmt "\n", variable(1), constantTerm);
+ else
+ fprintf(outputFile, " :::=> [%s <= " coef_fmt "]\n", variable(1), constantTerm);
+ }
+ }
+ }
+ else if (Lc > 0) {
+ for (Ue = topEqn; Ue >= 0; Ue--)
+ if (GEQs[Ue].coef[i] < 0) {
+ if (GEQs[Le].key != -GEQs[Ue].key) {
+ coef_t Uc = -GEQs[Ue].coef[i];
+ coefficient = check_mul(GEQs[Ue].coef[1], Lc) + check_mul(GEQs[Le].coef[1], Uc);
+ constantTerm = check_mul(GEQs[Ue].coef[0], Lc) + check_mul(GEQs[Le].coef[0], Uc);
+ if (DEBUG) {
+ printGEQextra(&(GEQs[Ue]));
+ fprintf(outputFile, "\n");
+ printGEQextra(&(GEQs[Le]));
+ fprintf(outputFile, "\n");
+ }
+ if (coefficient > 0) {
+ constantTerm = -(int_div(constantTerm, coefficient));
+ /* assert(black == 0) */
+ if (constantTerm > lowerBound ||
+ (constantTerm == lowerBound &&
+ (desiredResult != OC_SOLVE_SIMPLIFY || (GEQs[Ue].color == EQ_BLACK && GEQs[Le].color == EQ_BLACK)))) {
+ lowerBound = constantTerm;
+ lb_color = GEQs[Ue].color || GEQs[Le].color;
+ }
+ if (DEBUG) {
+ if (GEQs[Ue].color || GEQs[Le].color)
+ fprintf(outputFile, " ::=> [%s >= " coef_fmt "]\n", variable(1), constantTerm);
+ else
+ fprintf(outputFile, " ::=> %s >= " coef_fmt "\n", variable(1), constantTerm);
+ }
+ }
+ else if (coefficient < 0) {
+ constantTerm = (int_div(constantTerm, -coefficient));
+ if (constantTerm < upperBound ||
+ (constantTerm == upperBound && GEQs[Ue].color == EQ_BLACK && GEQs[Le].color == EQ_BLACK)) {
+ upperBound = constantTerm;
+ ub_color = GEQs[Ue].color || GEQs[Le].color;
+ }
+ if (DEBUG) {
+ if (GEQs[Ue].color || GEQs[Le].color)
+ fprintf(outputFile, " ::=> [%s <= " coef_fmt "]\n", variable(1), constantTerm);
+ else
+ fprintf(outputFile, " ::=> %s <= " coef_fmt "\n", variable(1), constantTerm);
+ }
+ }
+ }
+ }
+ }
+ nGEQs = 0;
+ if (DEBUG)
+ fprintf(outputFile, " therefore, %c" coef_fmt " <= %c%s%c <= " coef_fmt "%c\n", lb_color ? '[' : ' ', lowerBound, (lb_color && !ub_color) ? ']' : ' ', variable(1), (!lb_color && ub_color) ? '[' : ' ', upperBound, ub_color ? ']' : ' ');
+ if (lowerBound > upperBound)
+ return (false);
+
+ if (upperBound == lowerBound) {
+ int e = newEQ();
+ assert(e == 0);
+ EQs[e].coef[1] = -1;
+ EQs[e].coef[0] = upperBound;
+ EQs[e].color = ub_color | lb_color;
+ addingEqualityConstraint(0);
+ }
+ else if (safeVars == 1) {
+ if (upperBound != posInfinity) {
+ int e = newGEQ();
+ assert(e == 0);
+ GEQs[e].coef[1] = -1;
+ GEQs[e].coef[0] = upperBound;
+ GEQs[e].color = ub_color;
+ GEQs[e].key = -1;
+ GEQs[e].touched = 0;
+ }
+ if (lowerBound != negInfinity) {
+ int e = newGEQ();
+ GEQs[e].coef[1] = 1;
+ GEQs[e].coef[0] = -lowerBound;
+ GEQs[e].color = lb_color;
+ GEQs[e].key = 1;
+ GEQs[e].touched = 0;
+ }
+ }
+ if (safeVars == 0)
+ nVars = 0;
+ return (true);
+ }
+ eliminateAgain = 1;
+
+ {
+ int deadEqns[maxDead];
+ int numDead = 0;
+ int topEqn = nGEQs - 1;
+ int lowerBoundCount = 0;
+ for (Le = topEqn; Le >= 0; Le--)
+ if (GEQs[Le].coef[i] > 0)
+ lowerBoundCount++;
+ if (DEBUG)
+ fprintf(outputFile, "lower bound count = %d\n", lowerBoundCount);
+ if (lowerBoundCount == 0) {
+ if (desiredResult != OC_SOLVE_SIMPLIFY) fv = 0;
+ else fv = safeVars;
+ nVars++;
+ freeEliminations(fv);
+ continue;
+ }
+ for (Le = topEqn; Le >= 0; Le--)
+ if (GEQs[Le].coef[i] > 0) {
+ coef_t Lc = GEQs[Le].coef[i];
+ for (Ue = topEqn; Ue >= 0; Ue--)
+ if (GEQs[Ue].coef[i] < 0) {
+ if (GEQs[Le].key != -GEQs[Ue].key) {
+ coef_t Uc = -GEQs[Ue].coef[i];
+ int e2;
+ if (numDead == 0) {
+ /*( Big kludge warning ) */
+ /* this code is still using location nVars+1 */
+ /* but newGEQ, if it reallocates, only copies*/
+ /* locations up to nVars. This fixes that. */
+ nVars++;
+ e2 = newGEQ();
+ nVars--;
+ }
+ else {
+ e2 = deadEqns[--numDead];
+ }
+ if (DEBUG)
+ fprintf(outputFile, "Le = %d, Ue = %d, gen = %d\n", Le, Ue, e2);
+ if (DEBUG) {
+ printGEQextra(&(GEQs[Le]));
+ fprintf(outputFile, "\n");
+ printGEQextra(&(GEQs[Ue]));
+ fprintf(outputFile, "\n");
+ }
+ eliminateAgain = 0;
+ coef_t g = gcd(Lc,Uc);
+ coef_t Lc_over_g = Lc/g;
+ coef_t Uc_over_g = Uc/g;
+
+ for (k = nVars; k >= 0; k--)
+ GEQs[e2].coef[k] =
+ check_mul(GEQs[Ue].coef[k], Lc_over_g) + check_mul(GEQs[Le].coef[k], Uc_over_g);
+
+ GEQs[e2].coef[nVars + 1] = 0;
+ GEQs[e2].touched = true;
+ GEQs[e2].color = GEQs[Ue].color | GEQs[Le].color;
+
+ if (DEBUG) {
+ printGEQ(&(GEQs[e2]));
+ fprintf(outputFile, "\n");
+ }
+ }
+ if (lowerBoundCount == 1) {
+ deadEqns[numDead++] = Ue;
+ if (DEBUG)
+ fprintf(outputFile, "Killed %d\n", Ue);
+ }
+ }
+ lowerBoundCount--;
+ deadEqns[numDead++] = Le;
+ if (DEBUG)
+ fprintf(outputFile, "Killed %d\n", Le);
+ }
+
+ {
+ int isDead[maxmaxGEQs];
+ for (e = nGEQs - 1; e >= 0; e--)
+ isDead[e] = false;
+ while (numDead > 0) {
+ e = deadEqns[--numDead];
+ isDead[e] = true;
+ }
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (isDead[e]) {
+ nVars++;
+ deleteGEQ(e);
+ nVars--;
+ }
+ }
+ continue;
+ }
+ }
+ else {
+ Problem *rS, *iS;
+
+ rS = new Problem;
+ iS = new Problem;
+
+ iS->nVars = rS->nVars = nVars; // do this immed.; in case of reallocation, we
+ // need to know how much to copy
+ rS->get_var_name = get_var_name;
+ rS->getVarNameArgs = getVarNameArgs;
+ iS->get_var_name = get_var_name;
+ iS->getVarNameArgs = getVarNameArgs;
+
+ for (e = 0; e < nGEQs; e++)
+ if (GEQs[e].coef[i] == 0) {
+ int re2 = rS->newGEQ();
+ int ie2 = iS->newGEQ();
+ eqnncpy(&(rS->GEQs[re2]), &GEQs[e], nVars);
+ eqnncpy(&(iS->GEQs[ie2]), &GEQs[e], nVars);
+ if (DEBUG) {
+ int t;
+ fprintf(outputFile, "Copying (%d, " coef_fmt "): ", i, GEQs[e].coef[i]);
+ printGEQextra(&GEQs[e]);
+ fprintf(outputFile, "\n");
+ for (t = 0; t <= nVars + 1; t++)
+ fprintf(outputFile, coef_fmt " ", GEQs[e].coef[t]);
+ fprintf(outputFile, "\n");
+ }
+ }
+ for (Le = nGEQs - 1; Le >= 0; Le--)
+ if (GEQs[Le].coef[i] > 0) {
+ coef_t Lc = GEQs[Le].coef[i];
+ for (Ue = nGEQs - 1; Ue >= 0; Ue--)
+ if (GEQs[Ue].coef[i] < 0)
+ if (GEQs[Le].key != -GEQs[Ue].key) {
+ coef_t Uc = -GEQs[Ue].coef[i];
+ coef_t g = gcd(Lc,Uc);
+ coef_t Lc_over_g = Lc/g;
+ coef_t Uc_over_g = Uc/g;
+ int re2 = rS->newGEQ();
+ int ie2 = iS->newGEQ();
+ rS->GEQs[re2].touched = iS->GEQs[ie2].touched = true;
+ if (DEBUG) {
+ fprintf(outputFile, "---\n");
+ fprintf(outputFile, "Le(Lc) = %d(" coef_fmt "), Ue(Uc) = %d(" coef_fmt "), gen = %d\n", Le, Lc, Ue, Uc, ie2);
+ printGEQextra(&GEQs[Le]);
+ fprintf(outputFile, "\n");
+ printGEQextra(&GEQs[Ue]);
+ fprintf(outputFile, "\n");
+ }
+
+ if (Uc == Lc) {
+ for (k = nVars; k >= 0; k--)
+ iS->GEQs[ie2].coef[k] = rS->GEQs[re2].coef[k] =
+ GEQs[Ue].coef[k] + GEQs[Le].coef[k];
+ iS->GEQs[ie2].coef[0] -= (Uc - 1);
+ }
+ else {
+ for (k = nVars; k >= 0; k--)
+ iS->GEQs[ie2].coef[k] = rS->GEQs[re2].coef[k] =
+ check_mul(GEQs[Ue].coef[k], Lc_over_g) + check_mul(GEQs[Le].coef[k], Uc_over_g);
+ iS->GEQs[ie2].coef[0] -= check_mul(Uc_over_g-1, Lc_over_g-1);
+ }
+
+ iS->GEQs[ie2].color = rS->GEQs[re2].color
+ = GEQs[Ue].color || GEQs[Le].color;
+
+ if (DEBUG) {
+ printGEQ(&(rS->GEQs[re2]));
+ fprintf(outputFile, "\n");
+ }
+ // ie2 = iS->newGEQ();
+ // re2 = rS->newGEQ();
+ }
+
+ }
+ iS->variablesInitialized = rS->variablesInitialized = 1;
+ iS->nEQs = rS->nEQs = 0;
+ assert(desiredResult != OC_SOLVE_SIMPLIFY);
+ assert(nSUBs == 0);
+ iS->nSUBs = rS->nSUBs = nSUBs;
+ iS->safeVars = rS->safeVars = safeVars;
+ int t;
+ for (t = nVars; t >= 0; t--)
+ rS->var[t] = var[t];
+ for (t = nVars; t >= 0; t--)
+ iS->var[t] = var[t];
+ nVars++;
+ if (desiredResult != true) {
+ int t = trace;
+ if (TRACE)
+ fprintf(outputFile, "\nreal solution(%d):\n", depth);
+ depth++;
+ trace = 0;
+ if (originalProblem == noProblem) {
+ originalProblem = this;
+ result = rS->solveGEQ(false);
+ originalProblem = noProblem;
+ }
+ else
+ result = rS->solveGEQ(false);
+ trace = t;
+ depth--;
+ if (result == false) {
+ delete rS;
+ delete iS;
+ return (result);
+ }
+
+ if (nEQs > 0) {
+ /* An equality constraint must have been found */
+ delete rS;
+ delete iS;
+ return (solve(desiredResult));
+ }
+ }
+ if (desiredResult != false) {
+ if (darkShadowFeasible) {
+ if (TRACE)
+ fprintf(outputFile, "\ninteger solution(%d):\n", depth);
+ depth++;
+ conservative++;
+ result = iS->solveGEQ(desiredResult);
+ conservative--;
+ depth--;
+ if (result != false) {
+ delete rS;
+ delete iS;
+ return (result);
+ }
+ }
+ if (TRACE)
+ fprintf(outputFile, "have to do exact analysis\n");
+
+ {
+ coef_t worstLowerBoundConstant=1;
+ int lowerBounds = 0;
+ int lowerBound[maxmaxGEQs];
+ int smallest;
+ int t;
+ conservative++;
+ for (e = 0; e < nGEQs; e++)
+ if (GEQs[e].coef[i] < -1) {
+ set_max(worstLowerBoundConstant,
+ -GEQs[e].coef[i]);
+ }
+ else if (GEQs[e].coef[i] > 1)
+ lowerBound[lowerBounds++] = e;
+ /* sort array */
+ for (j = 0; j < lowerBounds; j++) {
+ smallest = j;
+ for (k = j + 1; k < lowerBounds; k++)
+ if (GEQs[lowerBound[smallest]].coef[i] > GEQs[lowerBound[k]].coef[i])
+ smallest = k;
+ t = lowerBound[smallest];
+ lowerBound[smallest] = lowerBound[j];
+ lowerBound[j] = t;
+ }
+ if (DEBUG) {
+ fprintf(outputFile, "lower bound coeeficients = ");
+ for (j = 0; j < lowerBounds; j++)
+ fprintf(outputFile, " " coef_fmt, GEQs[lowerBound[j]].coef[i]);
+ fprintf(outputFile, "\n");
+ }
+
+
+ for (j = 0; j < lowerBounds; j++) {
+ coef_t maxIncr;
+ coef_t c;
+ e = lowerBound[j];
+ maxIncr = (check_mul(GEQs[e].coef[i]-1, worstLowerBoundConstant-1) - 1) / worstLowerBoundConstant;
+
+ /* maxIncr += 2; */
+ if ((doTrace && desiredResult == OC_SOLVE_SIMPLIFY) || DBUG) {
+ fprintf(outputFile, "for equation ");
+ printGEQ(&GEQs[e]);
+ fprintf(outputFile, "\ntry decrements from 0 to " coef_fmt "\n", maxIncr);
+ printProblem();
+ }
+ if (maxIncr > 50) {
+ if (!smoothed && smoothWeirdEquations()) {
+ conservative--;
+ delete rS;
+ delete iS;
+ smoothed = 1;
+ goto solveGEQstart;
+ }
+ }
+ int neweq = newEQ();
+ assert(neweq == 0);
+ eqnncpy(&EQs[neweq], &GEQs[e], nVars);
+ /*
+ * if (GEQs[e].color) fprintf(outputFile,"warning: adding black equality constraint
+ * based on red inequality\n");
+ */
+ EQs[neweq].color = EQ_BLACK;
+ eqnnzero(&GEQs[e], nVars);
+ GEQs[e].touched = true;
+ for (c = maxIncr; c >= 0; c--) {
+ if (DBUG)
+ fprintf(outputFile, "trying next decrement of " coef_fmt "\n", maxIncr - c);
+ if (DBUG)
+ printProblem();
+ *rS = *this;
+ if (DEBUG)
+ rS->printProblem();
+ result = rS->solve(desiredResult);
+ if (result == true) {
+ delete rS;
+ delete iS;
+ conservative--;
+ return (true);
+ }
+ EQs[0].coef[0]--;
+ }
+ if (j + 1 < lowerBounds) {
+ nEQs = 0;
+ eqnncpy(&GEQs[e], &EQs[0], nVars);
+ GEQs[e].touched = 1;
+ GEQs[e].color = EQ_BLACK;
+ *rS = *this;
+ if (DEBUG)
+ fprintf(outputFile, "exhausted lower bound, checking if still feasible ");
+ result = rS->solve(false);
+ if (result == false)
+ break;
+ }
+ }
+ if ((doTrace && desiredResult == OC_SOLVE_SIMPLIFY) || DBUG)
+ fprintf(outputFile, "fall-off the end\n");
+ delete rS;
+ delete iS;
+
+ conservative--;
+ return (false);
+ }
+ }
+ delete rS;
+ delete iS;
+ }
+ return (OC_SOLVE_UNKNOWN);
+ }
+ }
+ while (eliminateAgain);
+ }
+}
+
+
+int Problem::parallelSplinter(int e, int diff, int desiredResult) {
+ Problem *tmpProblem;
+ int i;
+ if (DBUG) {
+ fprintf(outputFile, "Using parallel splintering\n");
+ printProblem();
+ }
+ tmpProblem = new Problem;
+ int neweq = newEQ();
+ assert(neweq == 0);
+ eqnncpy(&EQs[0], &GEQs[e], nVars);
+ for (i = 0; i <= diff; i++) {
+ *tmpProblem = * this;
+ tmpProblem->isTemporary = true;
+ if (DBUG) {
+ fprintf(outputFile, "Splinter # %i\n", i);
+ printProblem();
+ }
+ if (tmpProblem->solve(desiredResult)) {
+ delete tmpProblem;
+ return true;
+ }
+ EQs[0].coef[0]--;
+ }
+ delete tmpProblem;
+ return false;
+}
+
+
+int Problem::verifyProblem() {
+ int result;
+ int e;
+ int areRed;
+ check();
+ Problem tmpProblem(*this);
+ tmpProblem.varsOfInterest = 0;
+ tmpProblem.safeVars = 0;
+ tmpProblem.nSUBs = 0;
+ tmpProblem.nMemories = 0;
+ tmpProblem.isTemporary = true;
+ areRed = 0;
+ if (mayBeRed) {
+ for(e=0; e<nEQs; e++) if (EQs[e].color) areRed = 1;
+ for(e=0; e<nGEQs; e++) if (GEQs[e].color) areRed = 1;
+ if (areRed) tmpProblem.turnRedBlack();
+ }
+ originalProblem = this;
+ assert(!outerColor);
+ outerColor = areRed;
+ if (TRACE) {
+ fprintf(outputFile, "verifying problem: [\n");
+ printProblem();
+ }
+ tmpProblem.check();
+ tmpProblem.freeEliminations(0);
+ result = tmpProblem.solve(OC_SOLVE_UNKNOWN);
+ originalProblem = noProblem;
+ outerColor = 0;
+ if (TRACE) {
+ if (result)
+ fprintf(outputFile, "] verified problem\n");
+ else
+ fprintf(outputFile, "] disproved problem\n");
+ printProblem();
+ }
+ check();
+ return result;
+}
+
+
+void Problem:: freeEliminations(int fv) {
+ int tryAgain = 1;
+ int i, e, e2;
+ while (tryAgain) {
+ tryAgain = 0;
+ for (i = nVars; i > fv; i--) {
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (GEQs[e].coef[i])
+ break;
+ if (e < 0)
+ e2 = e;
+ else if (GEQs[e].coef[i] > 0) {
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (GEQs[e2].coef[i] < 0)
+ break;
+ }
+ else {
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (GEQs[e2].coef[i] > 0)
+ break;
+ }
+ if (e2 < 0) {
+ int e3;
+ for (e3 = nSUBs - 1; e3 >= 0; e3--)
+ if (SUBs[e3].coef[i])
+ break;
+ if (e3 >= 0)
+ continue;
+ for (e3 = nEQs - 1; e3 >= 0; e3--)
+ if (EQs[e3].coef[i])
+ break;
+ if (e3 >= 0)
+ continue;
+ if (DBUG)
+ fprintf(outputFile, "a free elimination of %s (%d)\n", variable(i),e);
+ if (e >= 0) {
+ deleteGEQ(e);
+ for (e--; e >= 0; e--)
+ if (GEQs[e].coef[i]) {
+ deleteGEQ(e);
+ }
+ tryAgain = (i < nVars);
+ }
+ deleteVariable(i);
+ }
+ }
+ }
+
+ if (DEBUG) {
+ fprintf(outputFile, "\nafter free eliminations:\n");
+ printProblem();
+ fprintf(outputFile, "\n");
+ }
+}
+
+
+void Problem::freeRedEliminations() {
+ int tryAgain = 1;
+ int i, e, e2;
+ int isRedVar[maxVars];
+ int isDeadVar[maxVars];
+ int isDeadGEQ[maxmaxGEQs];
+ for (i = nVars; i > 0; i--) {
+ isRedVar[i] = 0;
+ isDeadVar[i] = 0;
+ }
+ for (e = nGEQs - 1; e >= 0; e--) {
+ isDeadGEQ[e] = 0;
+ if (GEQs[e].color)
+ for (i = nVars; i > 0; i--)
+ if (GEQs[e].coef[i] != 0)
+ isRedVar[i] = 1;
+ }
+
+ while (tryAgain) {
+ tryAgain = 0;
+ for (i = nVars; i > 0; i--)
+ if (!isRedVar[i] && !isDeadVar[i]) {
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (!isDeadGEQ[e] && GEQs[e].coef[i])
+ break;
+ if (e < 0)
+ e2 = e;
+ else if (GEQs[e].coef[i] > 0) {
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (!isDeadGEQ[e2] && GEQs[e2].coef[i] < 0)
+ break;
+ }
+ else {
+ for (e2 = e - 1; e2 >= 0; e2--)
+ if (!isDeadGEQ[e2] && GEQs[e2].coef[i] > 0)
+ break;
+ }
+ if (e2 < 0) {
+ int e3;
+ for (e3 = nSUBs - 1; e3 >= 0; e3--)
+ if (SUBs[e3].coef[i])
+ break;
+ if (e3 >= 0)
+ continue;
+ for (e3 = nEQs - 1; e3 >= 0; e3--)
+ if (EQs[e3].coef[i])
+ break;
+ if (e3 >= 0)
+ continue;
+ if (DBUG)
+ fprintf(outputFile, "a free red elimination of %s\n", variable(i));
+ for (; e >= 0; e--)
+ if (GEQs[e].coef[i])
+ isDeadGEQ[e] = 1;
+ tryAgain = 1;
+ isDeadVar[i] = 1;
+ }
+ }
+ }
+
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (isDeadGEQ[e])
+ deleteGEQ(e);
+
+ for (i = nVars; i > safeVars; i--)
+ if (isDeadVar[i])
+ deleteVariable(i);
+
+
+ if (DEBUG) {
+ fprintf(outputFile, "\nafter free red eliminations:\n");
+ printProblem();
+ fprintf(outputFile, "\n");
+ }
+}
+
+} // namespace
diff --git a/omegalib/omega/src/omega_core/oc_util.cc b/omegalib/omega/src/omega_core/oc_util.cc
new file mode 100644
index 0000000..a7d21be
--- /dev/null
+++ b/omegalib/omega/src/omega_core/oc_util.cc
@@ -0,0 +1,327 @@
+/*****************************************************************************
+ Copyright (C) 1994-2000 the Omega Project Team
+ Copyright (C) 2005-2011 Chun Chen
+ All Rights Reserved.
+
+ Purpose:
+
+ Notes:
+
+ History:
+*****************************************************************************/
+
+#include <omega/omega_core/oc_i.h>
+#include <algorithm>
+
+namespace omega {
+
+void Problem:: problem_merge(Problem &p2) {
+ int newLocation[maxVars];
+ int i,e2;
+
+ resurrectSubs();
+ p2.resurrectSubs();
+ setExternals();
+ p2.setExternals();
+
+ assert(safeVars == p2.safeVars);
+ if(DBUG) {
+ fprintf(outputFile,"Merging:\n");
+ printProblem();
+ fprintf(outputFile,"and\n");
+ p2.printProblem();
+ }
+ for(i=1; i<= p2.safeVars; i++) {
+ assert(p2.var[i] > 0) ;
+ newLocation[i] = forwardingAddress[p2.var[i]];
+ }
+ for(; i<= p2.nVars; i++) {
+ int j = ++(nVars);
+ newLocation[i] = j;
+ zeroVariable(j);
+ var[j] = -1;
+ }
+ newLocation[0] = 0;
+
+ for (e2 = p2.nEQs - 1; e2 >= 0; e2--) {
+ int e1 = newEQ();
+ eqnnzero(&(EQs[e1]), nVars);
+ for(i=0;i<=p2.nVars;i++)
+ EQs[e1].coef[newLocation[i]] = p2.EQs[e2].coef[i];
+ }
+ for (e2 = p2.nGEQs - 1; e2 >= 0; e2--) {
+ int e1 = newGEQ();
+ eqnnzero(&(GEQs[e1]), nVars);
+ GEQs[e1].touched = 1;
+ for(i=0;i<=p2.nVars;i++)
+ GEQs[e1].coef[newLocation[i]] = p2.GEQs[e2].coef[i];
+ }
+ int w = -1;
+ for (i = 1; i <= nVars; i++)
+ if (var[i] < 0) var[i] = w--;
+ if(DBUG) {
+ fprintf(outputFile,"to get:\n");
+ printProblem();
+ }
+}
+
+
+
+void Problem::chainUnprotect() {
+ int i, e;
+ int unprotect[maxVars];
+ int any = 0;
+ for (i = 1; i <= safeVars; i++) {
+ unprotect[i] = (var[i] < 0);
+ for (e = nSUBs - 1; e >= 0; e--)
+ if (SUBs[e].coef[i])
+ unprotect[i] = 0;
+ }
+ for (i = 1; i <= safeVars; i++) if (unprotect[i]) any=1;
+ if (!any) return;
+
+ if (DBUG) {
+ fprintf(outputFile, "Doing chain reaction unprotection\n");
+ printProblem();
+ for (i = 1; i <= safeVars; i++)
+ if (unprotect[i])
+ fprintf(outputFile, "unprotecting %s\n", variable(i));
+ }
+ for (i = 1; i <= safeVars; i++)
+ if (unprotect[i]) {
+ /* wild card */
+ if (i < safeVars) {
+ int j = safeVars;
+ std::swap(var[i], var[j]);
+ for (e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].touched = 1;
+ std::swap(GEQs[e].coef[i], GEQs[e].coef[j]);
+ }
+ for (e = nEQs - 1; e >= 0; e--)
+ std::swap(EQs[e].coef[i], EQs[e].coef[j]);
+ for (e = nSUBs - 1; e >= 0; e--)
+ std::swap(SUBs[e].coef[i], SUBs[e].coef[j]);
+ std::swap(unprotect[i], unprotect[j]);
+ i--;
+ }
+ safeVars--;
+ }
+ if (DBUG) {
+ fprintf(outputFile, "After chain reactions\n");
+ printProblem();
+ }
+}
+
+void Problem::resurrectSubs() {
+ if (nSUBs > 0 && !pleaseNoEqualitiesInSimplifiedProblems) {
+ int i, e, n, m,mbr;
+ mbr = 0;
+ for (e = nGEQs - 1; e >= 0; e--) if (GEQs[e].color) mbr=1;
+ for (e = nEQs - 1; e >= 0; e--) if (EQs[e].color) mbr=1;
+ if (nMemories) mbr = 1;
+
+ assert(!mbr || mayBeRed);
+
+ if (DBUG) {
+ fprintf(outputFile, "problem reduced, bringing variables back to life\n");
+ if(mbr && !mayBeRed) fprintf(outputFile, "Red equations we don't expect\n");
+ printProblem();
+ }
+ if (DBUG && nEQs > 0)
+ fprintf(outputFile,"This is wierd: problem has equalities\n");
+
+ for (i = 1; i <= safeVars; i++)
+ if (var[i] < 0) {
+ /* wild card */
+ if (i < safeVars) {
+ int j = safeVars;
+ std::swap(var[i], var[j]);
+ for (e = nGEQs - 1; e >= 0; e--) {
+ GEQs[e].touched = 1;
+ std::swap(GEQs[e].coef[i], GEQs[e].coef[j]);
+ }
+ for (e = nEQs - 1; e >= 0; e--)
+ std::swap(EQs[e].coef[i], EQs[e].coef[j]);
+ for (e = nSUBs - 1; e >= 0; e--)
+ std::swap(SUBs[e].coef[i], SUBs[e].coef[j]);
+ i--;
+ }
+ safeVars--;
+ }
+
+ m = nSUBs;
+ n = nVars;
+ if (n < safeVars + m)
+ n = safeVars + m;
+ for (e = nGEQs - 1; e >= 0; e--) {
+ if (singleVarGEQ(&GEQs[e])) {
+ i = abs(GEQs[e].key);
+ if (i >= safeVars + 1)
+ GEQs[e].key += (GEQs[e].key > 0 ? m : -m);
+ }
+ else {
+ GEQs[e].touched = true;
+ GEQs[e].key = 0;
+ }
+ }
+ for (i = nVars; i >= safeVars + 1; i--) {
+ var[i + m] = var[i];
+ for (e = nGEQs - 1; e >= 0; e--)
+ GEQs[e].coef[i + m] = GEQs[e].coef[i];
+ for (e = nEQs - 1; e >= 0; e--)
+ EQs[e].coef[i + m] = EQs[e].coef[i];
+ for (e = nSUBs - 1; e >= 0; e--)
+ SUBs[e].coef[i + m] = SUBs[e].coef[i];
+ }
+ for (i = safeVars + m; i >= safeVars + 1; i--) {
+ for (e = nGEQs - 1; e >= 0; e--) GEQs[e].coef[i] = 0;
+ for (e = nEQs - 1; e >= 0; e--) EQs[e].coef[i] = 0;
+ for (e = nSUBs - 1; e >= 0; e--) SUBs[e].coef[i] = 0;
+ }
+ nVars += m;
+ safeVars += m;
+ for (e = nSUBs - 1; e >= 0; e--)
+ var[safeVars -m + 1 + e] = SUBs[e].key;
+ for (i = 1; i <= safeVars; i++)
+ forwardingAddress[var[i]] = i;
+ if (DBUG) {
+ fprintf(outputFile,"Ready to wake substitutions\n");
+ printProblem();
+ }
+ for (e = nSUBs - 1; e >= 0; e--) {
+ int neweq = newEQ();
+ eqnncpy(&(EQs[neweq]), &(SUBs[e]), nVars);
+ EQs[neweq].coef[safeVars -m + 1 + e] = -1;
+ EQs[neweq].color = EQ_BLACK;
+ if (DBUG) {
+ fprintf(outputFile, "brought back: ");
+ printEQ(&EQs[neweq]);
+ fprintf(outputFile, "\n");
+ }
+ }
+ nSUBs = 0;
+
+ if (DBUG) {
+ fprintf(outputFile, "variables brought back to life\n");
+ printProblem();
+ }
+ }
+
+ coalesce();
+ recallRedMemories();
+ cleanoutWildcards();
+}
+
+
+void Problem::reverseProtectedVariables() {
+ int v1,v2,e,i;
+ coef_t t;
+ for (v1 = 1; v1 <= safeVars; v1++) {
+ v2 = safeVars+1-v1;
+ if (v2>=v1) break;
+ for(e=0;e<nEQs;e++) {
+ t = EQs[e].coef[v1];
+ EQs[e].coef[v1] = EQs[e].coef[v2];
+ EQs[e].coef[v2] = t;
+ }
+ for(e=0;e<nGEQs;e++) {
+ t = GEQs[e].coef[v1];
+ GEQs[e].coef[v1] = GEQs[e].coef[v2];
+ GEQs[e].coef[v2] = t;
+ GEQs[e].touched = 1;
+ }
+ for(e=0;e<nSUBs;e++) {
+ t = SUBs[e].coef[v1];
+ SUBs[e].coef[v1] = SUBs[e].coef[v2];
+ SUBs[e].coef[v2] = t;
+ }
+ }
+
+ for (i = 1; i <= safeVars; i++)
+ forwardingAddress[var[i]] = i;
+ for (i = 0; i < nSUBs; i++)
+ forwardingAddress[SUBs[i].key] = -i - 1;
+}
+
+void Problem::ordered_elimination(int symbolic) {
+ int i,j,e;
+ int isDead[maxmaxGEQs];
+ for(e=0;e<nEQs;e++) isDead[e] = 0;
+
+ if (!variablesInitialized) {
+ initializeVariables();
+ }
+
+ for(i=nVars;i>symbolic;i--)
+ for(e=0;e<nEQs;e++)
+ if (EQs[e].coef[i] == 1 || EQs[e].coef[i] == -1) {
+ for(j=nVars;j>i;j--) if (EQs[e].coef[j]) break;
+ if (i==j) {
+ doElimination(e, i);
+ isDead[e] = 1;
+ break;
+ }
+ }
+
+ for(e=nEQs-1;e>=0;e--)
+ if (isDead[e]) {
+ nEQs--;
+ if (e < nEQs) eqnncpy(&EQs[e], &EQs[nEQs], nVars);
+ }
+
+ for (i = 1; i <= safeVars; i++)
+ forwardingAddress[var[i]] = i;
+ for (i = 0; i < nSUBs; i++)
+ forwardingAddress[SUBs[i].key] = -i - 1;
+}
+
+
+coef_t Problem::checkSum() const {
+ coef_t cs;
+ int e;
+ cs = 0;
+ for(e=0;e<nGEQs;e++) {
+ coef_t c = GEQs[e].coef[0];
+ cs += c*c*c;
+ }
+ return cs;
+}
+
+
+void Problem::coalesce() {
+ int e, e2, colors;
+ int isDead[maxmaxGEQs];
+ int foundSomething = 0;
+
+
+ colors = 0;
+ for (e = 0; e < nGEQs; e++)
+ if (GEQs[e].color)
+ colors++;
+ if (colors < 2)
+ return;
+ for (e = 0; e < nGEQs; e++)
+ isDead[e] = 0;
+ for (e = 0; e < nGEQs; e++)
+ if (!GEQs[e].touched)
+ for (e2 = e + 1; e2 < nGEQs; e2++)
+ if (!GEQs[e2].touched && GEQs[e].key == -GEQs[e2].key
+ && GEQs[e].coef[0] == -GEQs[e2].coef[0]) {
+ int neweq = newEQ();
+ eqnncpy(&EQs[neweq], &GEQs[e], nVars);
+ EQs[neweq].color = GEQs[e].color || GEQs[e2].color;
+ foundSomething++;
+ isDead[e] = 1;
+ isDead[e2] = 1;
+ }
+ for (e = nGEQs - 1; e >= 0; e--)
+ if (isDead[e]) {
+ deleteGEQ(e);
+ }
+ if (DEBUG && foundSomething) {
+ fprintf(outputFile, "Coalesced GEQs into %d EQ's:\n", foundSomething);
+ printProblem();
+ }
+}
+
+} // namespace
generated by cgit v1.2.3-70-g09d2 (git 2.48.1) at 2025-02-02 03:48:21 +0000