Merge branch 'master' into polysat

src/math/lp/.clang-format

@@ -0,0 +1,4 @@
+BasedOnStyle: Google
+IndentWidth: 4
+ColumnLimit: 0
+NamespaceIndentation: All

src/math/lp/bound_analyzer_on_row.h

@@ -19,6 +19,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include "util/vector.h"

src/math/lp/column_info.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 
 #pragma once
 #include "util/vector.h"

src/math/lp/column_namer.h

@@ -18,6 +18,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include <string>
 #include "math/lp/static_matrix.h"
 namespace lp {

src/math/lp/core_solver_pretty_printer.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include "math/lp/numeric_pair.h"
 #include "math/lp/core_solver_pretty_printer_def.h"
 template lp::core_solver_pretty_printer<lp::mpq, lp::mpq>::core_solver_pretty_printer(const lp::lp_core_solver_base<lp::mpq, lp::mpq> &, std::ostream & out);

src/math/lp/core_solver_pretty_printer.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include <limits>
 #include <string>

src/math/lp/core_solver_pretty_printer_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <limits>
@@ -279,9 +280,9 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print()
         print_row(i);
     }
     m_out << std::endl;
-    if (m_core_solver.inf_set().size()) {
+    if (m_core_solver.inf_heap().size()) {
         m_out << "inf columns: ";
-        print_vector(m_core_solver.inf_set(), m_out);
+        print_vector(m_core_solver.inf_heap(), m_out);
         m_out << std::endl;
     }
 }

src/math/lp/cross_nested.h

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 #include <functional>
 #include "math/lp/nex.h"

src/math/lp/dense_matrix.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include "math/lp/lp_settings.h"
 #include "math/lp/dense_matrix_def.h"
 #ifdef Z3DEBUG

src/math/lp/dense_matrix.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #ifdef Z3DEBUG
 #include "util/vector.h"

src/math/lp/dense_matrix_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/lp_settings.h"

src/math/lp/emonics.cpp

@@ -18,6 +18,7 @@
   replaced rooted_mons.h and rooted_mon, rooted_mon_tabled
 
 --*/
+// clang-format off
 
 #include "math/lp/emonics.h"
 #include "math/lp/nla_defs.h"

src/math/lp/emonics.h

@@ -18,6 +18,7 @@
     to replace rooted_mons.h and rooted_mon, rooted_mon_tabled
 
   --*/
+// clang-format off
 
 #pragma once
 #include "math/lp/lp_utils.h"

src/math/lp/explanation.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/lp_utils.h"
 #include "util/map.h"
@@ -48,6 +49,15 @@ public:
         SASSERT(m_vector.empty());
         m_set.insert(j);
     }
+
+    void remove(constraint_index j) {
+        m_set.remove(j);
+        unsigned i = 0;
+        for (auto& p : m_vector) 
+            if (p.first != j)
+                m_vector[i++] = p;
+        m_vector.shrink(i);
+    }
     
     void add_expl(const explanation& e) {
         if (e.m_vector.empty()) {

src/math/lp/factorization.h

@@ -18,6 +18,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 #include "util/rational.h"
 #include "math/lp/monic.h"

src/math/lp/factorization_factory_imp.cpp

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #include "math/lp/factorization_factory_imp.h"
 #include "math/lp/nla_core.h"
 namespace nla {

src/math/lp/factorization_factory_imp.h

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 #include "math/lp/factorization.h"
 namespace nla {

src/math/lp/general_matrix.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include <functional>
 namespace lp {

src/math/lp/gomory.cpp

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #include "math/lp/gomory.h"
 #include "math/lp/int_solver.h"
 #include "math/lp/lar_solver.h"

src/math/lp/gomory.h

@@ -15,6 +15,7 @@ Author:
 
 Revision History:
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/lar_term.h"
 #include "math/lp/lia_move.h"

src/math/lp/hnf.h

@@ -26,6 +26,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/numeric_pair.h"
 #include "util/ext_gcd.h"

src/math/lp/hnf_cutter.cpp

@@ -9,6 +9,7 @@ Author:
     Lev Nachmanson (levnach)
 
 --*/
+// clang-format off
 #include "math/lp/int_solver.h"
 #include "math/lp/lar_solver.h"
 #include "math/lp/hnf_cutter.h"

src/math/lp/hnf_cutter.h

@@ -14,6 +14,7 @@ Author:
     Lev Nachmanson (levnach)
 
 --*/
+// clang-format off
 
 #pragma once
 #include "math/lp/lar_term.h"

src/math/lp/horner.cpp

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 
 #include "math/lp/horner.h"
 #include "math/lp/nla_core.h"

src/math/lp/horner.h

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/nla_common.h"

src/math/lp/implied_bound.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/lp_settings.h"
 #include "math/lp/lar_constraints.h"

src/math/lp/incremental_vector.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 
 #pragma once
 #include "util/vector.h"

src/math/lp/indexed_value.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 namespace lp {

src/math/lp/indexed_vector.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include "util/vector.h"
 #include "math/lp/indexed_vector_def.h"
 namespace lp {

src/math/lp/indexed_vector.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 
 #pragma once
 #include "util/vector.h"

src/math/lp/indexed_vector_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include "util/vector.h"

src/math/lp/int_branch.cpp

@@ -15,7 +15,7 @@
 
   Revision History:
   --*/
-
+// clang-format off
 #include "math/lp/int_solver.h"
 #include "math/lp/lar_solver.h"
 #include "math/lp/int_branch.h"
@@ -63,7 +63,7 @@ int int_branch::find_inf_int_base_column() {
     mpq small_value(1024);
     unsigned n = 0;
     lar_core_solver & lcs = lra.m_mpq_lar_core_solver;
-    unsigned prev_usage = 0; // to quiet down the compile
+    unsigned prev_usage = 0; // to quiet down the compiler
     unsigned k = 0;
     unsigned usage;
     unsigned j;

src/math/lp/int_branch.h

@@ -15,6 +15,7 @@ Author:
 
 Revision History:
 --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/lar_solver.h"

src/math/lp/int_cube.cpp

@@ -15,6 +15,7 @@ Author:
 
 Revision History:
 --*/
+// clang-format off
 
 #include "math/lp/int_solver.h"
 #include "math/lp/lar_solver.h"

src/math/lp/int_cube.h

@@ -19,6 +19,7 @@ Author:
 
 Revision History:
 --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/lia_move.h"

src/math/lp/int_gcd_test.cpp

@@ -45,6 +45,7 @@ Accumulative:
  
 
 --*/
+// clang-format off
 
 #include "math/lp/int_solver.h"
 #include "math/lp/lar_solver.h"

src/math/lp/int_gcd_test.h

@@ -24,6 +24,7 @@ Author:
 
 Revision History:
 --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/lia_move.h"

src/math/lp/int_solver.cpp

@@ -2,8 +2,7 @@
   Copyright (c) 2017 Microsoft Corporation
   Author: Lev Nachmanson
 */
-
-#include <utility>
+// clang-format off
 #include "math/lp/int_solver.h"
 #include "math/lp/lar_solver.h"
 #include "math/lp/lp_utils.h"
@@ -14,56 +13,204 @@
 
 namespace lp {
 
-int_solver::patcher::patcher(int_solver& lia):
-    lia(lia),
-    lra(lia.lra),
-    lrac(lia.lrac),
-    m_num_nbasic_patches(0),
-    m_patch_cost(0),
-    m_next_patch(0),
-    m_delay(0)
-{}
+    int_solver::patcher::patcher(int_solver& lia):
+        lia(lia),
+        lra(lia.lra),
+        lrac(lia.lrac)
+    {}
+    
+    void int_solver::patcher::remove_fixed_vars_from_base() {
+        unsigned num = lra.A_r().column_count();
+        for (unsigned v = 0; v < num; v++) {
+            if (!lia.is_base(v) || !lia.is_fixed(v))
+                continue;
+            auto const & r = lra.basic2row(v);
+            for (auto const& c : r) {
+                if (c.var() != v && !lia.is_fixed(c.var())) {
+                    lra.pivot(c.var(), v); 
+                    break;
+                }
+            }        
+        }
+    }
+
+
+    unsigned int_solver::patcher::count_non_int() {
+        unsigned non_int = 0;
+        for (auto j : lra.r_basis()) 
+            if (lra.column_is_int(j) && !lra.column_value_is_int(j))
+                ++non_int;
+        return non_int;
+    }
+
+    lia_move int_solver::patcher::patch_basic_columns() {
+        remove_fixed_vars_from_base();
+        lia.settings().stats().m_patches++;
+        lp_assert(lia.is_feasible());
+        
+       // unsigned non_int_before, non_int_after;
+
+       // non_int_before = count_non_int();
+
+
+        // unsigned num = lra.A_r().column_count();
+        for (unsigned j : lra.r_basis()) 
+            if (!lra.get_value(j).is_int())
+                patch_basic_column(j);
+        // non_int_after = count_non_int();
+        // verbose_stream() << non_int_before << " -> " << non_int_after << "\n";
+
+        if (!lia.has_inf_int()) {
+            lia.settings().stats().m_patches_success++;
+            return lia_move::sat;
+        }
+        return lia_move::undef;
+    }
+
+    // clang-format on
+    /**
+     * \brief find integral and minimal, in the absolute values, deltas such that x - alpha*delta is integral too.
+     */
+    bool get_patching_deltas(const rational& x, const rational& alpha,
+                             rational& delta_plus, rational& delta_minus) {
+        auto a1 = numerator(alpha);
+        auto a2 = denominator(alpha);
+        auto x1 = numerator(x);
+        auto x2 = denominator(x);
+        if (!divides(x2, a2))
+            return false;
+
+        // delta has to be integral.
+        // We need to find delta such that x1/x2 + (a1/a2)*delta is integral (we are going to flip the delta sign later).
+        // Then a2*x1/x2 + a1*delta is integral, but x2 and x1 are coprime:
+        // that means that t = a2/x2 is
+        // integral. We established that a2 = x2*t Then x1 + a1*delta*(x2/a2)  = x1
+        // + a1*(delta/t)  is integral. Taking into account that t and a1 are
+        // coprime we have delta = t*k, where k is an integer.
+        rational t = a2 / x2;
+        // std::cout << "t = " << t << std::endl;
+        // Now we have x1/x2 + (a1/x2)*k is integral, or (x1 + a1*k)/x2 is integral.
+        // It is equivalent to x1 + a1*k = x2*m, where m is an integer
+        // We know that a2 and a1 are coprime, and x2 divides a2, so x2 and a1 are
+        // coprime. We can find u and v such that u*a1 + v*x2 = 1.
+        rational u, v;
+        gcd(a1, x2, u, v);
+        lp_assert(gcd(a1, x2, u, v).is_one());
+        // std::cout << "u = " << u << ", v = " << v << std::endl;
+        // std::cout << "x= " << (x1 / x2) << std::endl;
+        // std::cout << "x + (a1 / a2) * (-u * t) * x1 = "
+        //           << x + (a1 / a2) * (-u * t) * x1 << std::endl;
+        lp_assert((x + (a1 / a2) * (-u * t) * x1).is_int());
+        // 1 = (u- l*x2 ) * a1 + (v + l*a1)*x2, for every integer l.
+        rational d = u * t * x1;
+        // We can prove that x+alpha*d is integral,
+        // and any other delta, satisfying x+alpha*delta, is equal to d modulo a2.
+        delta_plus = mod(d, a2);
+        lp_assert(delta_plus > 0);
+        delta_minus = delta_plus - a2;
+        lp_assert(delta_minus < 0);
 
-bool int_solver::patcher::should_apply() {
-#if 1
-    return true;
-#else
-    if (m_delay == 0) {
         return true;
     }
-    --m_delay;
-    return false;
-#endif
-}
+    /**
+     * \brief try to patch the basic column v
+     */
+    bool int_solver::patcher::patch_basic_column_on_row_cell(unsigned v, row_cell<mpq> const& c) {
+        if (v == c.var())
+            return false;
+        if (!lra.column_is_int(c.var()))  // could use real to patch integer
+            return false;
+        if (c.coeff().is_int())
+            return false;
+        mpq a = fractional_part(c.coeff());
+        mpq r = fractional_part(lra.get_value(v));
+        lp_assert(0 < r && r < 1);
+        lp_assert(0 < a && a < 1);
+        mpq delta_plus, delta_minus;
+        if (!get_patching_deltas(r, a, delta_plus, delta_minus))
+            return false;
 
-lia_move int_solver::patcher::operator()() {
-    return patch_nbasic_columns();
-}
+        if (lia.random() % 2) {
+            return try_patch_column(v, c.var(), delta_plus) ||
+                   try_patch_column(v, c.var(), delta_minus);
+        } else {
+            return try_patch_column(v, c.var(), delta_minus) ||
+                   try_patch_column(v, c.var(), delta_plus);
+        }
+    }
+    // clang-format off
+    
+    bool int_solver::patcher::try_patch_column(unsigned v, unsigned j, mpq const& delta) {
+        const auto & A = lra.A_r();
+        if (delta < 0) {
+            if (lia.has_lower(j) && lia.get_value(j) + impq(delta) < lra.get_lower_bound(j))
+                return false;
+        }
+        else {
+            if (lia.has_upper(j) && lia.get_value(j) + impq(delta) > lra.get_upper_bound(j))
+                return false;
+        }
+        for (auto const& c : A.column(j)) {
+            unsigned row_index = c.var();
+            unsigned i = lrac.m_r_basis[row_index];
+            auto old_val = lia.get_value(i);
+            auto new_val = old_val - impq(c.coeff()*delta);
+            if (lia.has_lower(i) && new_val < lra.get_lower_bound(i))
+                return false;
+            if (lia.has_upper(i) && new_val > lra.get_upper_bound(i))
+                return false;
+            if (old_val.is_int() && !new_val.is_int()){
+                return false; // do not waste resources on this case
+            }
+            lp_assert(i != v || new_val.is_int())
+        }
+        
+        lra.set_value_for_nbasic_column(j, lia.get_value(j) + impq(delta));
+         
+        return true;
+    }
+    
+    void int_solver::patcher::patch_basic_column(unsigned v) {
+        SASSERT(!lia.is_fixed(v));
+        for (auto const& c : lra.basic2row(v))
+            if (patch_basic_column_on_row_cell(v, c))
+                return;                                       
+    }
 
-lia_move int_solver::patcher::patch_nbasic_columns() {
-    lia.settings().stats().m_patches++;
-    lp_assert(lia.is_feasible());
-    m_num_nbasic_patches = 0;
-    m_patch_cost = 0;
-    for (unsigned j : lia.lrac.m_r_nbasis) {
-        patch_nbasic_column(j);
+    lia_move int_solver::patcher::patch_nbasic_columns() {
+        remove_fixed_vars_from_base();
+        lia.settings().stats().m_patches++;
+        lp_assert(lia.is_feasible());
+        m_patch_success = 0;
+        m_patch_fail = 0;
+        m_num_ones = 0;
+        m_num_divides = 0;
+        //unsigned non_int_before = count_non_int();
+
+        unsigned num = lra.A_r().column_count();
+        for (unsigned v = 0; v < num; v++) {
+            if (lia.is_base(v))
+                continue;
+            patch_nbasic_column(v);
+        }
+        unsigned num_fixed = 0;
+        for (unsigned v = 0; v < num; v++) 
+            if (lia.is_fixed(v))
+                ++num_fixed;
+            
+        lp_assert(lia.is_feasible());
+        //verbose_stream() << "patch " << m_patch_success << " fails " << m_patch_fail << " ones " << m_num_ones << " divides " << m_num_divides << " num fixed " << num_fixed << "\n";
+        //lra.display(verbose_stream());
+        //exit(0);
+        //unsigned non_int_after = count_non_int();
+
+        // verbose_stream() << non_int_before << " -> " << non_int_after << "\n";
+        if (!lia.has_inf_int()) {
+            lia.settings().stats().m_patches_success++;
+            return lia_move::sat;
+        }
+        return lia_move::undef;
     }
-    lp_assert(lia.is_feasible());
-    if (!lia.has_inf_int()) {
-        lia.settings().stats().m_patches_success++;
-        m_delay = 0;
-        m_next_patch = 0;
-        return lia_move::sat;
-    }
-    if (m_patch_cost > 0 && m_num_nbasic_patches * 10 < m_patch_cost) {
-        m_delay = std::min(20u, m_next_patch++);
-    }
-    else {
-        m_delay = 0;
-        m_next_patch = 0;
-    }
-    return lia_move::undef;
-}
 
 void int_solver::patcher::patch_nbasic_column(unsigned j) {
     impq & val = lrac.m_r_x[j];
@@ -71,17 +218,48 @@ void int_solver::patcher::patch_nbasic_column(unsigned j) {
     impq l, u;
     mpq m;
     bool has_free = lia.get_freedom_interval_for_column(j, inf_l, l, inf_u, u, m);
-    m_patch_cost += lra.A_r().number_of_non_zeroes_in_column(j);
-    if (!has_free) {
+    if (!has_free) 
         return;
-    }
     bool m_is_one = m.is_one();
     bool val_is_int = lia.value_is_int(j);
 
+#if 0
+const auto & A = lra.A_r();
+#endif
     // check whether value of j is already a multiple of m.
     if (val_is_int && (m_is_one || (val.x / m).is_int())) {
+        if (m_is_one)
+            ++m_num_ones;
+        else
+            ++m_num_divides;
+#if 0
+        for (auto c : A.column(j)) {
+            unsigned row_index = c.var();
+            unsigned i = lrac.m_r_basis[row_index];
+            if (!lia.get_value(i).is_int() ||
+                (lia.has_lower(i) && lia.get_value(i) < lra.get_lower_bound(i)) ||
+                (lia.has_upper(i) && lia.get_value(i) > lra.get_upper_bound(i))) {
+                verbose_stream() << "skip " << j << " " << m << ": ";
+                lia.display_row(verbose_stream(), A.m_rows[row_index]);
+                verbose_stream() << "\n";                    
+            }
+        }
+#endif
         return;
     }
+
+#if 0
+    if (!m_is_one) {
+        // lia.display_column(verbose_stream(), j);
+        for (auto c : A.column(j)) {
+            continue;
+            unsigned row_index = c.var();
+            lia.display_row(verbose_stream(), A.m_rows[row_index]);
+            verbose_stream() << "\n";
+        }
+    }
+#endif
+
     TRACE("patch_int",
           tout << "TARGET j" << j << " -> [";
           if (inf_l) tout << "-oo"; else tout << l;
@@ -89,9 +267,33 @@ void int_solver::patcher::patch_nbasic_column(unsigned j) {
           if (inf_u) tout << "oo"; else tout << u;
           tout << "]";
           tout << ", m: " << m << ", val: " << val << ", is_int: " << lra.column_is_int(j) << "\n";);
+
+#if 0
+    verbose_stream() << "path " << m << " ";
+    if (!inf_l) verbose_stream() << "infl " << l.x << " ";
+    if (!inf_u) verbose_stream() << "infu " << u.x << " ";
+    verbose_stream() << "\n";
     if (m.is_big() || (!inf_l && l.x.is_big()) || (!inf_u && u.x.is_big())) {
         return;
     }
+#endif
+
+#if 0
+    verbose_stream() << "TARGET v" << j << " -> [";
+    if (inf_l) verbose_stream() << "-oo"; else verbose_stream() << ceil(l.x) << " " << l << "\n";
+    verbose_stream() << ", ";
+    if (inf_u) verbose_stream() << "oo"; else verbose_stream() << floor(u.x) << " " << u << "\n";
+    verbose_stream() << "]";
+    verbose_stream() << ", m: " << m << ", val: " << val << ", is_int: " << lra.column_is_int(j) << "\n";
+#endif
+
+#if 0
+    if (!inf_l)
+        l = impq(ceil(l));
+    if (!inf_u)
+        u = impq(floor(u));
+#endif
+    
     if (!inf_l) {
         l = impq(m_is_one ? ceil(l) : m * ceil(l / m));
         if (inf_u || l <= u) {
@@ -99,8 +301,23 @@ void int_solver::patcher::patch_nbasic_column(unsigned j) {
             lra.set_value_for_nbasic_column(j, l);
         }
         else {
-            --m_num_nbasic_patches;
+            //verbose_stream() << "fail: " << j << " " << m << "\n";
+            ++m_patch_fail;
             TRACE("patch_int", tout << "not patching " << l << "\n";);
+#if 0
+            verbose_stream() << "not patched\n";
+            for (auto c : A.column(j)) {
+                unsigned row_index = c.var();
+                unsigned i = lrac.m_r_basis[row_index];
+                if (!lia.get_value(i).is_int() ||
+                    (lia.has_lower(i) && lia.get_value(i) < lra.get_lower_bound(i)) ||
+                    (lia.has_upper(i) && lia.get_value(i) > lra.get_upper_bound(i))) {
+                    lia.display_row(verbose_stream(), A.m_rows[row_index]);
+                    verbose_stream() << "\n";                    
+                }
+            }
+#endif
+            return;
         }
     }
     else if (!inf_u) {
@@ -112,7 +329,21 @@ void int_solver::patcher::patch_nbasic_column(unsigned j) {
         lra.set_value_for_nbasic_column(j, impq(0));
         TRACE("patch_int", tout << "patching with 0\n";);
     }
-    ++m_num_nbasic_patches;
+    ++m_patch_success;
+#if 0
+    verbose_stream() << "patched " << j << "\n";
+    for (auto c : A.column(j)) {
+        unsigned row_index = c.var();
+        unsigned i = lrac.m_r_basis[row_index];
+        if (!lia.get_value(i).is_int() ||
+            (lia.has_lower(i) && lia.get_value(i) < lra.get_lower_bound(i)) ||
+            (lia.has_upper(i) && lia.get_value(i) > lra.get_upper_bound(i))) {
+            lia.display_row(verbose_stream(), A.m_rows[row_index]);
+            verbose_stream() << "\n";                    
+        }
+    }
+#endif
+    
 }
 
 int_solver::int_solver(lar_solver& lar_slv) :
@@ -326,7 +557,7 @@ static void set_upper(impq & u, bool & inf_u, impq const & v) {
 
 // this function assumes that all basic columns dependend on j are feasible
 bool int_solver::get_freedom_interval_for_column(unsigned j, bool & inf_l, impq & l, bool & inf_u, impq & u, mpq & m) {
-    if (lrac.m_r_heading[j] >= 0) // the basic var
+    if (lrac.m_r_heading[j] >= 0 || is_fixed(j)) // basic or fixed var
         return false;
 
     TRACE("random_update", display_column(tout, j) << ", is_int = " << column_is_int(j) << "\n";);
@@ -361,10 +592,9 @@ bool int_solver::get_freedom_interval_for_column(unsigned j, bool & inf_l, impq
         unsigned i = lrac.m_r_basis[row_index];
         impq const & xi = get_value(i);
         lp_assert(lrac.m_r_solver.column_is_feasible(i));
-        if (column_is_int(i) && !a.is_int())
+        if (column_is_int(i) && !a.is_int() && xi.is_int()) 
             m = lcm(m, denominator(a));
-
-
+        
         if (!inf_l && !inf_u) {
             if (l == u) 
                 continue;            
@@ -372,15 +602,15 @@ bool int_solver::get_freedom_interval_for_column(unsigned j, bool & inf_l, impq
 
         if (a.is_neg()) {
             if (has_lower(i)) 
-                set_lower(l, inf_l, delta(a, xi, lrac.m_r_lower_bounds()[i]));
+                set_lower(l, inf_l, delta(a, xi, lra.get_lower_bound(i)));
             if (has_upper(i)) 
-                set_upper(u, inf_u, delta(a, xi, lrac.m_r_upper_bounds()[i]));
+                set_upper(u, inf_u, delta(a, xi, lra.get_upper_bound(i)));
         }
         else {
             if (has_upper(i)) 
-                set_lower(l, inf_l, delta(a, xi, lrac.m_r_upper_bounds()[i]));
+                set_lower(l, inf_l, delta(a, xi, lra.get_upper_bound(i)));
             if (has_lower(i)) 
-                set_upper(u, inf_u, delta(a, xi, lrac.m_r_lower_bounds()[i]));
+                set_upper(u, inf_u, delta(a, xi, lra.get_lower_bound(i)));
         }
     }
 
@@ -479,14 +709,11 @@ bool int_solver::at_upper(unsigned j) const {
 }
 
 std::ostream & int_solver::display_row(std::ostream & out, lp::row_strip<rational> const & row) const {
-bool first = true;
+    bool first = true;
     auto & rslv = lrac.m_r_solver;
-for (const auto &c : row)
-    {
-        if (is_fixed(c.var()))
-        {
-            if (!get_value(c.var()).is_zero())
-            {
+    for (const auto &c : row) {
+        if (is_fixed(c.var())) {
+            if (!get_value(c.var()).is_zero()) {
                 impq val = get_value(c.var()) * c.coeff();
                 if (!first && val.is_pos())
                     out << "+";
@@ -505,17 +732,11 @@ for (const auto &c : row)
         }
         else if (c.coeff().is_minus_one())
             out << "-";
-        else
-        {
-            if (c.coeff().is_pos())
-            {
-                if (!first)
-                    out << "+";
-            }
+        else {
+            if (c.coeff().is_pos() && !first)
+                out << "+";
             if (c.coeff().is_big())
-            {
                 out << " b*";
-            }
             else
                 out << c.coeff();
         }
@@ -523,8 +744,7 @@ for (const auto &c : row)
         first = false;
     }
     out << "\n";
-    for (const auto &c : row)
-    {
+    for (const auto &c : row) {
         if (is_fixed(c.var()))
             continue;
         rslv.print_column_info(c.var(), out);
@@ -533,14 +753,13 @@ for (const auto &c : row)
     }
     return out;
 }
+    
 std::ostream& int_solver::display_row_info(std::ostream & out, unsigned row_index) const  {    
     auto & rslv = lrac.m_r_solver;
     auto const& row = rslv.m_A.m_rows[row_index];
     return display_row(out, row);
 }
 
-
-
 bool int_solver::shift_var(unsigned j, unsigned range) {
     if (is_fixed(j) || is_base(j))
         return false;
@@ -549,11 +768,13 @@ bool int_solver::shift_var(unsigned j, unsigned range) {
     bool inf_l = false, inf_u = false;
     impq l, u;
     mpq m;
-    VERIFY(get_freedom_interval_for_column(j, inf_l, l, inf_u, u, m) || settings().get_cancel_flag());
+    if (!get_freedom_interval_for_column(j, inf_l, l, inf_u, u, m))
+        return false;
     if (settings().get_cancel_flag())
         return false;
     const impq & x = get_value(j);
     // x, the value of j column, might be shifted on a multiple of m
+
     if (inf_l && inf_u) {
         impq new_val = m * impq(random() % (range + 1)) + x;
         lra.set_value_for_nbasic_column(j, new_val);
@@ -570,6 +791,7 @@ bool int_solver::shift_var(unsigned j, unsigned range) {
     if (!inf_l && !inf_u && l >= u)
         return false;
 
+
     if (inf_u) {
         SASSERT(!inf_l);
         impq new_val = x + m * impq(random() % (range + 1));
@@ -640,21 +862,14 @@ int int_solver::select_int_infeasible_var() {
     unsigned n = 0;
     lar_core_solver & lcs = lra.m_mpq_lar_core_solver;
     unsigned prev_usage = 0; // to quiet down the compile
-    unsigned k = 0;
-    unsigned usage;
-    unsigned j;
 
     enum state { small_box, is_small_value,  any_value, not_found };
     state st = not_found;
 
-    // 1. small box
-    // 2. small value
-    // 3. any value
-    for (; k < lra.r_basis().size(); k++) {
-        j = lra.r_basis()[k];
+    for (unsigned j : lra.r_basis()) {
         if (!column_is_int_inf(j))
             continue;
-        usage = lra.usage_in_terms(j);
+        unsigned usage = lra.usage_in_terms(j);
         if (is_boxed(j) && (new_range = lcs.m_r_upper_bounds()[j].x - lcs.m_r_lower_bounds()[j].x - rational(2*usage)) <= small_value) {
             SASSERT(!is_fixed(j));
             if (st != small_box) {
@@ -688,12 +903,12 @@ int int_solver::select_int_infeasible_var() {
             continue;
         SASSERT(st == not_found || st == any_value);
         st = any_value;
-        if (n == 0 /*|| usage > prev_usage*/) {
+        if (n == 0 || usage > prev_usage) {
             result = j;
             prev_usage = usage;
             n = 1;
         } 
-        else if (usage > 0 && /*usage == prev_usage && */ (random() % (++n) == 0)) 
+        else if (usage > 0 && usage == prev_usage && (random() % (++n) == 0)) 
             result = j;
     }
     

src/math/lp/int_solver.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/lp_settings.h"
 #include "math/lp/static_matrix.h"
@@ -44,17 +45,23 @@ class int_solver {
         int_solver&         lia;
         lar_solver&         lra;
         lar_core_solver&    lrac;
-        unsigned            m_num_nbasic_patches;
-        unsigned            m_patch_cost;
-        unsigned            m_next_patch;
-        unsigned            m_delay;
+        unsigned            m_patch_success = 0;
+        unsigned            m_patch_fail = 0;
+        unsigned            m_num_ones = 0;
+        unsigned            m_num_divides = 0;
     public:
         patcher(int_solver& lia);
-        bool should_apply();
-        lia_move operator()();
+        bool should_apply() const { return true; }
+        lia_move operator()() { return patch_basic_columns(); }
         void patch_nbasic_column(unsigned j);
+        bool patch_basic_column_on_row_cell(unsigned v, row_cell<mpq> const& c);
+        void patch_basic_column(unsigned j);
+        bool try_patch_column(unsigned v, unsigned j, mpq const& delta);
+        unsigned count_non_int();
     private:
+        void remove_fixed_vars_from_base();
         lia_move patch_nbasic_columns();
+        lia_move patch_basic_columns();
     };
 
     lar_solver&         lra;

src/math/lp/lar_constraints.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 
 #pragma once
 #include <utility>

src/math/lp/lar_core_solver.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include <utility>
 #include <memory>
 #include <string>

src/math/lp/lar_core_solver.h

@@ -5,6 +5,7 @@ Author:
     Lev Nachmanson (levnach)
 
 --*/
+// clang-format off
 #pragma once
 #include "util/vector.h"
 #include <string>
@@ -16,7 +17,7 @@ Author:
 #include "math/lp/stacked_vector.h"
 #include "util/stacked_value.h"
 namespace lp {
-
+// clang-format off
 class lar_core_solver  {
     vector<std::pair<mpq, unsigned>> m_infeasible_linear_combination;
     int m_infeasible_sum_sign; // todo: get rid of this field
@@ -93,6 +94,8 @@ public:
 
     void solve();
 
+    void pivot(int entering, int leaving) { m_r_solver.pivot(entering, leaving); }
+    
     bool lower_bounds_are_set() const { return true; }
 
     const indexed_vector<mpq> & get_pivot_row() const {

src/math/lp/lar_core_solver_def.h

@@ -9,6 +9,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <string>
@@ -85,8 +86,8 @@ unsigned lar_core_solver::get_number_of_non_ints() const {
 void lar_core_solver::solve() {
     TRACE("lar_solver", tout << m_r_solver.get_status() << "\n";);
     lp_assert(m_r_solver.non_basic_columns_are_set_correctly());
-    lp_assert(m_r_solver.inf_set_is_correct());
-	TRACE("find_feas_stats", tout << "infeasibles = " << m_r_solver.inf_set_size() << ", int_infs = " << get_number_of_non_ints() << std::endl;);
+    lp_assert(m_r_solver.inf_heap_is_correct());
+	TRACE("find_feas_stats", tout << "infeasibles = " << m_r_solver.inf_heap_size() << ", int_infs = " << get_number_of_non_ints() << std::endl;);
 	if (m_r_solver.current_x_is_feasible() && m_r_solver.m_look_for_feasible_solution_only) {
             m_r_solver.set_status(lp_status::OPTIMAL);
             TRACE("lar_solver", tout << m_r_solver.get_status() << "\n";);
@@ -117,11 +118,9 @@ void lar_core_solver::solve() {
     }
     lp_assert(r_basis_is_OK());
     lp_assert(m_r_solver.non_basic_columns_are_set_correctly());
-    lp_assert(m_r_solver.inf_set_is_correct());
-
-    TRACE("lar_solver", tout << m_r_solver.get_status() << "\n";);
-}
-
+    lp_assert(m_r_solver.inf_heap_is_correct());
 
+  TRACE("lar_solver", tout << m_r_solver.get_status() << "\n";);
 }
 
+} // namespace lp

src/math/lp/lar_solver.cpp

@@ -2,7 +2,7 @@
   Copyright (c) 2017 Microsoft Corporation
   Author: Nikolaj Bjorner, Lev Nachmanson
 */
-
+// clang-format off
 #include "math/lp/lar_solver.h"
 #include "smt/params/smt_params_helper.hpp"
 
@@ -41,7 +41,6 @@ namespace lp {
         for (auto t : m_terms)
             delete t;
     }
-
     
     bool lar_solver::sizes_are_correct() const {
         lp_assert(A_r().column_count() == m_mpq_lar_core_solver.m_r_solver.m_column_types.size());
@@ -50,7 +49,6 @@ namespace lp {
         return true;
     }
 
-
     std::ostream& lar_solver::print_implied_bound(const implied_bound& be, std::ostream& out) const {
         out << "implied bound\n";
         unsigned v = be.m_j;
@@ -215,7 +213,7 @@ namespace lp {
 
     void lar_solver::fill_explanation_from_crossed_bounds_column(explanation& evidence) const {
         lp_assert(static_cast<int>(get_column_type(m_crossed_bounds_column)) >= static_cast<int>(column_type::boxed));
-        lp_assert(!m_mpq_lar_core_solver.m_r_solver.column_is_feasible(m_crossed_bounds_column));
+        lp_assert(!column_is_feasible(m_crossed_bounds_column));
 
         // this is the case when the lower bound is in conflict with the upper one
         const ul_pair& ul = m_columns_to_ul_pairs[m_crossed_bounds_column];
@@ -244,6 +242,14 @@ namespace lp {
             set.erase(j);
     }
 
+    void lar_solver::clean_popped_elements_for_heap(unsigned n, lpvar_heap& heap) {
+        vector<int> to_remove;
+        for (unsigned j : heap)
+            if (j >= n)
+                to_remove.push_back(j);
+        for (unsigned j : to_remove)
+            heap.erase(j);
+    }
     
 
     void lar_solver::pop(unsigned k) {
@@ -271,7 +277,7 @@ namespace lp {
 
         unsigned m = A_r().row_count();
         clean_popped_elements(m, m_rows_with_changed_bounds);
-        clean_inf_set_of_r_solver_after_pop();
+        clean_inf_heap_of_r_solver_after_pop();
         lp_assert(
             m_settings.simplex_strategy() == simplex_strategy_enum::undecided ||
             m_mpq_lar_core_solver.m_r_solver.reduced_costs_are_correct_tableau());
@@ -328,7 +334,7 @@ namespace lp {
 
     void lar_solver::set_costs_to_zero(const lar_term& term) {
         auto& rslv = m_mpq_lar_core_solver.m_r_solver;
-        auto& jset = m_mpq_lar_core_solver.m_r_solver.inf_set(); // hijack this set that should be empty right now
+        auto& jset = m_mpq_lar_core_solver.m_r_solver.inf_heap(); // hijack this set that should be empty right now
         lp_assert(jset.empty());
 
         for (lar_term::ival p : term) {
@@ -667,16 +673,16 @@ namespace lp {
            m_mpq_lar_core_solver.m_r_solver.add_delta_to_x_and_track_feasibility(bj, -A_r().get_val(c) * delta);
            TRACE("change_x_del",
               tout << "changed basis column " << bj << ", it is " <<
-              (m_mpq_lar_core_solver.m_r_solver.column_is_feasible(bj) ? "feas" : "inf") << std::endl;);
+              (column_is_feasible(bj) ? "feas" : "inf") << std::endl;);
         }
     }
 
     void lar_solver::update_x_and_inf_costs_for_column_with_changed_bounds(unsigned j) {
         if (m_mpq_lar_core_solver.m_r_heading[j] >= 0) {
             if (costs_are_used()) {
-                bool was_infeas = m_mpq_lar_core_solver.m_r_solver.inf_set_contains(j);
+                bool was_infeas = m_mpq_lar_core_solver.m_r_solver.inf_heap_contains(j);
                 m_mpq_lar_core_solver.m_r_solver.track_column_feasibility(j);
-                if (was_infeas != m_mpq_lar_core_solver.m_r_solver.inf_set_contains(j))
+                if (was_infeas != m_mpq_lar_core_solver.m_r_solver.inf_heap_contains(j))
                     m_basic_columns_with_changed_cost.insert(j);
             }
             else {
@@ -1293,12 +1299,12 @@ namespace lp {
         lp_assert(m_mpq_lar_core_solver.m_r_solver.basis_heading_is_correct());
     }
 
-    void lar_solver::clean_inf_set_of_r_solver_after_pop() {
+    void lar_solver::clean_inf_heap_of_r_solver_after_pop() {
         vector<unsigned> became_feas;
-        clean_popped_elements(A_r().column_count(), m_mpq_lar_core_solver.m_r_solver.inf_set());
+        clean_popped_elements_for_heap(A_r().column_count(), m_mpq_lar_core_solver.m_r_solver.inf_heap());
         std::unordered_set<unsigned> basic_columns_with_changed_cost;
         m_inf_index_copy.reset();
-        for (auto j : m_mpq_lar_core_solver.m_r_solver.inf_set())
+        for (auto j : m_mpq_lar_core_solver.m_r_solver.inf_heap())
             m_inf_index_copy.push_back(j);
         for (auto j : m_inf_index_copy) {
             if (m_mpq_lar_core_solver.m_r_heading[j] >= 0) {
@@ -1316,16 +1322,16 @@ namespace lp {
             lp_assert(m_mpq_lar_core_solver.m_r_solver.m_basis_heading[j] < 0);
             m_mpq_lar_core_solver.m_r_solver.m_d[j] -= m_mpq_lar_core_solver.m_r_solver.m_costs[j];
             m_mpq_lar_core_solver.m_r_solver.m_costs[j] = zero_of_type<mpq>();
-            m_mpq_lar_core_solver.m_r_solver.remove_column_from_inf_set(j);
+            m_mpq_lar_core_solver.m_r_solver.remove_column_from_inf_heap(j);
         }
         became_feas.clear();
-        for (unsigned j : m_mpq_lar_core_solver.m_r_solver.inf_set()) {
+        for (unsigned j : m_mpq_lar_core_solver.m_r_solver.inf_heap()) {
             lp_assert(m_mpq_lar_core_solver.m_r_heading[j] >= 0);
-            if (m_mpq_lar_core_solver.m_r_solver.column_is_feasible(j))
+            if (column_is_feasible(j))
                 became_feas.push_back(j);
         }
         for (unsigned j : became_feas)
-            m_mpq_lar_core_solver.m_r_solver.remove_column_from_inf_set(j);
+            m_mpq_lar_core_solver.m_r_solver.remove_column_from_inf_heap(j);
 
     }
 
@@ -1346,8 +1352,8 @@ namespace lp {
     }
 
     bool lar_solver::term_is_int(const vector<std::pair<mpq, unsigned int>>& coeffs) const {
-        for (auto const& p : coeffs)
-            if (!(column_is_int(p.second) && p.first.is_int()))
+        for (auto const& [coeff, v] : coeffs)
+            if (!(column_is_int(v) && coeff.is_int()))
                 return false;
         return true;
     }
@@ -1374,66 +1380,6 @@ namespace lp {
         return m_mpq_lar_core_solver.column_is_free(j);
     }
 
-    // column is at lower or upper bound, lower and upper bound are different.
-    // the lower/upper bound is not strict.
-    // the LP obtained by making the bound strict is infeasible
-    // -> the column has to be fixed
-    bool lar_solver::is_fixed_at_bound(column_index const& j) {
-        if (column_is_fixed(j))
-            return false;
-        mpq val;
-        if (!has_value(j, val))
-            return false;
-        lp::lconstraint_kind k;
-        if (column_has_upper_bound(j) &&
-            get_upper_bound(j).x == val) {
-            verbose_stream() << "check upper " << j << "\n";
-            push();
-            if (column_is_int(j))
-                k = LE, val -= 1;
-            else
-                k = LT;
-            auto ci = mk_var_bound(j, k, val);
-            update_column_type_and_bound(j, k, val, ci);
-            auto st = find_feasible_solution();
-            pop(1);
-            return st == lp_status::INFEASIBLE;            
-        }
-        if (column_has_lower_bound(j) &&
-            get_lower_bound(j).x == val) {
-            verbose_stream() << "check lower " << j << "\n";
-            push();
-            if (column_is_int(j))
-                k = GE, val += 1;
-            else
-                k = GT;
-            auto ci = mk_var_bound(j, k, val);
-            update_column_type_and_bound(j, k, val, ci);
-            auto st = find_feasible_solution();
-            pop(1);
-            return st == lp_status::INFEASIBLE;                        
-        }
-        
-        return false;
-    }
-
-    bool lar_solver::has_fixed_at_bound() {
-        verbose_stream() << "has-fixed-at-bound\n";
-        unsigned num_fixed = 0;
-        for (unsigned j = 0; j < A_r().m_columns.size(); ++j) {
-            auto ci = column_index(j);
-            if (is_fixed_at_bound(ci)) {
-                ++num_fixed;
-                verbose_stream() << "fixed " << j << "\n";
-            }
-        }
-        verbose_stream() << "num fixed " << num_fixed << "\n";
-        if (num_fixed > 0)
-            find_feasible_solution();
-        return num_fixed > 0;
-    }
-
-
     // below is the initialization functionality of lar_solver
 
     bool lar_solver::strategy_is_undecided() const {
@@ -1504,7 +1450,7 @@ namespace lp {
         m_mpq_lar_core_solver.m_r_x.resize(j + 1);
         m_mpq_lar_core_solver.m_r_lower_bounds.increase_size_by_one();
         m_mpq_lar_core_solver.m_r_upper_bounds.increase_size_by_one();
-        m_mpq_lar_core_solver.m_r_solver.inf_set_increase_size_by_one();
+        m_mpq_lar_core_solver.m_r_solver.inf_heap_increase_size_by_one();
         m_mpq_lar_core_solver.m_r_solver.m_costs.resize(j + 1);
         m_mpq_lar_core_solver.m_r_solver.m_d.resize(j + 1);
         lp_assert(m_mpq_lar_core_solver.m_r_heading.size() == j); // as A().column_count() on the entry to the method
@@ -1609,6 +1555,18 @@ namespace lp {
         return ret;
     }
 
+    /**
+    * \brief ensure there is a column index corresponding to vi
+    * If vi is already a column, just return vi
+    * If vi is for a term, then create a row that uses the term.
+    */
+    var_index lar_solver::ensure_column(var_index vi) {
+        if (lp::tv::is_term(vi))
+            return to_column(vi);
+        else
+            return vi;
+    }
+
 
     void lar_solver::add_row_from_term_no_constraint(const lar_term* term, unsigned term_ext_index) {
         TRACE("dump_terms", print_term(*term, tout) << std::endl;);
@@ -1780,13 +1738,20 @@ namespace lp {
         lconstraint_kind kind,
         const mpq& right_side,
         constraint_index constr_index) {
+        TRACE("lar_solver_feas", tout << "j = " << j << " was " << (this->column_is_feasible(j)?"feas":"non-feas") << std::endl;);
         m_constraints.activate(constr_index);
         if (column_has_upper_bound(j))
             update_column_type_and_bound_with_ub(j, kind, right_side, constr_index);
         else
             update_column_type_and_bound_with_no_ub(j, kind, right_side, constr_index);
+        TRACE("lar_solver_feas", tout << "j = " << j << " became " << (this->column_is_feasible(j)?"feas":"non-feas") << ", and " << (this->column_is_bounded(j)? "bounded":"non-bounded") << std::endl;);    
     }
-
+    // clang-format on
+    void lar_solver::insert_to_columns_with_changed_bounds(unsigned j) {
+        m_columns_with_changed_bounds.insert(j);
+        TRACE("lar_solver", tout << "column " << j << (column_is_feasible(j) ? " feas" : " non-feas") << "\n";);
+    }
+    // clang-format off
     void lar_solver::update_column_type_and_bound_check_on_equal(unsigned j,
         lconstraint_kind kind,
         const mpq& right_side,
@@ -1876,117 +1841,111 @@ namespace lp {
         }
     }
 
-
+    // clang-format on
     void lar_solver::update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
         lp_assert(column_has_lower_bound(j) && column_has_upper_bound(j));
         lp_assert(m_mpq_lar_core_solver.m_column_types[j] == column_type::boxed ||
-            m_mpq_lar_core_solver.m_column_types[j] == column_type::fixed);
+                  m_mpq_lar_core_solver.m_column_types[j] == column_type::fixed);
 
         mpq y_of_bound(0);
         switch (kind) {
-        case LT:
-            y_of_bound = -1;
-        case LE:
-        {
-            auto up = numeric_pair<mpq>(right_side, y_of_bound);
-            if (up < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
-                set_infeasible_column(j);
+            case LT:
+                y_of_bound = -1;
+            case LE: {
+                auto up = numeric_pair<mpq>(right_side, y_of_bound);
+                if (up < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
+                    set_infeasible_column(j);
+                }
+                if (up >= m_mpq_lar_core_solver.m_r_upper_bounds[j]) return;
+                m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
+                set_upper_bound_witness(j, ci);
+                insert_to_columns_with_changed_bounds(j);
+                break;
             }
-            if (up >= m_mpq_lar_core_solver.m_r_upper_bounds[j]) return;
-            m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
-            set_upper_bound_witness(j, ci);
-            insert_to_columns_with_changed_bounds(j);
-        }
-        break;
-        case GT:
-            y_of_bound = 1;
-        case GE:
-        {
-            auto low = numeric_pair<mpq>(right_side, y_of_bound);
-            if (low > m_mpq_lar_core_solver.m_r_upper_bounds[j]) {
-                set_infeasible_column(j);
+            case GT:
+                y_of_bound = 1;
+            case GE: {
+                auto low = numeric_pair<mpq>(right_side, y_of_bound);
+                if (low > m_mpq_lar_core_solver.m_r_upper_bounds[j]) {
+                    set_infeasible_column(j);
+                }
+                if (low < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
+                    return;
+                }
+                m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
+                set_lower_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_column_types[j] = (low == m_mpq_lar_core_solver.m_r_upper_bounds[j] ? column_type::fixed : column_type::boxed);
+                insert_to_columns_with_changed_bounds(j);
+                break;
+            } 
+            case EQ: {
+                auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
+                if (v > m_mpq_lar_core_solver.m_r_upper_bounds[j] || v < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
+                    set_infeasible_column(j);
+                }
+                set_upper_bound_witness(j, ci);
+                set_lower_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
+                break;
             }
-            if (low < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
-                return;
-            }
-            m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
-            insert_to_columns_with_changed_bounds(j);
-            set_lower_bound_witness(j, ci);
-            m_mpq_lar_core_solver.m_column_types[j] = (low == m_mpq_lar_core_solver.m_r_upper_bounds[j] ? column_type::fixed : column_type::boxed);
-        }
-        break;
-        case EQ:
-        {
-            auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
-            if (v > m_mpq_lar_core_solver.m_r_upper_bounds[j] || v < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
-                set_infeasible_column(j);
-            }
-            set_upper_bound_witness(j, ci);
-            set_lower_bound_witness(j, ci);
-            m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
-            break;
-        }
 
-        default:
-            UNREACHABLE();
+            default:
+                UNREACHABLE();
         }
         if (m_mpq_lar_core_solver.m_r_upper_bounds[j] == m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
             m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
         }
     }
+    // clang-format off
     void lar_solver::update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
         lp_assert(column_has_lower_bound(j) && !column_has_upper_bound(j));
         lp_assert(m_mpq_lar_core_solver.m_column_types[j] == column_type::lower_bound);
 
         mpq y_of_bound(0);
         switch (kind) {
-        case LT:
-            y_of_bound = -1;
-        case LE:
-        {
-            auto up = numeric_pair<mpq>(right_side, y_of_bound);
-            if (up < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
-                set_infeasible_column(j);
-            }
-            m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
-            set_upper_bound_witness(j, ci);
-            insert_to_columns_with_changed_bounds(j);
-            m_mpq_lar_core_solver.m_column_types[j] = (up == m_mpq_lar_core_solver.m_r_lower_bounds[j] ? column_type::fixed : column_type::boxed);
-        }
-        break;
-        case GT:
-            y_of_bound = 1;
-        case GE:
-        {
-            auto low = numeric_pair<mpq>(right_side, y_of_bound);
-            if (low < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
-                return;
-            }
-            m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
-            insert_to_columns_with_changed_bounds(j);
-            set_lower_bound_witness(j, ci);
-        }
-        break;
-        case EQ:
-        {
-            auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
-            if (v < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
-                set_infeasible_column(j);
+            case LT:
+                y_of_bound = -1;
+            case LE: {
+                auto up = numeric_pair<mpq>(right_side, y_of_bound);
+                if (up < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
+                    set_infeasible_column(j);
+                }
+                m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
+                set_upper_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_column_types[j] = (up == m_mpq_lar_core_solver.m_r_lower_bounds[j] ? column_type::fixed : column_type::boxed);
+                insert_to_columns_with_changed_bounds(j);
+                break;
+            } 
+            case GT:
+                y_of_bound = 1;
+            case GE: {
+                auto low = numeric_pair<mpq>(right_side, y_of_bound);
+                if (low < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
+                    return;
+                }
+                m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
+                set_lower_bound_witness(j, ci);
+                insert_to_columns_with_changed_bounds(j);
+                break;
+            } 
+            case EQ: {
+                auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
+                if (v < m_mpq_lar_core_solver.m_r_lower_bounds[j]) {
+                    set_infeasible_column(j);
+                }
+
+                set_upper_bound_witness(j, ci);
+                set_lower_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
+                m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
+                break;
             }
 
-            set_upper_bound_witness(j, ci);
-            set_lower_bound_witness(j, ci);
-            m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
-            m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
-            break;
+            default:
+                UNREACHABLE();
         }
-
-        default:
-            UNREACHABLE();
-        }
-
     }
-
+    // clang-format off
     void lar_solver::update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
         lp_assert(!column_has_lower_bound(j) && column_has_upper_bound(j));
         lp_assert(m_mpq_lar_core_solver.m_column_types[j] == column_type::upper_bound);
@@ -2012,9 +1971,10 @@ namespace lp {
                 set_infeasible_column(j);
             }
             m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
-            insert_to_columns_with_changed_bounds(j);
             set_lower_bound_witness(j, ci);
             m_mpq_lar_core_solver.m_column_types[j] = (low == m_mpq_lar_core_solver.m_r_upper_bounds[j] ? column_type::fixed : column_type::boxed);
+            insert_to_columns_with_changed_bounds(j);
+            
         }
         break;
         case EQ:
@@ -2035,48 +1995,44 @@ namespace lp {
             UNREACHABLE();
         }
     }
+    // clang-format on
     void lar_solver::update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index ci) {
         lp_assert(!column_has_lower_bound(j) && !column_has_upper_bound(j));
-        insert_to_columns_with_changed_bounds(j);
 
         mpq y_of_bound(0);
         switch (kind) {
-        case LT:
-            y_of_bound = -1;
-        case LE:
-        {
-            auto up = numeric_pair<mpq>(right_side, y_of_bound);
-            m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
-            set_upper_bound_witness(j, ci);
-            m_mpq_lar_core_solver.m_column_types[j] = column_type::upper_bound;
-        }
-        break;
-        case GT:
-            y_of_bound = 1;
-        case GE:
-        {
-            auto low = numeric_pair<mpq>(right_side, y_of_bound);
-            m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
-            insert_to_columns_with_changed_bounds(j);
-            set_lower_bound_witness(j, ci);
-            m_mpq_lar_core_solver.m_column_types[j] = column_type::lower_bound;
-        }
-        break;
-        case EQ:
-        {
-            auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
-            set_upper_bound_witness(j, ci);
-            set_lower_bound_witness(j, ci);
-            m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
-            m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
-            break;
-        }
+            case LT:
+                y_of_bound = -1;
+            case LE: {
+                auto up = numeric_pair<mpq>(right_side, y_of_bound);
+                m_mpq_lar_core_solver.m_r_upper_bounds[j] = up;
+                set_upper_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_column_types[j] = column_type::upper_bound;
+            } break;
+            case GT:
+                y_of_bound = 1;
+            case GE: {
+                auto low = numeric_pair<mpq>(right_side, y_of_bound);
+                m_mpq_lar_core_solver.m_r_lower_bounds[j] = low;
+                set_lower_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_column_types[j] = column_type::lower_bound;
 
-        default:
-            UNREACHABLE();
+            } break;
+            case EQ: {
+                auto v = numeric_pair<mpq>(right_side, zero_of_type<mpq>());
+                set_upper_bound_witness(j, ci);
+                set_lower_bound_witness(j, ci);
+                m_mpq_lar_core_solver.m_r_upper_bounds[j] = m_mpq_lar_core_solver.m_r_lower_bounds[j] = v;
+                m_mpq_lar_core_solver.m_column_types[j] = column_type::fixed;
+                break;
+            }
+
+            default:
+                UNREACHABLE();
         }
+        insert_to_columns_with_changed_bounds(j);
     }
-
+    // clang-format off
     bool lar_solver::column_corresponds_to_term(unsigned j) const {
         return tv::is_term(m_var_register.local_to_external(j));
     }
@@ -2180,8 +2136,8 @@ namespace lp {
     }
 
     bool lar_solver::get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq& rs, constraint_index& ci, bool& upper_bound) const {
-        lp_assert(t.is_term())
-            unsigned j;
+        lp_assert(t.is_term());
+        unsigned j;
         bool is_int;
         if (!m_var_register.external_is_used(t.index(), j, is_int))
             return false; // the term does not have a bound because it does not correspond to a column

src/math/lp/lar_solver.h

@@ -17,30 +17,32 @@
 
 
   --*/
+// clang-format off
 #pragma once
-#include "util/vector.h"
-#include <utility>
-#include "util/debug.h"
-#include "util/buffer.h"
+#include <algorithm>
+#include <functional>
+#include <stack>
+#include <string>
 #include <unordered_map>
 #include <unordered_set>
-#include <string>
-#include <algorithm>
-#include <stack>
-#include <functional>
+#include <utility>
+
+#include "math/lp/bound_analyzer_on_row.h"
+#include "math/lp/implied_bound.h"
+#include "math/lp/int_solver.h"
 #include "math/lp/lar_constraints.h"
 #include "math/lp/lar_core_solver.h"
-#include "math/lp/numeric_pair.h"
-#include "math/lp/lp_primal_core_solver.h"
-#include "math/lp/random_updater.h"
-#include "util/stacked_value.h"
-#include "math/lp/stacked_vector.h"
-#include "math/lp/implied_bound.h"
-#include "math/lp/bound_analyzer_on_row.h"
-#include "math/lp/int_solver.h"
-#include "math/lp/nra_solver.h"
-#include "math/lp/lp_types.h"
 #include "math/lp/lp_bound_propagator.h"
+#include "math/lp/lp_primal_core_solver.h"
+#include "math/lp/lp_types.h"
+#include "math/lp/nra_solver.h"
+#include "math/lp/numeric_pair.h"
+#include "math/lp/random_updater.h"
+#include "math/lp/stacked_vector.h"
+#include "util/buffer.h"
+#include "util/debug.h"
+#include "util/stacked_value.h"
+#include "util/vector.h"
 
 namespace lp {
 
@@ -48,10 +50,9 @@ class int_branch;
 class int_solver;
 class lar_solver : public column_namer {
     struct term_hasher {
-        std::size_t operator()(const lar_term &t) const
-        {            
-            using std::size_t;
+        std::size_t operator()(const lar_term& t) const {
             using std::hash;
+            using std::size_t;
             using std::string;
             size_t seed = 0;
             int i = 0;
@@ -66,110 +67,106 @@ class lar_solver : public column_namer {
     };
 
     struct term_comparer {
-        bool operator()(const lar_term &a, const lar_term& b) const
-        {
-            return a == b;            
+        bool operator()(const lar_term& a, const lar_term& b) const {
+            return a == b;
         }
     };
-    
+
     //////////////////// fields //////////////////////////
-    lp_settings                                         m_settings;
-    lp_status                                           m_status = lp_status::UNKNOWN;
-    stacked_value<simplex_strategy_enum>                m_simplex_strategy;
+    lp_settings m_settings;
+    lp_status m_status = lp_status::UNKNOWN;
+    stacked_value<simplex_strategy_enum> m_simplex_strategy;
     // such can be found at the initialization step: u < l
-    stacked_value<int>                                  m_crossed_bounds_column; 
-    lar_core_solver                                     m_mpq_lar_core_solver;
-    int_solver *                                        m_int_solver = nullptr;
-    bool                                                m_need_register_terms = false;
-    var_register                                        m_var_register;
-    var_register                                        m_term_register;
-    stacked_vector<ul_pair>                             m_columns_to_ul_pairs;
-    constraint_set                                      m_constraints;
+    stacked_value<int> m_crossed_bounds_column;
+    lar_core_solver m_mpq_lar_core_solver;
+    int_solver* m_int_solver = nullptr;
+    bool m_need_register_terms = false;
+    var_register m_var_register;
+    var_register m_term_register;
+    stacked_vector<ul_pair> m_columns_to_ul_pairs;
+    constraint_set m_constraints;
     // the set of column indices j such that bounds have changed for j
-    u_set                                               m_columns_with_changed_bounds;
-    u_set                                               m_rows_with_changed_bounds;
-    unsigned_vector                                     m_row_bounds_to_replay;
-    
-    u_set                                               m_basic_columns_with_changed_cost;
+    u_set m_columns_with_changed_bounds;
+    u_set m_rows_with_changed_bounds;
+    unsigned_vector m_row_bounds_to_replay;
+
+    u_set m_basic_columns_with_changed_cost;
     // these are basic columns with the value changed, so the corresponding row in the tableau
     // does not sum to zero anymore
-    u_set                                               m_incorrect_columns;
-    // copy of m_r_solver.inf_set()
-    unsigned_vector                                     m_inf_index_copy;
-    stacked_value<unsigned>                             m_term_count;
-    vector<lar_term*>                                   m_terms;
-    indexed_vector<mpq>                                 m_column_buffer;
+    u_set m_incorrect_columns;
+    // copy of m_r_solver.inf_heap()
+    unsigned_vector m_inf_index_copy;
+    stacked_value<unsigned> m_term_count;
+    vector<lar_term*> m_terms;
+    indexed_vector<mpq> m_column_buffer;
     std::unordered_map<lar_term, std::pair<mpq, unsigned>, term_hasher, term_comparer>
-    m_normalized_terms_to_columns;
-    vector<impq>                                        m_backup_x;
-    stacked_vector<unsigned>                            m_usage_in_terms;
+        m_normalized_terms_to_columns;
+    vector<impq> m_backup_x;
+    stacked_vector<unsigned> m_usage_in_terms;
     // ((x[j], is_int(j))->j)  for fixed j, used in equalities propagation
     // maps values to integral fixed vars
-    map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>>  m_fixed_var_table_int;
+    map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>> m_fixed_var_table_int;
     // maps values to non-integral fixed vars
-    map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>>  m_fixed_var_table_real;
+    map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>> m_fixed_var_table_real;
     // end of fields
 
     ////////////////// methods ////////////////////////////////
-    
-    static bool valid_index(unsigned j) { return static_cast<int>(j) >= 0;}
-    const lar_term & get_term(unsigned j) const;
+
+    static bool valid_index(unsigned j) { return static_cast<int>(j) >= 0; }
+    const lar_term& get_term(unsigned j) const;
     bool row_has_a_big_num(unsigned i) const;
     // init region
     bool strategy_is_undecided() const;
     void register_new_ext_var_index(unsigned ext_v, bool is_int);
-    bool term_is_int(const lar_term * t) const;
-    bool term_is_int(const vector<std::pair<mpq, unsigned int>> & coeffs) const;
+    bool term_is_int(const lar_term* t) const;
+    bool term_is_int(const vector<std::pair<mpq, unsigned int>>& coeffs) const;
     void add_non_basic_var_to_core_fields(unsigned ext_j, bool is_int);
     void add_new_var_to_core_fields_for_mpq(bool register_in_basis);
     mpq adjust_bound_for_int(lpvar j, lconstraint_kind&, const mpq&);
 
     // terms
-    bool all_vars_are_registered(const vector<std::pair<mpq, var_index>> & coeffs);
-    var_index add_term_undecided(const vector<std::pair<mpq, var_index>> & coeffs);
-    bool term_coeffs_are_ok(const vector<std::pair<mpq, var_index>> & coeffs);
+    bool all_vars_are_registered(const vector<std::pair<mpq, var_index>>& coeffs);
+    var_index add_term_undecided(const vector<std::pair<mpq, var_index>>& coeffs);
+    bool term_coeffs_are_ok(const vector<std::pair<mpq, var_index>>& coeffs);
     void push_term(lar_term* t);
-    void add_row_from_term_no_constraint(const lar_term * term, unsigned term_ext_index);
+    void add_row_from_term_no_constraint(const lar_term* term, unsigned term_ext_index);
     void add_basic_var_to_core_fields();
-    bool compare_values(impq const& lhs, lconstraint_kind k, const mpq & rhs);
+    bool compare_values(impq const& lhs, lconstraint_kind k, const mpq& rhs);
 
     inline void clear_columns_with_changed_bounds() { m_columns_with_changed_bounds.clear(); }
     inline void increase_by_one_columns_with_changed_bounds() { m_columns_with_changed_bounds.increase_size_by_one(); }
-    inline void insert_to_columns_with_changed_bounds(unsigned j) { m_columns_with_changed_bounds.insert(j); }
-    
-    void update_column_type_and_bound_check_on_equal(unsigned j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index, unsigned&);
-    void update_column_type_and_bound(unsigned j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
-    void update_column_type_and_bound_with_ub(var_index j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
-    void update_column_type_and_bound_with_no_ub(var_index j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
-    void update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
-    void update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
-    void update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
-    void update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq & right_side, constraint_index constr_index);
+    void insert_to_columns_with_changed_bounds(unsigned j);
+    void update_column_type_and_bound_check_on_equal(unsigned j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index, unsigned&);
+    void update_column_type_and_bound(unsigned j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
+    void update_column_type_and_bound_with_ub(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
+    void update_column_type_and_bound_with_no_ub(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
+    void update_bound_with_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
+    void update_bound_with_no_ub_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
+    void update_bound_with_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
+    void update_bound_with_no_ub_no_lb(var_index j, lconstraint_kind kind, const mpq& right_side, constraint_index constr_index);
     void register_in_fixed_var_table(unsigned, unsigned&);
     void remove_non_fixed_from_fixed_var_table();
-    constraint_index add_var_bound_on_constraint_for_term(var_index j, lconstraint_kind kind, const mpq & right_side);
+    constraint_index add_var_bound_on_constraint_for_term(var_index j, lconstraint_kind kind, const mpq& right_side);
     inline void set_infeasible_column(unsigned j) {
         set_status(lp_status::INFEASIBLE);
         m_crossed_bounds_column = j;
     }
     constraint_index add_constraint_from_term_and_create_new_column_row(unsigned term_j, const lar_term* term,
-                                                                        lconstraint_kind kind, const mpq & right_side);
+                                                                        lconstraint_kind kind, const mpq& right_side);
     unsigned row_of_basic_column(unsigned) const;
     void decide_on_strategy_and_adjust_initial_state();
     void adjust_initial_state();
     void adjust_initial_state_for_tableau_rows();
     bool sizes_are_correct() const;
-    bool implied_bound_is_correctly_explained(implied_bound const & be, const vector<std::pair<mpq, unsigned>> & explanation) const;
-    
+    bool implied_bound_is_correctly_explained(implied_bound const& be, const vector<std::pair<mpq, unsigned>>& explanation) const;
 
     void substitute_basis_var_in_terms_for_row(unsigned i);
-    
+
     template <typename T>
-    unsigned calculate_implied_bounds_for_row(unsigned row_index, lp_bound_propagator<T> & bp) {
-        
-        if (A_r().m_rows[row_index].size() > settings().max_row_length_for_bound_propagation  || row_has_a_big_num(row_index))
+    unsigned calculate_implied_bounds_for_row(unsigned row_index, lp_bound_propagator<T>& bp) {
+        if (A_r().m_rows[row_index].size() > settings().max_row_length_for_bound_propagation || row_has_a_big_num(row_index))
             return 0;
-        
+
         return bound_analyzer_on_row<row_strip<mpq>, lp_bound_propagator<T>>::analyze_row(
             A_r().m_rows[row_index],
             null_ci,
@@ -177,27 +174,28 @@ class lar_solver : public column_namer {
             row_index,
             bp);
     }
-    
+
+    static void clean_popped_elements_for_heap(unsigned n, lpvar_heap& set);
     static void clean_popped_elements(unsigned n, u_set& set);
-    bool maximize_term_on_tableau(const lar_term & term,
-                                  impq &term_max);
+    bool maximize_term_on_tableau(const lar_term& term,
+                                  impq& term_max);
     bool costs_are_zeros_for_r_solver() const;
     bool reduced_costs_are_zeroes_for_r_solver() const;
-    void set_costs_to_zero(const lar_term & term);
-    void prepare_costs_for_r_solver(const lar_term & term);
-    bool maximize_term_on_corrected_r_solver(lar_term & term, impq &term_max);    
+    void set_costs_to_zero(const lar_term& term);
+    void prepare_costs_for_r_solver(const lar_term& term);
+    bool maximize_term_on_corrected_r_solver(lar_term& term, impq& term_max);
     void pop_core_solver_params();
     void pop_core_solver_params(unsigned k);
     void set_upper_bound_witness(var_index j, constraint_index ci);
     void set_lower_bound_witness(var_index j, constraint_index ci);
-    void substitute_terms_in_linear_expression( const vector<std::pair<mpq, var_index>>& left_side_with_terms,
-                                                vector<std::pair<mpq, var_index>> &left_side) const;
-    
+    void substitute_terms_in_linear_expression(const vector<std::pair<mpq, var_index>>& left_side_with_terms,
+                                               vector<std::pair<mpq, var_index>>& left_side) const;
+
     void detect_rows_of_bound_change_column_for_nbasic_column_tableau(unsigned j);
     bool use_tableau_costs() const;
     bool tableau_with_costs() const;
     bool costs_are_used() const;
-    void change_basic_columns_dependend_on_a_given_nb_column(unsigned j, const numeric_pair<mpq> & delta);
+    void change_basic_columns_dependend_on_a_given_nb_column(unsigned j, const numeric_pair<mpq>& delta);
     void update_x_and_inf_costs_for_column_with_changed_bounds(unsigned j);
     unsigned num_changed_bounds() const { return m_rows_with_changed_bounds.size(); }
     void insert_row_with_changed_bounds(unsigned rid);
@@ -209,19 +207,19 @@ class lar_solver : public column_namer {
     numeric_pair<mpq> get_basic_var_value_from_row(unsigned i);
     bool all_constrained_variables_are_registered(const vector<std::pair<mpq, var_index>>& left_side);
     bool all_constraints_hold() const;
-    bool constraint_holds(const lar_base_constraint & constr, std::unordered_map<var_index, mpq> & var_map) const;
-    static void register_in_map(std::unordered_map<var_index, mpq> & coeffs, const lar_base_constraint & cn, const mpq & a);
-    static void register_monoid_in_map(std::unordered_map<var_index, mpq> & coeffs, const mpq & a, unsigned j);
-    bool the_left_sides_sum_to_zero(const vector<std::pair<mpq, unsigned>> & evidence) const;
+    bool constraint_holds(const lar_base_constraint& constr, std::unordered_map<var_index, mpq>& var_map) const;
+    static void register_in_map(std::unordered_map<var_index, mpq>& coeffs, const lar_base_constraint& cn, const mpq& a);
+    static void register_monoid_in_map(std::unordered_map<var_index, mpq>& coeffs, const mpq& a, unsigned j);
+    bool the_left_sides_sum_to_zero(const vector<std::pair<mpq, unsigned>>& evidence) const;
     bool explanation_is_correct(explanation&) const;
     bool inf_explanation_is_correct() const;
-    mpq sum_of_right_sides_of_explanation(explanation &) const;
+    mpq sum_of_right_sides_of_explanation(explanation&) const;
     void get_infeasibility_explanation_for_inf_sign(
-        explanation & exp,
-        const vector<std::pair<mpq, unsigned>> & inf_row,
+        explanation& exp,
+        const vector<std::pair<mpq, unsigned>>& inf_row,
         int inf_sign) const;
-    mpq get_left_side_val(const lar_base_constraint &  cns, const std::unordered_map<var_index, mpq> & var_map) const;
-    void fill_var_set_for_random_update(unsigned sz, var_index const * vars, vector<unsigned>& column_list);
+    mpq get_left_side_val(const lar_base_constraint& cns, const std::unordered_map<var_index, mpq>& var_map) const;
+    void fill_var_set_for_random_update(unsigned sz, var_index const* vars, vector<unsigned>& column_list);
     bool column_represents_row_in_tableau(unsigned j);
     void make_sure_that_the_bottom_right_elem_not_zero_in_tableau(unsigned i, unsigned j);
     void remove_last_row_and_column_from_tableau(unsigned j);
@@ -230,25 +228,25 @@ class lar_solver : public column_namer {
     void remove_last_column_from_basis_tableau(unsigned j);
     void remove_last_column_from_tableau();
     void pop_tableau();
-    void clean_inf_set_of_r_solver_after_pop();
+    void clean_inf_heap_of_r_solver_after_pop();
     inline bool column_value_is_integer(unsigned j) const { return get_column_value(j).is_int(); }
     bool model_is_int_feasible() const;
-    
-    bool bound_is_integer_for_integer_column(unsigned j, const mpq & right_side) const;
-    inline lar_core_solver & get_core_solver() { return m_mpq_lar_core_solver; }
+
+    bool bound_is_integer_for_integer_column(unsigned j, const mpq& right_side) const;
+    inline lar_core_solver& get_core_solver() { return m_mpq_lar_core_solver; }
     var_index to_column(unsigned ext_j) const;
     void fix_terms_with_rounded_columns();
     bool remove_from_basis(unsigned);
     lar_term get_term_to_maximize(unsigned ext_j) const;
-    bool sum_first_coords(const lar_term& t, mpq & val) const;
+    bool sum_first_coords(const lar_term& t, mpq& val) const;
     void register_normalized_term(const lar_term&, lpvar);
     void deregister_normalized_term(const lar_term&);
 
     mutable std::unordered_set<impq> m_set_of_different_pairs;
-    mutable std::unordered_set<mpq>  m_set_of_different_singles;
+    mutable std::unordered_set<mpq> m_set_of_different_singles;
     mutable mpq m_delta;
 
-public:
+   public:
     // this function just looks at the status
     bool is_feasible() const;
 
@@ -256,7 +254,6 @@ public:
         return m_fixed_var_table_int;
     }
 
-    
     const map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>>& fixed_var_table_real() const {
         return m_fixed_var_table_real;
     }
@@ -265,11 +262,12 @@ public:
         return m_fixed_var_table_real;
     }
 
-    bool find_in_fixed_tables(const rational& mpq, bool is_int, unsigned & j) const {
-        return is_int? fixed_var_table_int().find(mpq, j) : fixed_var_table_real().find(mpq, j);
+    bool find_in_fixed_tables(const rational& mpq, bool is_int, unsigned& j) const {
+        return is_int ? fixed_var_table_int().find(mpq, j) : fixed_var_table_real().find(mpq, j);
     }
-    
-    template <typename T> void remove_non_fixed_from_table(T&);
+
+    template <typename T>
+    void remove_non_fixed_from_table(T&);
 
     unsigned external_to_column_index(unsigned) const;
 
@@ -282,38 +280,38 @@ public:
     inline void set_column_value_test(unsigned j, const impq& v) {
         set_column_value(j, v);
     }
-    
+
     var_index add_named_var(unsigned ext_j, bool is_integer, const std::string&);
 
-    lp_status maximize_term(unsigned j_or_term, impq &term_max);
+    lp_status maximize_term(unsigned j_or_term, impq& term_max);
 
-    inline core_solver_pretty_printer<lp::mpq, lp::impq> pp(std::ostream& out) const { 
-        return core_solver_pretty_printer<lp::mpq, lp::impq>(m_mpq_lar_core_solver.m_r_solver, out); 
+    inline core_solver_pretty_printer<lp::mpq, lp::impq> pp(std::ostream& out) const {
+        return core_solver_pretty_printer<lp::mpq, lp::impq>(m_mpq_lar_core_solver.m_r_solver, out);
     }
 
-    void get_infeasibility_explanation(explanation &) const;
+    void get_infeasibility_explanation(explanation&) const;
 
     inline void backup_x() { m_backup_x = m_mpq_lar_core_solver.m_r_x; }
 
     inline void restore_x() { m_mpq_lar_core_solver.m_r_x = m_backup_x; }
 
     template <typename T>
-    void explain_implied_bound(const implied_bound & ib, lp_bound_propagator<T> & bp) {
+    void explain_implied_bound(const implied_bound& ib, lp_bound_propagator<T>& bp) {
         unsigned i = ib.m_row_or_term_index;
         int bound_sign = (ib.m_is_lower_bound ? 1 : -1);
         int j_sign = (ib.m_coeff_before_j_is_pos ? 1 : -1) * bound_sign;
         unsigned bound_j = ib.m_j;
-        if (tv::is_term(bound_j)) 
+        if (tv::is_term(bound_j))
             bound_j = m_var_register.external_to_local(bound_j);
 
         for (auto const& r : get_row(i)) {
             unsigned j = r.var();
-            if (j == bound_j) 
+            if (j == bound_j)
                 continue;
             mpq const& a = r.coeff();
             int a_sign = is_pos(a) ? 1 : -1;
             int sign = j_sign * a_sign;
-            const ul_pair & ul =  m_columns_to_ul_pairs[j];
+            const ul_pair& ul = m_columns_to_ul_pairs[j];
             auto witness = sign > 0 ? ul.upper_bound_witness() : ul.lower_bound_witness();
             lp_assert(is_valid(witness));
             bp.consume(a, witness);
@@ -327,13 +325,13 @@ public:
     }
 
     // lp_assert(implied_bound_is_correctly_explained(ib, explanation)); }
-    constraint_index mk_var_bound(var_index j, lconstraint_kind kind, const mpq & right_side);
+    constraint_index mk_var_bound(var_index j, lconstraint_kind kind, const mpq& right_side);
     void activate_check_on_equal(constraint_index, var_index&);
     void activate(constraint_index);
-    void random_update(unsigned sz, var_index const * vars);
+    void random_update(unsigned sz, var_index const* vars);
     void mark_rows_for_bound_prop(lpvar j);
     template <typename T>
-    void propagate_bounds_for_touched_rows(lp_bound_propagator<T> & bp) {
+    void propagate_bounds_for_touched_rows(lp_bound_propagator<T>& bp) {
         unsigned num_prop = 0;
         for (unsigned i : m_rows_with_changed_bounds) {
             num_prop += calculate_implied_bounds_for_row(i, bp);
@@ -347,7 +345,7 @@ public:
             bp.clear_for_eq();
             for (unsigned i : m_rows_with_changed_bounds) {
                 unsigned offset_eqs = stats().m_offset_eqs;
-                bp.cheap_eq_tree(i);                
+                bp.cheap_eq_tree(i);
                 if (settings().get_cancel_flag())
                     return;
                 if (stats().m_offset_eqs > offset_eqs)
@@ -358,74 +356,74 @@ public:
     }
 
     template <typename T>
-    void check_missed_propagations(lp_bound_propagator<T> & bp) {
-        for (unsigned i = 0; i < A_r().row_count(); i++) 
-            if (!m_rows_with_changed_bounds.contains(i))                
+    void check_missed_propagations(lp_bound_propagator<T>& bp) {
+        for (unsigned i = 0; i < A_r().row_count(); i++)
+            if (!m_rows_with_changed_bounds.contains(i))
                 if (0 < calculate_implied_bounds_for_row(i, bp)) {
                     verbose_stream() << i << ": " << get_row(i) << "\n";
                 }
     }
 
-    bool is_fixed_at_bound(column_index const& j);
-    bool has_fixed_at_bound();
-    
-    bool is_fixed(column_index const& j) const { return column_is_fixed(j); }    
+    bool is_fixed(column_index const& j) const { return column_is_fixed(j); }
     inline column_index to_column_index(unsigned v) const { return column_index(external_to_column_index(v)); }
     bool external_is_used(unsigned) const;
     void pop(unsigned k);
-    bool compare_values(var_index j, lconstraint_kind kind, const mpq & right_side);
-    var_index add_term(const vector<std::pair<mpq, var_index>> & coeffs, unsigned ext_i);
+    bool compare_values(var_index j, lconstraint_kind kind, const mpq& right_side);
+    var_index add_term(const vector<std::pair<mpq, var_index>>& coeffs, unsigned ext_i);
     void register_existing_terms();
-    constraint_index add_var_bound(var_index, lconstraint_kind, const mpq &);
-    constraint_index add_var_bound_check_on_equal(var_index, lconstraint_kind, const mpq &, var_index&);
-    
+    var_index ensure_column(var_index vi);
+    constraint_index add_var_bound(var_index, lconstraint_kind, const mpq&);
+    constraint_index add_var_bound_check_on_equal(var_index, lconstraint_kind, const mpq&, var_index&);
+
     var_index add_var(unsigned ext_j, bool is_integer);
     void set_cut_strategy(unsigned cut_frequency);
     inline unsigned column_count() const { return A_r().column_count(); }
     inline var_index local_to_external(var_index idx) const {
-        return tv::is_term(idx)?
-            m_term_register.local_to_external(idx) : m_var_register.local_to_external(idx);
+        return tv::is_term(idx) ? m_term_register.local_to_external(idx) : m_var_register.local_to_external(idx);
     }
     bool column_corresponds_to_term(unsigned) const;
+    const lar_term& column_to_term(unsigned j) const {
+        SASSERT(column_corresponds_to_term(j));
+        return get_term(column2tv(to_column_index(j)));
+    }
+
     inline unsigned row_count() const { return A_r().row_count(); }
     bool var_is_registered(var_index vj) const;
-    void clear_inf_set() {
-        m_mpq_lar_core_solver.m_r_solver.inf_set().clear();        
+    void clear_inf_heap() {
+        m_mpq_lar_core_solver.m_r_solver.inf_heap().clear();
     }
-    inline void remove_column_from_inf_set(unsigned j) {
-        m_mpq_lar_core_solver.m_r_solver.remove_column_from_inf_set(j);
+
+    void pivot(int entering, int leaving) {
+        m_mpq_lar_core_solver.pivot(entering, leaving);
     }
+
     template <typename ChangeReport>
     void change_basic_columns_dependend_on_a_given_nb_column_report(unsigned j,
-                                                                    const numeric_pair<mpq> & delta,
+                                                                    const numeric_pair<mpq>& delta,
                                                                     const ChangeReport& after) {
-        
-      for (const auto & c : A_r().m_columns[j]) {
-                unsigned bj = m_mpq_lar_core_solver.m_r_basis[c.var()];
-                if (tableau_with_costs()) {
-                    m_basic_columns_with_changed_cost.insert(bj);
-                }
-                m_mpq_lar_core_solver.m_r_solver.add_delta_to_x_and_track_feasibility(bj, - A_r().get_val(c) * delta);
-                after(bj);
-                TRACE("change_x_del",
-                      tout << "changed basis column " << bj << ", it is " <<
-                      ( m_mpq_lar_core_solver.m_r_solver.column_is_feasible(bj)?  "feas":"inf") << std::endl;);
-            }
-   }
+        for (const auto& c : A_r().m_columns[j]) {
+            unsigned bj = m_mpq_lar_core_solver.m_r_basis[c.var()];
+            if (tableau_with_costs())
+                m_basic_columns_with_changed_cost.insert(bj);
+            m_mpq_lar_core_solver.m_r_solver.add_delta_to_x_and_track_feasibility(bj, -A_r().get_val(c) * delta);
+            after(bj);
+            TRACE("change_x_del",
+                  tout << "changed basis column " << bj << ", it is " << (m_mpq_lar_core_solver.m_r_solver.column_is_feasible(bj) ? "feas" : "inf") << std::endl;);
+        }
+    }
 
     template <typename ChangeReport>
     void set_value_for_nbasic_column_report(unsigned j,
-                                            const impq & new_val,
+                                            const impq& new_val,
                                             const ChangeReport& after) {
-
         lp_assert(!is_base(j));
-        auto & x = m_mpq_lar_core_solver.m_r_x[j];        
+        auto& x = m_mpq_lar_core_solver.m_r_x[j];
         auto delta = new_val - x;
         x = new_val;
         after(j);
         change_basic_columns_dependend_on_a_given_nb_column_report(j, delta, after);
     }
-    
+
     template <typename Blocker, typename ChangeReport>
     bool try_to_patch(lpvar j, const mpq& val,
                       const Blocker& is_blocked,
@@ -442,8 +440,8 @@ public:
         impq delta = get_column_value(j) - ival;
         for (auto c : A_r().column(j)) {
             unsigned row_index = c.var();
-            const mpq & a = c.coeff();        
-            unsigned rj = m_mpq_lar_core_solver.m_r_basis[row_index];      
+            const mpq& a = c.coeff();
+            unsigned rj = m_mpq_lar_core_solver.m_r_basis[row_index];
             impq rj_new_val = a * delta + get_column_value(rj);
             // if (column_is_int(rj) && !rj_new_val.is_int())
             //     return false;
@@ -460,65 +458,63 @@ public:
 
     inline bool column_has_lower_bound(unsigned j) const {
         return m_mpq_lar_core_solver.m_r_solver.column_has_lower_bound(j);
-    }    
+    }
 
-    inline
-    constraint_index get_column_upper_bound_witness(unsigned j) const {
+    inline constraint_index get_column_upper_bound_witness(unsigned j) const {
         if (tv::is_term(j)) {
             j = m_var_register.external_to_local(j);
         }
         return m_columns_to_ul_pairs()[j].upper_bound_witness();
     }
 
-    inline
-    const impq& get_upper_bound(column_index j) const {
+    inline const impq& get_upper_bound(column_index j) const {
         return m_mpq_lar_core_solver.m_r_solver.m_upper_bounds[j];
     }
 
-    inline
-    const impq& get_lower_bound(column_index j) const {
+    inline const impq& get_lower_bound(column_index j) const {
         return m_mpq_lar_core_solver.m_r_solver.m_lower_bounds[j];
     }
-    bool has_lower_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const;    
+    bool has_lower_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const;
     bool has_upper_bound(var_index var, constraint_index& ci, mpq& value, bool& is_strict) const;
     bool has_value(var_index var, mpq& value) const;
-    bool fetch_normalized_term_column(const lar_term& t, std::pair<mpq, lpvar>& ) const;
+    bool fetch_normalized_term_column(const lar_term& t, std::pair<mpq, lpvar>&) const;
     unsigned map_term_index_to_column_index(unsigned j) const;
     bool column_is_fixed(unsigned j) const;
     bool column_is_free(unsigned j) const;
+    bool column_is_feasible(unsigned j) const { return m_mpq_lar_core_solver.m_r_solver.column_is_feasible(j);}
     unsigned column_to_reported_index(unsigned j) const;
-    lp_settings & settings();
-    lp_settings const & settings() const;
+    lp_settings& settings();
+    lp_settings const& settings() const;
     statistics& stats();
- 
+
     void updt_params(params_ref const& p);
     column_type get_column_type(unsigned j) const { return m_mpq_lar_core_solver.m_column_types()[j]; }
-    const impq & get_lower_bound(unsigned j) const { return m_mpq_lar_core_solver.m_r_lower_bounds()[j]; }
-    const impq & get_upper_bound(unsigned j) const { return m_mpq_lar_core_solver.m_r_upper_bounds()[j]; }
+    const impq& get_lower_bound(unsigned j) const { return m_mpq_lar_core_solver.m_r_lower_bounds()[j]; }
+    const impq& get_upper_bound(unsigned j) const { return m_mpq_lar_core_solver.m_r_upper_bounds()[j]; }
     std::ostream& print_terms(std::ostream& out) const;
-    std::ostream& print_term(lar_term const& term, std::ostream & out) const;
-    static std::ostream& print_term_as_indices(lar_term const& term, std::ostream & out);
-    std::ostream& print_constraint_indices_only(const lar_base_constraint * c, std::ostream & out) const;
-    std::ostream& print_implied_bound(const implied_bound& be, std::ostream & out) const;
+    std::ostream& print_term(lar_term const& term, std::ostream& out) const;
+    static std::ostream& print_term_as_indices(lar_term const& term, std::ostream& out);
+    std::ostream& print_constraint_indices_only(const lar_base_constraint* c, std::ostream& out) const;
+    std::ostream& print_implied_bound(const implied_bound& be, std::ostream& out) const;
     std::ostream& print_values(std::ostream& out) const;
     std::ostream& display(std::ostream& out) const;
 
     bool init_model() const;
     mpq get_value(column_index const& j) const;
     mpq get_tv_value(tv const& t) const;
-    const impq & get_tv_ivalue(tv const& t) const;
-    void get_model(std::unordered_map<var_index, mpq> & variable_values) const;
+    const impq& get_tv_ivalue(tv const& t) const;
+    void get_model(std::unordered_map<var_index, mpq>& variable_values) const;
     void get_rid_of_inf_eps();
-    void get_model_do_not_care_about_diff_vars(std::unordered_map<var_index, mpq> & variable_values) const;
+    void get_model_do_not_care_about_diff_vars(std::unordered_map<var_index, mpq>& variable_values) const;
     std::string get_variable_name(var_index vi) const override;
     void set_variable_name(var_index vi, std::string);
     inline unsigned number_of_vars() const { return m_var_register.size(); }
     inline bool is_base(unsigned j) const { return m_mpq_lar_core_solver.m_r_heading[j] >= 0; }
-    inline const impq & column_lower_bound(unsigned j) const {
+    inline const impq& column_lower_bound(unsigned j) const {
         return m_mpq_lar_core_solver.lower_bound(j);
     }
 
-    inline const impq & column_upper_bound(unsigned j) const {
+    inline const impq& column_upper_bound(unsigned j) const {
         return m_mpq_lar_core_solver.upper_bound(j);
     }
 
@@ -531,9 +527,9 @@ public:
     }
 
     std::pair<constraint_index, constraint_index> add_equality(lpvar j, lpvar k);
-    
-    inline void get_bound_constraint_witnesses_for_column(unsigned j, constraint_index & lc, constraint_index & uc) const {
-        const ul_pair & ul = m_columns_to_ul_pairs[j];
+
+    inline void get_bound_constraint_witnesses_for_column(unsigned j, constraint_index& lc, constraint_index& uc) const {
+        const ul_pair& ul = m_columns_to_ul_pairs[j];
         lc = ul.lower_bound_witness();
         uc = ul.upper_bound_witness();
     }
@@ -550,19 +546,19 @@ public:
     inline tv column2tv(column_index const& c) const {
         return tv::raw(column_to_reported_index(c));
     }
-    
+
     inline std::ostream& print_column_info(unsigned j, std::ostream& out) const {
         m_mpq_lar_core_solver.m_r_solver.print_column_info(j, out);
         if (tv::is_term(j)) {
             print_term_as_indices(get_term(j), out) << "\n";
-            
-        } else if (column_corresponds_to_term(j)) { 
+
+        } else if (column_corresponds_to_term(j)) {
             const lar_term& t = get_term(m_var_register.local_to_external(j));
             print_term_as_indices(t, out) << "\n";
         }
         return out;
     }
-    
+
     void subst_known_terms(lar_term*);
 
     inline std::ostream& print_column_bound_info(unsigned j, std::ostream& out) const {
@@ -573,67 +569,68 @@ public:
 
     inline bool has_inf_int() const {
         for (unsigned j = 0; j < column_count(); j++) {
-            if (column_is_int(j) && ! column_value_is_int(j))
+            if (column_is_int(j) && !column_value_is_int(j))
                 return true;
         }
         return false;
     }
-    inline const vector<lar_term*> & terms() const { return m_terms; }
+    inline const vector<lar_term*>& terms() const { return m_terms; }
     inline lar_term const& term(unsigned i) const { return *m_terms[i]; }
-    inline void set_int_solver(int_solver * int_slv) { m_int_solver = int_slv; }
-    inline int_solver * get_int_solver() { return m_int_solver; }
-    inline const int_solver * get_int_solver() const { return m_int_solver; }
-    inline const lar_term & get_term(tv const& t) const { lp_assert(t.is_term()); return *m_terms[t.id()]; }
-    lp_status find_feasible_solution();   
+    inline void set_int_solver(int_solver* int_slv) { m_int_solver = int_slv; }
+    inline int_solver* get_int_solver() { return m_int_solver; }
+    inline const int_solver* get_int_solver() const { return m_int_solver; }
+    inline const lar_term& get_term(tv const& t) const {
+        lp_assert(t.is_term());
+        return *m_terms[t.id()];
+    }
+    lp_status find_feasible_solution();
     void move_non_basic_columns_to_bounds(bool);
     bool move_non_basic_column_to_bounds(unsigned j, bool);
     inline bool r_basis_has_inf_int() const {
         for (unsigned j : r_basis()) {
-            if (column_is_int(j) && ! column_value_is_int(j))
+            if (column_is_int(j) && !column_value_is_int(j))
                 return true;
         }
         return false;
     }
     void round_to_integer_solution();
-    inline const row_strip<mpq> &  get_row(unsigned i) const { return A_r().m_rows[i]; }
-    inline const row_strip<mpq> &  basic2row(unsigned i) const { return A_r().m_rows[row_of_basic_column(i)]; }
-    inline const column_strip &  get_column(unsigned i) const { return A_r().m_columns[i]; }
+    inline const row_strip<mpq>& get_row(unsigned i) const { return A_r().m_rows[i]; }
+    inline const row_strip<mpq>& basic2row(unsigned i) const { return A_r().m_rows[row_of_basic_column(i)]; }
+    inline const column_strip& get_column(unsigned i) const { return A_r().m_columns[i]; }
     bool row_is_correct(unsigned i) const;
     bool ax_is_correct() const;
-    bool get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq &rs, constraint_index& ci, bool &upper_bound) const;
+    bool get_equality_and_right_side_for_term_on_current_x(tv const& t, mpq& rs, constraint_index& ci, bool& upper_bound) const;
     bool var_is_int(var_index v) const;
-    inline const vector<int> & r_heading() const { return m_mpq_lar_core_solver.m_r_heading; }
-    inline const vector<unsigned> & r_basis() const { return m_mpq_lar_core_solver.r_basis(); }
-    inline const vector<unsigned> & r_nbasis() const { return m_mpq_lar_core_solver.r_nbasis(); }
-    inline bool column_is_real(unsigned j) const { return !column_is_int(j); }	
+    inline const vector<int>& r_heading() const { return m_mpq_lar_core_solver.m_r_heading; }
+    inline const vector<unsigned>& r_basis() const { return m_mpq_lar_core_solver.r_basis(); }
+    inline const vector<unsigned>& r_nbasis() const { return m_mpq_lar_core_solver.r_nbasis(); }
+    inline bool column_is_real(unsigned j) const { return !column_is_int(j); }
     lp_status get_status() const;
-    bool has_changed_columns() const { return !m_columns_with_changed_bounds.empty();  }
+    bool has_changed_columns() const { return !m_columns_with_changed_bounds.empty(); }
     void set_status(lp_status s);
     lp_status solve();
-    void fill_explanation_from_crossed_bounds_column(explanation & evidence) const;
+    void fill_explanation_from_crossed_bounds_column(explanation& evidence) const;
     bool term_is_used_as_row(unsigned term) const;
     bool tighten_term_bounds_by_delta(tv const& t, const impq&);
     lar_solver();
     void set_track_pivoted_rows(bool v);
-    bool get_track_pivoted_rows() const;    
+    bool get_track_pivoted_rows() const;
     ~lar_solver() override;
     const vector<impq>& r_x() const { return m_mpq_lar_core_solver.m_r_x; }
     bool column_is_int(unsigned j) const;
     inline bool column_value_is_int(unsigned j) const { return m_mpq_lar_core_solver.m_r_x[j].is_int(); }
-    inline static_matrix<mpq, impq> & A_r() { return m_mpq_lar_core_solver.m_r_A; }
-    inline const static_matrix<mpq, impq> & A_r() const { return m_mpq_lar_core_solver.m_r_A; }
+    inline static_matrix<mpq, impq>& A_r() { return m_mpq_lar_core_solver.m_r_A; }
+    inline const static_matrix<mpq, impq>& A_r() const { return m_mpq_lar_core_solver.m_r_A; }
     // columns
     bool column_is_int(column_index const& j) const { return column_is_int((unsigned)j); }
-//    const impq& get_ivalue(column_index const& j) const { return get_column_value(j); }
+    //    const impq& get_ivalue(column_index const& j) const { return get_column_value(j); }
     const impq& get_column_value(column_index const& j) const { return m_mpq_lar_core_solver.m_r_x[j]; }
-    inline
-    var_index external_to_local(unsigned j) const {
+    inline var_index external_to_local(unsigned j) const {
         var_index local_j;
         if (m_var_register.external_is_used(j, local_j) ||
             m_term_register.external_is_used(j, local_j)) {
             return local_j;
-        }
-        else {
+        } else {
             return -1;
         }
     }
@@ -644,6 +641,5 @@ public:
     }
     friend int_solver;
     friend int_branch;
-    
 };
-}
+}  // namespace lp

src/math/lp/lar_term.h

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 #include "math/lp/indexed_vector.h"
 #include "util/map.h"

src/math/lp/lia_move.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 namespace lp {
 enum class lia_move {

src/math/lp/lp_api.h

@@ -7,6 +7,7 @@ Author:
     Nikolaj Bjorner (nbjorner)
 
 --*/
+// clang-format off
 #pragma once
 
 #include "util/inf_rational.h"

src/math/lp/lp_bound_propagator.h

@@ -4,39 +4,40 @@
   Nikolaj Bjorner (nbjorner)
   Lev Nachmanson   (levnach)
 */
+//clang-format off
 #pragma once
-#include "math/lp/lp_settings.h"
 #include <utility>
+
+#include "math/lp/lp_settings.h"
+#include "util/uint_set.h"
 namespace lp {
 template <typename T>
 class lp_bound_propagator {
-    class edge; // forward definition
+    class edge;  // forward definition
     // vertex represents a column
-    // The set of vertices is organised in a tree.
+    // The set of vertices is organized in a tree.
     // The edges of the tree are rows,
     // Vertices with m_neg set to false grow with the same rate as the root.
     // Vertices with m_neq set to true diminish with the same rate as the roow grows.
     // When two vertices with the same m_neg have the same value of columns
-    // then we have an equality betweet the columns.
-    class vertex {        
-        unsigned           m_column;
-        vector<edge>       m_edges;
-        edge               m_edge_from_parent;
-        unsigned           m_level; // the distance in hops to the root;
-                                    // it is handy to find the common ancestor
-    public:
+    // then we have an equality between the columns.
+    class vertex {
+        unsigned m_column;
+        vector<edge> m_edges;
+        edge m_edge_from_parent;
+        unsigned m_level;  // the distance in hops to the root;
+                           // it is handy to find the common ancestor
+       public:
         vertex() {}
-        vertex(unsigned column) :
-            m_column(column),
-            m_level(0)
-        {}
+        vertex(unsigned column) : m_column(column),
+                                  m_level(0) {}
         unsigned column() const { return m_column; }
         const vertex* parent() const { return m_edge_from_parent.source(); }
         vertex* parent() { return m_edge_from_parent.source(); }
         unsigned level() const { return m_level; }
-        void set_edge_from_parent(edge &e) { m_edge_from_parent = e; }
+        void set_edge_from_parent(edge& e) { m_edge_from_parent = e; }
         const edge& edge_from_parent() const { return m_edge_from_parent; }
-        
+
         void add_child(int row, vertex* child) {
             SASSERT(*this != *child);
             SASSERT(child->parent() == nullptr);
@@ -45,7 +46,7 @@ class lp_bound_propagator {
             child->set_edge_from_parent(e);
             child->m_level = m_level + 1;
         }
-        const vector<edge> & edges() const { return m_edges; }
+        const vector<edge>& edges() const { return m_edges; }
         bool operator==(const vertex& o) const {
             return m_column == o.m_column;
         }
@@ -58,7 +59,8 @@ class lp_bound_propagator {
         vertex* m_source;
         vertex* m_target;
         int m_row;
-    public:
+
+       public:
         edge(vertex* source, vertex* target, int row) : m_source(source), m_target(target), m_row(row) {}
         edge() : m_source(nullptr), m_target(nullptr), m_row(-1) {}
         const vertex* source() const { return m_source; }
@@ -69,57 +71,58 @@ class lp_bound_propagator {
         edge reverse() const { return edge(m_target, m_source, m_row); }
     };
 
-    static int other(int x, int y, int z) { SASSERT(x == z || y == z); return x == z ? y : x; }
-    std::ostream& print_vert(std::ostream & out, const vertex* v) const {
+    static int other(int x, int y, int z) {
+        SASSERT(x == z || y == z);
+        return x == z ? y : x;
+    }
+    std::ostream& print_vert(std::ostream& out, const vertex* v) const {
         out << "(c = " << v->column() << ", parent = {";
         if (v->parent())
             out << "(" << v->parent()->column() << ")";
         else
-            out << "null"; 
-        out <<  "} , lvl = " << v->level();
-        if (m_pol.contains(v->column())) 
-            out << (pol(v) == -1? " -":" +");
+            out << "null";
+        out << "} , lvl = " << v->level();
+        if (m_pol.contains(v->column()))
+            out << (pol(v) == -1 ? " -" : " +");
         else
             out << " not in m_pol";
         out << ')';
         return out;
     }
 
-    hashtable<unsigned, u_hash, u_eq>         m_visited_rows;
-    hashtable<unsigned, u_hash, u_eq>         m_visited_columns;
-    u_map<vertex*>                            m_vertices;
-    vertex*                                   m_root = nullptr;
+    hashtable<unsigned, u_hash, u_eq> m_visited_rows;
+    hashtable<unsigned, u_hash, u_eq> m_visited_columns;
+    u_map<vertex*> m_vertices;
+    vertex* m_root = nullptr;
     // At some point we can find a row with a single vertex non fixed vertex
     // then we can fix the whole tree,
     // by adjusting the vertices offsets, so they become absolute.
     // If the tree is fixed then in addition to checking with the m_vals_to_verts
     // we are going to check with the m_fixed_var_tables.
-    const vertex*                             m_fixed_vertex = nullptr;
-    explanation                               m_fixed_vertex_explanation;
+    const vertex* m_fixed_vertex = nullptr;
+    explanation m_fixed_vertex_explanation;
     // a pair (o, j) belongs to m_vals_to_verts iff x[j] = x[m_root->column()] + o
-    map<mpq, const vertex*, obj_hash<mpq>, default_eq<mpq>>  m_vals_to_verts;
+    map<mpq, const vertex*, obj_hash<mpq>, default_eq<mpq>> m_vals_to_verts;
     // a pair (o, j) belongs to m_vals_to_verts_neg iff -x[j] = x[m_root->column()] + o
-    map<mpq, const vertex*, obj_hash<mpq>, default_eq<mpq>>  m_vals_to_verts_neg;
+    map<mpq, const vertex*, obj_hash<mpq>, default_eq<mpq>> m_vals_to_verts_neg;
     // x[m_root->column()] - m_pol[j].pol()*x[j] == const;
     // to bind polarity and the vertex in the table
-    u_map<int>                                m_pol; 
+    u_map<int> m_pol;
     // if m_pos.contains(j) then  x[j] = x[m_root->column()] + o
-    uint_set                                  m_pos;
-    
+    uint_set m_pos;
+
     // these maps map a column index to the corresponding index in ibounds
-    std::unordered_map<unsigned, unsigned>    m_improved_lower_bounds;
-    std::unordered_map<unsigned, unsigned>    m_improved_upper_bounds;
-    
-    T&                                        m_imp;
-    vector<implied_bound>                     m_ibounds;
+    std::unordered_map<unsigned, unsigned> m_improved_lower_bounds;
+    std::unordered_map<unsigned, unsigned> m_improved_upper_bounds;
 
+    T& m_imp;
+    vector<implied_bound> m_ibounds;
 
+    map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>> m_val2fixed_row;
 
-    map<mpq, unsigned, obj_hash<mpq>, default_eq<mpq>>  m_val2fixed_row;
-
-    bool is_fixed_row(unsigned r, unsigned & x) {
+    bool is_fixed_row(unsigned r, unsigned& x) {
         x = UINT_MAX;
-        const auto & row = lp().get_row(r);
+        const auto& row = lp().get_row(r);
         for (unsigned k = 0; k < row.size(); k++) {
             const auto& c = row[k];
             if (column_is_fixed(c.var()))
@@ -130,7 +133,7 @@ class lp_bound_propagator {
         }
         return x != UINT_MAX;
     }
-    
+
     void try_add_equation_with_internal_fixed_tables(unsigned r1) {
         SASSERT(m_fixed_vertex);
         unsigned v1, v2;
@@ -154,8 +157,8 @@ class lp_bound_propagator {
         TRACE("eq", print_row(tout, r1); print_row(tout, r2); tout << v1 << " == " << v2 << " = " << val(v1) << "\n");
         add_eq_on_columns(ex, v1, v2, true);
     }
-    
-    void try_add_equation_with_lp_fixed_tables(unsigned row_index, const vertex *v) {
+
+    void try_add_equation_with_lp_fixed_tables(unsigned row_index, const vertex* v) {
         SASSERT(m_fixed_vertex);
         unsigned v_j = v->column();
         unsigned j = null_lpvar;
@@ -163,23 +166,24 @@ class lp_bound_propagator {
             try_add_equation_with_internal_fixed_tables(row_index);
             return;
         }
-       
+
         TRACE("cheap_eq",
-              tout << "v_j = "; lp().print_column_info(v_j, tout) << std::endl;
+              tout << "v_j = ";
+              lp().print_column_info(v_j, tout) << std::endl;
               tout << "v = "; print_vert(tout, v) << std::endl;
-              tout << "found j " << j << std::endl; lp().print_column_info(j, tout)<< std::endl;
+              tout << "found j " << j << std::endl; lp().print_column_info(j, tout) << std::endl;
               tout << "found j = " << j << std::endl;);
         vector<edge> path = connect_in_tree(v, m_fixed_vertex);
         explanation ex = get_explanation_from_path(path);
         ex.add_expl(m_fixed_vertex_explanation);
         explain_fixed_column(j, ex);
-        add_eq_on_columns(ex, j, v_j, true);         
+        add_eq_on_columns(ex, j, v_j, true);
     }
 
-    void try_add_equation_with_val_table(const vertex *v) {
+    void try_add_equation_with_val_table(const vertex* v) {
         SASSERT(m_fixed_vertex);
         unsigned v_j = v->column();
-        const vertex *u = nullptr;
+        const vertex* u = nullptr;
         if (!m_vals_to_verts.find(val(v_j), u)) {
             m_vals_to_verts.insert(val(v_j), v);
             return;
@@ -187,7 +191,7 @@ class lp_bound_propagator {
         unsigned j = u->column();
         if (j == v_j || is_int(j) != is_int(v_j))
             return;
-    
+
         TRACE("cheap_eq", tout << "found j=" << j << " for v=";
               print_vert(tout, v) << "\n in m_vals_to_verts\n";);
         vector<edge> path = connect_in_tree(u, v);
@@ -198,7 +202,7 @@ class lp_bound_propagator {
 
     static bool not_set(unsigned j) { return j == UINT_MAX; }
     static bool is_set(unsigned j) { return j != UINT_MAX; }
-    
+
     void create_root(unsigned row_index) {
         SASSERT(!m_root && !m_fixed_vertex);
         unsigned x, y;
@@ -210,31 +214,30 @@ class lp_bound_propagator {
         }
         TRACE("cheap_eq", print_row(tout, row_index););
         m_root = alloc_v(x);
-        set_polarity(m_root, 1); // keep m_root in the positive table
+        set_polarity(m_root, 1);  // keep m_root in the positive table
         if (not_set(y)) {
             set_fixed_vertex(m_root);
             explain_fixed_in_row(row_index, m_fixed_vertex_explanation);
-        }
-        else {            
-            vertex *v = add_child_with_check(row_index, y, m_root, polarity);
+        } else {
+            vertex* v = add_child_with_check(row_index, y, m_root, polarity);
             if (v)
                 explore_under(v);
         }
         explore_under(m_root);
     }
 
-        void explore_under(vertex * v) {
+    void explore_under(vertex* v) {
         check_for_eq_and_add_to_val_tables(v);
         go_over_vertex_column(v);
     }
 
     // In case of only one non fixed column, and the function returns true,
     // this column would be represened by x.
-    bool is_tree_offset_row(unsigned row_index, unsigned & x, unsigned & y, int & polarity) const {
-        x = y =  UINT_MAX;
+    bool is_tree_offset_row(unsigned row_index, unsigned& x, unsigned& y, int& polarity) const {
+        x = y = UINT_MAX;
         const row_cell<mpq>* x_cell = nullptr;
         const row_cell<mpq>* y_cell = nullptr;
-        const auto & row = lp().get_row(row_index);
+        const auto& row = lp().get_row(row_index);
         for (unsigned k = 0; k < row.size(); k++) {
             const auto& c = row[k];
             if (column_is_fixed(c.var()))
@@ -242,25 +245,21 @@ class lp_bound_propagator {
             if (not_set(x)) {
                 if (c.coeff().is_one() || c.coeff().is_minus_one()) {
                     x = c.var();
-                    x_cell = & c;
-                }
-                else 
+                    x_cell = &c;
+                } else
                     return false;
-            }
-            else if (not_set(y)) {
+            } else if (not_set(y)) {
                 if (c.coeff().is_one() || c.coeff().is_minus_one()) {
                     y = c.var();
-                    y_cell = & c;
-                }
-                else
+                    y_cell = &c;
+                } else
                     return false;
-            }
-            else
+            } else
                 return false;
         }
         if (is_set(x)) {
             if (is_set(y))
-                polarity = x_cell->coeff().is_pos() == y_cell->coeff().is_pos()? -1 : 1;
+                polarity = x_cell->coeff().is_pos() == y_cell->coeff().is_pos() ? -1 : 1;
             else
                 polarity = 1;
             return true;
@@ -268,18 +267,18 @@ class lp_bound_propagator {
         return false;
     }
 
-    void go_over_vertex_column(vertex * v) {
+    void go_over_vertex_column(vertex* v) {
         lpvar j = v->column();
         if (!check_insert(m_visited_columns, j))
             return;
-        
-        for (const auto & c : lp().get_column(j)) {
+
+        for (const auto& c : lp().get_column(j)) {
             unsigned row_index = c.var();
             if (!check_insert(m_visited_rows, row_index))
                 continue;
             vertex* u = get_child_from_row(row_index, v);
-            if (u) 
-                explore_under(u);            
+            if (u)
+                explore_under(u);
         }
     }
 
@@ -295,21 +294,18 @@ class lp_bound_propagator {
         m_pol.reset();
         m_vertices.reset();
     }
-    
+
     struct reset_cheap_eq {
         lp_bound_propagator& p;
-        reset_cheap_eq(lp_bound_propagator& p):p(p) {}
+        reset_cheap_eq(lp_bound_propagator& p) : p(p) {}
         ~reset_cheap_eq() { p.reset_cheap_eq_eh(); }
     };
 
-
-public:
-    
-    lp_bound_propagator(T& imp):
-        m_imp(imp) {}
+   public:
+    lp_bound_propagator(T& imp) : m_imp(imp) {}
 
     const vector<implied_bound>& ibounds() const { return m_ibounds; }
-    
+
     void init() {
         m_improved_upper_bounds.clear();
         m_improved_lower_bounds.clear();
@@ -318,24 +314,24 @@ public:
 
     const lar_solver& lp() const { return m_imp.lp(); }
     lar_solver& lp() { return m_imp.lp(); }
-    
+
     column_type get_column_type(unsigned j) const {
         return m_imp.lp().get_column_type(j);
     }
-    
-    const impq & get_lower_bound(unsigned j) const {
+
+    const impq& get_lower_bound(unsigned j) const {
         return m_imp.lp().get_lower_bound(j);
     }
 
-    const mpq & get_lower_bound_rational(unsigned j) const {
+    const mpq& get_lower_bound_rational(unsigned j) const {
         return m_imp.lp().get_lower_bound(j).x;
     }
-    
-    const impq & get_upper_bound(unsigned j) const {
+
+    const impq& get_upper_bound(unsigned j) const {
         return m_imp.lp().get_upper_bound(j);
     }
 
-    const mpq & get_upper_bound_rational(unsigned j) const {
+    const mpq& get_upper_bound_rational(unsigned j) const {
         return m_imp.lp().get_upper_bound(j).x;
     }
 
@@ -343,40 +339,37 @@ public:
     bool column_is_fixed(lpvar j) const {
         return lp().column_is_fixed(j) && get_lower_bound(j).y.is_zero();
     }
-    
-    void try_add_bound(mpq const& v, unsigned j, bool is_low, bool coeff_before_j_is_pos, unsigned row_or_term_index, bool strict) {
-        j = m_imp.lp().column_to_reported_index(j);    
 
-        lconstraint_kind kind = is_low? GE : LE;
+    void try_add_bound(mpq const& v, unsigned j, bool is_low, bool coeff_before_j_is_pos, unsigned row_or_term_index, bool strict) {
+        j = m_imp.lp().column_to_reported_index(j);
+
+        lconstraint_kind kind = is_low ? GE : LE;
         if (strict)
             kind = static_cast<lconstraint_kind>(kind / 2);
-    
+
         if (!m_imp.bound_is_interesting(j, kind, v))
             return;
-        unsigned k; // index to ibounds
+        unsigned k;  // index to ibounds
         if (is_low) {
             if (try_get_value(m_improved_lower_bounds, j, k)) {
-                auto & found_bound = m_ibounds[k];
+                auto& found_bound = m_ibounds[k];
                 if (v > found_bound.m_bound || (v == found_bound.m_bound && !found_bound.m_strict && strict)) {
                     found_bound = implied_bound(v, j, is_low, coeff_before_j_is_pos, row_or_term_index, strict);
                     TRACE("try_add_bound", m_imp.lp().print_implied_bound(found_bound, tout););
                 }
-            }
-            else {
+            } else {
                 m_improved_lower_bounds[j] = m_ibounds.size();
                 m_ibounds.push_back(implied_bound(v, j, is_low, coeff_before_j_is_pos, row_or_term_index, strict));
                 TRACE("try_add_bound", m_imp.lp().print_implied_bound(m_ibounds.back(), tout););
             }
-        }
-        else { // the upper bound case
+        } else {  // the upper bound case
             if (try_get_value(m_improved_upper_bounds, j, k)) {
-                auto & found_bound = m_ibounds[k];
+                auto& found_bound = m_ibounds[k];
                 if (v < found_bound.m_bound || (v == found_bound.m_bound && !found_bound.m_strict && strict)) {
                     found_bound = implied_bound(v, j, is_low, coeff_before_j_is_pos, row_or_term_index, strict);
                     TRACE("try_add_bound", m_imp.lp().print_implied_bound(found_bound, tout););
                 }
-            }
-            else {
+            } else {
                 m_improved_upper_bounds[j] = m_ibounds.size();
                 m_ibounds.push_back(implied_bound(v, j, is_low, coeff_before_j_is_pos, row_or_term_index, strict));
                 TRACE("try_add_bound", m_imp.lp().print_implied_bound(m_ibounds.back(), tout););
@@ -388,18 +381,18 @@ public:
         m_imp.consume(a, ci);
     }
 
-    const mpq& val(unsigned j) const { 
-        return  lp().get_column_value(j).x;
+    const mpq& val(unsigned j) const {
+        return lp().get_column_value(j).x;
     }
-    
+
     const mpq& val(const vertex* v) const {
         return val(v->column());
     }
-    
-    bool tree_contains_r(vertex* root, vertex *v) const {
+
+    bool tree_contains_r(vertex* root, vertex* v) const {
         if (*root == *v)
             return true;
-        for (auto e : root->edges()) 
+        for (auto e : root->edges())
             if (tree_contains_r(e.target(), v))
                 return true;
         return false;
@@ -415,7 +408,7 @@ public:
         m_pol.insert(j, p);
     }
 
-    void check_and_set_polarity(vertex* v, int polarity, unsigned row_index, vertex*v_parent) {
+    void check_and_set_polarity(vertex* v, int polarity, unsigned row_index, vertex* v_parent) {
         int prev_pol;
         if (!m_pol.find(v->column(), prev_pol)) {
             set_polarity(v, polarity);
@@ -429,106 +422,107 @@ public:
         m_fixed_vertex_explanation = get_explanation_from_path(path);
         explain_fixed_in_row(row_index, m_fixed_vertex_explanation);
         set_fixed_vertex(v);
-        TRACE("cheap_eq", 
-              tout << "polarity switch: " << polarity << "\nv = "; print_vert(tout , v) << "\nu = ";        tout << "fixed vertex explanation\n";
-              for (auto p : m_fixed_vertex_explanation) 
-                  lp().constraints().display(tout, [this](lpvar j) { return lp().get_variable_name(j);}, p.ci()););
-
+        TRACE("cheap_eq",
+              tout << "polarity switch: " << polarity << "\nv = ";
+              print_vert(tout, v) << "\nu = "; tout << "fixed vertex explanation\n";
+              for (auto p
+                   : m_fixed_vertex_explanation)
+                  lp()
+                      .constraints()
+                      .display(
+                          tout, [this](lpvar j) { return lp().get_variable_name(j); }, p.ci()););
     }
-    
-    bool tree_contains(vertex *v) const {
+
+    bool tree_contains(vertex* v) const {
         if (!m_root)
             return false;
         return tree_contains_r(m_root, v);
     }
-    
-    vertex * alloc_v(unsigned column) {
-        vertex * v = alloc(vertex, column);
+
+    vertex* alloc_v(unsigned column) {
+        vertex* v = alloc(vertex, column);
         m_vertices.insert(column, v);
         SASSERT(!tree_contains(v));
         return v;
     }
 
-
     unsigned column(unsigned row, unsigned index) {
         return lp().get_row(row)[index].var();
     }
 
     bool fixed_phase() const { return m_fixed_vertex; }
-    
-    
+
     // Returns the vertex to start exploration from, or nullptr.
     // It is assumed that parent->column() is present in the row
     vertex* get_child_from_row(unsigned row_index, vertex* parent) {
         TRACE("cheap_eq_det", print_row(tout, row_index););
-        unsigned x, y; int row_polarity;
+        unsigned x, y;
+        int row_polarity;
         if (!is_tree_offset_row(row_index, x, y, row_polarity)) {
-            TRACE("cheap_eq_det", tout << "not an offset row\n"; );
+            TRACE("cheap_eq_det", tout << "not an offset row\n";);
             return nullptr;
         }
-        if (not_set(y)) { // there is only one fixed variable in the row
+        if (not_set(y)) {  // there is only one fixed variable in the row
             if (!fixed_phase()) {
                 set_fixed_vertex(parent);
                 explain_fixed_in_row(row_index, m_fixed_vertex_explanation);
             }
             return nullptr;
         }
-        
+
         SASSERT(is_set(x) && is_set(y));
         unsigned col = other(x, y, parent->column());
         return add_child_with_check(row_index, col, parent, row_polarity);
     }
 
-    vertex * add_child_with_check(unsigned row_index, unsigned col, vertex* parent, int row_polarity) {
+    vertex* add_child_with_check(unsigned row_index, unsigned col, vertex* parent, int row_polarity) {
         vertex* vy;
         if (m_vertices.find(col, vy)) {
             SASSERT(vy != nullptr);
             if (!fixed_phase()) {
-                check_and_set_polarity(vy, pol(parent) * row_polarity, row_index, parent); 
+                check_and_set_polarity(vy, pol(parent) * row_polarity, row_index, parent);
             }
-            return nullptr; // it is not a new vertex
-        } 
+            return nullptr;  // it is not a new vertex
+        }
         vy = alloc_v(col);
         parent->add_child(row_index, vy);
         if (!fixed_phase())
             check_and_set_polarity(vy, row_polarity * pol(parent), row_index, parent);
-        return vy;        
+        return vy;
     }
 
-    bool is_equal(lpvar j, lpvar k) const {        
+    bool is_equal(lpvar j, lpvar k) const {
         return m_imp.is_equal(col_to_imp(j), col_to_imp(k));
     }
 
-    void check_for_eq_and_add_to_val_table(vertex* v,  map<mpq, const vertex*, obj_hash<mpq>, default_eq<mpq>>& table) {
+    void check_for_eq_and_add_to_val_table(vertex* v, map<mpq, const vertex*, obj_hash<mpq>, default_eq<mpq>>& table) {
         TRACE("cheap_eq", tout << "v = "; print_vert(tout, v) << "\n";);
-        const vertex *k; // the other vertex
+        const vertex* k;  // the other vertex
         if (table.find(val(v), k)) {
-            TRACE("cheap_eq", tout << "found k " ; print_vert(tout, k) << "\n";);
+            TRACE("cheap_eq", tout << "found k "; print_vert(tout, k) << "\n";);
             if (k->column() != v->column() &&
                 is_int(k->column()) == is_int(v->column()) &&
                 !is_equal(k->column(), v->column())) {
                 report_eq(k, v);
-            }
-            else {
+            } else {
                 TRACE("cheap_eq", tout << "no report\n";);
             }
-        }
-        else {
+        } else {
             TRACE("cheap_eq", tout << "registered: " << val(v) << " -> { "; print_vert(tout, v) << "} \n";);
             table.insert(val(v), v);
         }
     }
-    
+
     void check_for_eq_and_add_to_val_tables(vertex* v) {
         TRACE("cheap_eq_det", print_vert(tout, v) << "\n";);
         if (!fixed_phase()) {
             if (pol(v->column()) == -1)
                 check_for_eq_and_add_to_val_table(v, m_vals_to_verts_neg);
-            else 
+            else
                 check_for_eq_and_add_to_val_table(v, m_vals_to_verts);
         }
     }
-        
+
     void clear_for_eq() {
         m_visited_rows.reset();
         m_visited_columns.reset();
@@ -542,41 +536,40 @@ public:
 
     std::ostream& print_path(const vector<edge>& path, std::ostream& out) const {
         out << "path = \n";
-        for (const edge& k : path) 
+        for (const edge& k : path)
             print_edge(k, out) << "\n";
         return out;
     }
-    
+
     // we have v_i and v_j, indices of vertices at the same offsets
     void report_eq(const vertex* v_i, const vertex* v_j) {
         SASSERT(v_i != v_j);
         SASSERT(lp().get_column_value(v_i->column()) == lp().get_column_value(v_j->column()));
         TRACE("cheap_eq", tout << v_i->column() << " = " << v_j->column() << "\nu = ";
-              print_vert(tout, v_i) << "\nv = "; print_vert(tout, v_j) <<"\n");
-        
+              print_vert(tout, v_i) << "\nv = "; print_vert(tout, v_j) << "\n");
+
         vector<edge> path = connect_in_tree(v_i, v_j);
         lp::explanation exp = get_explanation_from_path(path);
         add_eq_on_columns(exp, v_i->column(), v_j->column(), false);
-        
     }
 
-    std::ostream& print_expl(std::ostream & out, const explanation& exp) const {
-        for (auto p : exp) 
-            lp().constraints().display(out, [this](lpvar j) { return lp().get_variable_name(j);}, p.ci());
+    std::ostream& print_expl(std::ostream& out, const explanation& exp) const {
+        for (auto p : exp)
+            lp().constraints().display(
+                out, [this](lpvar j) { return lp().get_variable_name(j); }, p.ci());
         return out;
     }
-    
+
     bool add_eq_on_columns(const explanation& exp, lpvar j, lpvar k, bool is_fixed) {
         SASSERT(j != k);
         unsigned je = lp().column_to_reported_index(j);
         unsigned ke = lp().column_to_reported_index(k);
-        TRACE("cheap_eq", 
-              tout << "reporting eq " << j  << ", " << k << "\n";
+        TRACE("cheap_eq",
+              tout << "reporting eq " << j << ", " << k << "\n";
               tout << "reported idx " << je << ", " << ke << "\n";
               print_expl(tout, exp);
-              tout << "theory_vars v" << lp().local_to_external(je) << " == v" << lp().local_to_external(ke) << "\n";
-              );
-        
+              tout << "theory_vars v" << lp().local_to_external(je) << " == v" << lp().local_to_external(ke) << "\n";);
+
         bool added = m_imp.add_eq(je, ke, exp, is_fixed);
         if (added) {
             if (is_fixed)
@@ -600,75 +593,75 @@ public:
     bool is_int(lpvar j) const {
         return lp().column_is_int(j);
     }
-    
+
     explanation get_explanation_from_path(vector<edge>& path) const {
         explanation ex;
-        for (edge &e : path)
+        for (edge& e : path)
             explain_fixed_in_row(e.row(), ex);
         return ex;
     }
 
     void explain_fixed_in_row(unsigned row, explanation& ex) const {
         TRACE("cheap_eq", tout << lp().get_row(row) << std::endl);
-        for (const auto & c : lp().get_row(row)) 
-            if (lp().is_fixed(c.var())) 
+        for (const auto& c : lp().get_row(row))
+            if (lp().is_fixed(c.var()))
                 explain_fixed_column(c.var(), ex);
     }
 
-    void explain_fixed_column(unsigned j, explanation & ex) const {
+    void explain_fixed_column(unsigned j, explanation& ex) const {
         SASSERT(column_is_fixed(j));
-        constraint_index lc, uc;            
+        constraint_index lc, uc;
         lp().get_bound_constraint_witnesses_for_column(j, lc, uc);
         ex.push_back(lc);
-        ex.push_back(uc);        
+        ex.push_back(uc);
     }
-    
+
     vector<edge> connect_in_tree(const vertex* u, const vertex* v) const {
         vector<edge> path;
-        TRACE("cheap_eq_details", tout << "u = " ; print_vert(tout, u); tout << "\nv = ";print_vert(tout, v) << "\n";); 
+        TRACE("cheap_eq_details", tout << "u = "; print_vert(tout, u); tout << "\nv = "; print_vert(tout, v) << "\n";);
         vector<edge> v_branch;
         // equalize the levels
         while (u->level() > v->level()) {
             path.push_back(u->edge_from_parent().reverse());
             u = u->parent();
         }
-        
+
         while (u->level() < v->level()) {
             v_branch.push_back(v->edge_from_parent());
             v = v->parent();
         }
         SASSERT(u->level() == v->level());
-        TRACE("cheap_eq_details", tout << "u = " ; print_vert(tout, u); tout << "\nv = "; print_vert(tout, v) << "\n";); 
+        TRACE("cheap_eq_details", tout << "u = "; print_vert(tout, u); tout << "\nv = "; print_vert(tout, v) << "\n";);
         while (u != v) {
             path.push_back(u->edge_from_parent().reverse());
             v_branch.push_back(v->edge_from_parent());
             u = u->parent();
             v = v->parent();
         }
-        for (unsigned i = v_branch.size(); i--; ) {
-           path.push_back(v_branch[i]);
+        for (unsigned i = v_branch.size(); i--;) {
+            path.push_back(v_branch[i]);
         }
         TRACE("cheap_eq", print_path(path, tout););
         return path;
     }
-    
+
     bool tree_is_correct() const {
         std::unordered_set<int> vs;
         return tree_is_correct(m_root, vs);
     }
-    
+
     bool tree_is_correct(vertex* v, std::unordered_set<int>& visited_verts) const {
         if (fixed_phase())
             return true;
-        if (visited_verts.find(v->column()) != visited_verts.end()) 
+        if (visited_verts.find(v->column()) != visited_verts.end())
             return false;
         visited_verts.insert(v->column());
-        for (auto e : v->edges()) 
+        for (auto e : v->edges())
             if (!tree_is_correct(e.target(), visited_verts))
                 return false;
         return true;
     }
-    std::ostream& print_tree(std::ostream & out, vertex* v) const {
+    std::ostream& print_tree(std::ostream& out, vertex* v) const {
         print_vert(out, v);
         out << "\nchildren :\n";
         for (auto c : v->edges()) {
@@ -683,75 +676,74 @@ public:
         try_add_equation_with_lp_fixed_tables(row_index, v);
         try_add_equation_with_val_table(v);
     }
-    
+
     void handle_fixed_phase(unsigned row_index) {
         if (!fixed_phase())
             return;
         const vertex* v = m_root;
         try_add_equation_with_fixed_tables(row_index, v);
-        for (auto e: v->edges())
+        for (auto e : v->edges())
             try_add_equation_with_fixed_tables(row_index, e.target());
     }
 
-    
     void cheap_eq_tree(unsigned row_index) {
         reset_cheap_eq _reset(*this);
         TRACE("cheap_eq_det", tout << "row_index = " << row_index << "\n";);
-        if (!check_insert(m_visited_rows, row_index)) 
+        if (!check_insert(m_visited_rows, row_index))
             return;
         create_root(row_index);
         if (!m_root)
             return;
-        
-        TRACE("cheap_eq", tout << "tree = "; print_tree(tout, m_root) << "\n";);        
+
+        TRACE("cheap_eq", tout << "tree = "; print_tree(tout, m_root) << "\n";);
         SASSERT(tree_is_correct());
         handle_fixed_phase(row_index);
-        
+
         TRACE("cheap_eq",
               tout << "done for row_index " << row_index << "\n";
               tout << "tree size = " << verts_size(););
     }
 
-    std::ostream& print_row(std::ostream & out, unsigned row_index) const  {
-        unsigned x, y; int polarity;
+    std::ostream& print_row(std::ostream& out, unsigned row_index) const {
+        unsigned x, y;
+        int polarity;
         if (true || !is_tree_offset_row(row_index, x, y, polarity))
             return lp().get_int_solver()->display_row_info(out, row_index);
-        
+
         bool first = true;
-        for (const auto &c: lp().A_r().m_rows[row_index]) {
+        for (const auto& c : lp().A_r().m_rows[row_index]) {
             if (lp().column_is_fixed(c.var()))
                 continue;
             if (c.coeff().is_one()) {
-                    if (!first)
-                        out << "+";
-                }
-                else if (c.coeff().is_minus_one())
-                    out << "-";                             
+                if (!first)
+                    out << "+";
+            } else if (c.coeff().is_minus_one())
+                out << "-";
             out << lp().get_variable_name(c.var()) << " ";
-            first = false;        
+            first = false;
         }
         out << "\n";
         return out;
     }
-    
-    void set_fixed_vertex(vertex *v) {        
+
+    void set_fixed_vertex(vertex* v) {
         TRACE("cheap_eq", if (v) print_vert(tout, v); else tout << "set m_fixed_vertex to nullptr"; tout << "\n";);
         SASSERT(!m_fixed_vertex || v == nullptr);
         m_fixed_vertex = v;
     }
-    
+
     unsigned verts_size() const {
         return subtree_size(m_root);
     }
 
     unsigned subtree_size(vertex* v) const {
-        unsigned r = 1; // 1 for v
+        unsigned r = 1;  // 1 for v
         for (auto e : v->edges())
             r += subtree_size(e.target());
         return r;
     }
-            
-    void delete_tree(vertex * v) {
+
+    void delete_tree(vertex* v) {
         for (auto p : v->edges())
             delete_tree(p.target());
         dealloc(v);
@@ -763,7 +755,6 @@ public:
             return false;
         table.insert(j);
         return true;
-    }        
-
+    }
 };
-}
+}  // namespace lp

src/math/lp/lp_core_solver_base.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include <utility>
 #include <memory>
 #include <string>
@@ -63,8 +64,8 @@ template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::colu
 template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq>>::pivot_column_tableau(unsigned int, unsigned int);
 template bool lp::lp_core_solver_base<lp::mpq, lp::mpq>::pivot_column_tableau(unsigned int, unsigned int);
 template void lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::transpose_rows_tableau(unsigned int, unsigned int);
-template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::inf_set_is_correct() const;
-template bool lp::lp_core_solver_base<lp::mpq, lp::mpq>::inf_set_is_correct() const;
+template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::inf_heap_is_correct() const;
+template bool lp::lp_core_solver_base<lp::mpq, lp::mpq>::inf_heap_is_correct() const;
 template bool lp::lp_core_solver_base<lp::mpq, lp::numeric_pair<lp::mpq> >::remove_from_basis(unsigned int);
 
 

src/math/lp/lp_core_solver_base.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include <set>
 #include "util/vector.h"
@@ -28,9 +29,13 @@ Revision History:
 #include "math/lp/permutation_matrix.h"
 #include "math/lp/column_namer.h"
 #include "math/lp/u_set.h"
-
+#include "util/heap.h"
 
 namespace lp {
+struct lpvar_lt {
+  bool operator()(lpvar v1, lpvar v2) const { return v1 < v2; }
+};
+typedef heap<lpvar_lt> lpvar_heap;
 template <typename T, typename X>
 X dot_product(const vector<T> & a, const vector<X> & b) {
     lp_assert(a.size() == b.size());
@@ -50,25 +55,25 @@ private:
     lp_status m_status;
 public:
     bool current_x_is_feasible() const {
-        TRACE("feas",
-              if (m_inf_set.size()) {
-                  tout << "column " << *m_inf_set.begin() << " is infeasible" << std::endl;
-                  print_column_info(*m_inf_set.begin(), tout);
+        TRACE("feas_bug",
+              if (!m_inf_heap.empty()) {
+                  tout << "column " << *m_inf_heap.begin() << " is infeasible" << std::endl;
+                  print_column_info(*m_inf_heap.begin(), tout);
               } else {
                   tout << "x is feasible\n";
               }
               );
-        return m_inf_set.size() == 0;
+        return m_inf_heap.empty();
     }
-    bool current_x_is_infeasible() const { return m_inf_set.size() != 0; }
+    bool current_x_is_infeasible() const { return m_inf_heap.size() != 0; }
 private:
-    u_set m_inf_set;
+    lpvar_heap m_inf_heap;
 public:
-    const u_set& inf_set() const { return m_inf_set; }
-    u_set& inf_set() { return m_inf_set; }
-    void inf_set_increase_size_by_one() { m_inf_set.increase_size_by_one(); }
-    bool inf_set_contains(unsigned j) const { return m_inf_set.contains(j); }
-    unsigned inf_set_size() const { return m_inf_set.size(); }    
+    const lpvar_heap& inf_heap() const { return m_inf_heap; }
+    lpvar_heap& inf_heap() { return m_inf_heap; }
+    void inf_heap_increase_size_by_one() { m_inf_heap.reserve(m_inf_heap.size() + 1); }
+    bool inf_heap_contains(unsigned j) const { return m_inf_heap.contains(j); }
+    unsigned inf_heap_size() const { return m_inf_heap.size(); }    
     indexed_vector<T>     m_pivot_row; // this is the real pivot row of the simplex tableu
     static_matrix<T, X> & m_A; // the matrix A
     // vector<X> const &           m_b; // the right side
@@ -255,7 +260,7 @@ public:
 
     bool calc_current_x_is_feasible_include_non_basis() const;
 
-    bool inf_set_is_correct() const;
+    bool inf_heap_is_correct() const;
     
     bool column_is_dual_feasible(unsigned j) const;
 
@@ -304,9 +309,9 @@ public:
         case column_type::boxed:
             if (x < m_lower_bounds[j]) {
                 delta = m_lower_bounds[j] - x;
-                ret = true;;
+                ret = true;
             }
-            if (x > m_upper_bounds[j]) {
+            else if (x > m_upper_bounds[j]) {
                 delta = m_upper_bounds[j] - x;
                 ret = true;
             }
@@ -526,63 +531,55 @@ public:
         swap(m_basis_heading, m_basis[i], m_basis[ii]);
     }
 
-    bool column_is_in_inf_set(unsigned j) const {
-        return m_inf_set.contains(j);
+    bool column_is_in_inf_heap(unsigned j) const {
+        return m_inf_heap.contains(j);
     }
 
     bool column_is_base(unsigned j) const {
         return m_basis_heading[j] >= 0;
     }
 
-    
-    void update_x_with_feasibility_tracking(unsigned j, const X & v) {
-        TRACE("lar_solver", tout << "j = " << j << ", v = " << v << "\n";);
-        m_x[j] = v;
-        track_column_feasibility(j);
-    }
-
     void add_delta_to_x_and_track_feasibility(unsigned j, const X & del) {
-        TRACE("lar_solver", tout << "del = " << del << ", was x[" << j << "] = " << m_x[j] << "\n";);
+        TRACE("lar_solver_feas_bug", tout << "del = " << del << ", was x[" << j << "] = " << m_x[j] << "\n";);
         m_x[j] += del;
-        TRACE("lar_solver", tout << "became x[" << j << "] = " << m_x[j] << "\n";);
+        TRACE("lar_solver_feas_bug", tout << "became x[" << j << "] = " << m_x[j] << "\n";);
         track_column_feasibility(j);
     }
 
     void update_x(unsigned j, const X & v) {
-        TRACE("lar_solver", tout << "j = " << j << ", v = " << v << "\n";);
         m_x[j] = v;
+        TRACE("lar_solver_feas", tout << "not tracking feas j = " << j << ", v = " << v << (column_is_feasible(j)? " feas":" non-feas") << "\n";);
     }
 
-    void add_delta_to_x(unsigned j, const X & delta) {
-        TRACE("lar_solver", tout << "j = " << j << ", delta = " << delta << "\n";);
+    void add_delta_to_x(unsigned j, const X& delta) {
         m_x[j] += delta;
+        TRACE("lar_solver_feas", tout << "not tracking feas j = " << j << " v = " << m_x[j] << " delta = " << delta << (column_is_feasible(j) ? " feas" : " non-feas") << "\n";);
     }
-   
+        
     void track_column_feasibility(unsigned j) {
         if (column_is_feasible(j))
-            remove_column_from_inf_set(j);
+            remove_column_from_inf_heap(j);
         else
-            insert_column_into_inf_set(j);
+            insert_column_into_inf_heap(j);
     }
-    void insert_column_into_inf_set(unsigned j) {
-        TRACE("lar_solver", tout << "j = " << j << "\n";);
-        m_inf_set.insert(j);
+    void insert_column_into_inf_heap(unsigned j) {        
+		if (!m_inf_heap.contains(j)) {
+	        m_inf_heap.insert(j);
+            TRACE("lar_solver_inf_heap", tout << "insert into inf_heap j = " << j << "\n";);
+        }
         lp_assert(!column_is_feasible(j));
     }
-    void remove_column_from_inf_set(unsigned j) {
-        TRACE("lar_solver", tout << "j = " << j << "\n";);
-        m_inf_set.erase(j);
+    void remove_column_from_inf_heap(unsigned j) {
+		if (m_inf_heap.contains(j)) {
+            TRACE("lar_solver_inf_heap", tout << "insert into heap j = " << j << "\n";);
+        	m_inf_heap.erase(j);
+        }
         lp_assert(column_is_feasible(j));
     }
 
-    void resize_inf_set(unsigned size) {
-        TRACE("lar_solver",);
-        m_inf_set.resize(size);
-    }
-
-    void clear_inf_set() {
-        TRACE("lar_solver",);
-        m_inf_set.clear();
+    void clear_inf_heap() {
+        TRACE("lar_solver_feas",);
+        m_inf_heap.clear();
     }
     
     bool costs_on_nbasis_are_zeros() const {

src/math/lp/lp_core_solver_base_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <set>
@@ -42,7 +43,7 @@ lp_core_solver_base(static_matrix<T, X> & A,
     m_total_iterations(0),
     m_iters_with_no_cost_growing(0),
     m_status(lp_status::FEASIBLE),
-    m_inf_set(A.column_count()),
+    m_inf_heap(std::max(static_cast<unsigned>(1024), A.column_count())),
     m_pivot_row(A.column_count()),
     m_A(A),
     m_basis(basis),
@@ -250,9 +251,9 @@ template <typename T, typename X> bool lp_core_solver_base<T, X>::calc_current_x
     return true;
 }
 
-template <typename T, typename X> bool lp_core_solver_base<T, X>::inf_set_is_correct() const {
+template <typename T, typename X> bool lp_core_solver_base<T, X>::inf_heap_is_correct() const {
     for (unsigned j = 0; j < this->m_n(); j++) {
-        bool belongs_to_set = m_inf_set.contains(j);
+        bool belongs_to_set = m_inf_heap.contains(j);
         bool is_feas = column_is_feasible(j);
         if (is_feas == belongs_to_set) {
             TRACE("lp_core", tout << "incorrectly set column in inf set "; print_column_info(j, tout) << "\n";);

src/math/lp/lp_primal_core_solver.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include <utility>
 #include <memory>
 #include <string>

src/math/lp/lp_primal_core_solver_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <list>

src/math/lp/lp_primal_core_solver_tableau_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 // this is a part of lp_primal_core_solver that deals with the tableau
@@ -30,7 +31,7 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::one_iteratio
     else {
         advance_on_entering_tableau(entering);
     }
-    lp_assert(this->inf_set_is_correct());
+    lp_assert(this->inf_heap_is_correct());
 }
 
 template <typename T, typename X> void lp_primal_core_solver<T, X>::advance_on_entering_tableau(int entering) {
@@ -59,6 +60,7 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::advance_on_e
     }
     unsigned j_nz = this->m_m() + 1; // this number is greater than the max column size
     std::list<unsigned>::iterator entering_iter = m_non_basis_list.end();
+    unsigned n = 0;
     for (auto non_basis_iter = m_non_basis_list.begin(); number_of_benefitial_columns_to_go_over && non_basis_iter != m_non_basis_list.end(); ++non_basis_iter) {
         unsigned j = *non_basis_iter;
         if (!column_is_benefitial_for_entering_basis(j))
@@ -71,8 +73,9 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::advance_on_e
             entering_iter = non_basis_iter;
             if (number_of_benefitial_columns_to_go_over)
                 number_of_benefitial_columns_to_go_over--;
+            n = 1;
         }
-        else if (t == j_nz && this->m_settings.random_next() % 2 == 0) {
+        else if (t == j_nz && this->m_settings.random_next(++n) == 0) {
             entering_iter = non_basis_iter;
         }
     }// while (number_of_benefitial_columns_to_go_over && initial_offset_in_non_basis != offset_in_nb);
@@ -166,7 +169,8 @@ template <typename T, typename X>void lp_primal_core_solver<T, X>::advance_on_en
         }
         this->update_basis_and_x_tableau(entering, leaving, t);
         this->iters_with_no_cost_growing() = 0;
-    } else {
+    }
+    else {
         this->pivot_column_tableau(entering, this->m_basis_heading[leaving]);
         this->change_basis(entering, leaving);
     }
@@ -256,7 +260,7 @@ template <typename T, typename X> void lp_primal_core_solver<T, X>::init_run_tab
         this->m_basis_sort_counter = 0; // to initiate the sort of the basis
         //  this->set_total_iterations(0);
         this->iters_with_no_cost_growing() = 0;
-		lp_assert(this->inf_set_is_correct());
+		lp_assert(this->inf_heap_is_correct());
         if (this->current_x_is_feasible() && this->m_look_for_feasible_solution_only)
             return;
         if (this->m_settings.backup_costs)

src/math/lp/lp_settings.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include <memory>
 #include "util/vector.h"
 #include "smt/params/smt_params_helper.hpp"

src/math/lp/lp_settings.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 
 #pragma once
 #include "util/vector.h"
@@ -218,6 +219,8 @@ public:
     unsigned hnf_cut_period() const { return m_hnf_cut_period; }
     void set_hnf_cut_period(unsigned period) { m_hnf_cut_period = period;  }
     unsigned random_next() { return m_rand(); }
+    unsigned random_next(unsigned u ) { return m_rand(u); }
+    
     void set_random_seed(unsigned s) { m_rand.set_seed(s); }
 
     bool bound_progation() const { 

src/math/lp/lp_settings_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <cmath>

src/math/lp/lp_types.h

@@ -17,11 +17,13 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <sstream>
-
-
+// clang-format off
+#include <limits.h>
+#include "util/debug.h"
 namespace nla {
     class core;
 }

src/math/lp/lp_utils.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include <string>
 #include "math/lp/numeric_pair.h"

src/math/lp/matrix.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include "math/lp/lp_settings.h"
 #include "math/lp/matrix_def.h"
 #include "math/lp/static_matrix.h"

src/math/lp/matrix.h

@@ -6,6 +6,7 @@ Author:
     Lev Nachmanson (levnach)
 
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/numeric_pair.h"
 #include "util/vector.h"

src/math/lp/matrix_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include <cmath>

src/math/lp/monomial_bounds.cpp

@@ -6,6 +6,7 @@
   Lev Nachmanson (levnach)
 
   --*/
+// clang-format off
 
 #include "math/lp/monomial_bounds.h"
 #include "math/lp/nla_core.h"

src/math/lp/monomial_bounds.h

@@ -6,6 +6,7 @@
   Lev Nachmanson (levnach)
 
   --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/nla_common.h"

src/math/lp/nex.h

@@ -4,6 +4,7 @@
   Author:
     Lev Nachmanson (levnach)
 --*/
+// clang-format off
 
 #pragma once
 #include <initializer_list>

src/math/lp/nex_creator.cpp

@@ -4,6 +4,7 @@
   Author:
    Lev Nachmanson (levnach)
   --*/
+// clang-format off
 #include "util/lbool.h"
 #include "math/lp/nex_creator.h"
 #include <map>

src/math/lp/nex_creator.h

@@ -5,6 +5,7 @@
    Lev Nachmanson (levnach)
    Nikolaj Bjorner (nbjorner)
   --*/
+// clang-format off
 #pragma once
 #include <map>
 #include <set>

src/math/lp/nla_basics_lemmas.cpp

@@ -6,6 +6,7 @@
     Nikolaj Bjorner (nbjorner)
 
 --*/
+// clang-format off
 
 #include "math/lp/nla_basics_lemmas.h"
 #include "math/lp/nla_core.h"

src/math/lp/nla_basics_lemmas.h

@@ -6,6 +6,7 @@
     Nikolaj Bjorner (nbjorner)
 
   --*/
+// clang-format off
 #pragma once
 #include "math/lp/monic.h"    
 #include "math/lp/factorization.h"    

src/math/lp/nla_common.cpp

@@ -6,6 +6,7 @@
     Lev Nachmanson (levnach)
     Nikolaj Bjorner (nbjorner)
   --*/
+// clang-format off
 #include "math/lp/nla_common.h"
 #include "math/lp/nla_core.h"
 

src/math/lp/nla_common.h

@@ -6,6 +6,7 @@
    Nikolaj Bjorner (nbjorner)
 
   --*/
+// clang-format off
 #pragma once
 #include "util/rational.h"
 #include "math/lp/nla_defs.h"

src/math/lp/nla_core.cpp

@@ -10,6 +10,7 @@ Author:
     Nikolaj Bjorner (nbjorner)
 
 --*/
+// clang-format off
 #include "util/uint_set.h"
 #include "math/lp/nla_core.h"
 #include "math/lp/factorization_factory_imp.h"
@@ -166,7 +167,9 @@ bool core::check_monic(const monic& m) const {
     if (!is_relevant(m.var()))
         return true;
 #endif
-    SASSERT((!m_lar_solver.column_is_int(m.var())) || m_lar_solver.get_column_value(m.var()).is_int());
+    if (m_lar_solver.column_is_int(m.var()) && !m_lar_solver.get_column_value(m.var()).is_int())
+        return true;
+    
     bool ret = product_value(m) == m_lar_solver.get_column_value(m.var()).x; 
     CTRACE("nla_solver_check_monic", !ret, print_monic(m, tout) << '\n';);
     return ret;
@@ -1436,7 +1439,7 @@ void core::patch_monomials() {
     } else {
         m_lar_solver.pop();
         restore_tableau();
-        m_lar_solver.clear_inf_set();
+        m_lar_solver.clear_inf_heap();
     }
     SASSERT(m_lar_solver.ax_is_correct());
 }

src/math/lp/nla_core.h

@@ -10,6 +10,7 @@
     Nikolaj Bjorner (nbjorner)
 
 --*/
+// clang-format off
 #pragma once
 #include "math/lp/factorization.h"
 #include "math/lp/lp_types.h"

src/math/lp/nla_defs.h

@@ -8,6 +8,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 #include "math/lp/lp_types.h"
 #include "math/lp/column_info.h"

src/math/lp/nla_divisions.cpp

@@ -14,6 +14,7 @@ Description:
   Check divisions
 
 --*/
+// clang-format off
 #include "math/lp/nla_core.h"
 
 namespace nla {

src/math/lp/nla_divisions.h

@@ -13,6 +13,7 @@ Description:
   Check division constraints.
   
 --*/
+// clang-format off
 
 #include "math/lp/nla_types.h"
 

src/math/lp/nla_grobner.cpp

@@ -10,6 +10,7 @@ Author:
     Nikolaj Bjorner (nbjorner)
 
 --*/
+// clang-format off
 #include "util/uint_set.h"
 #include "math/lp/nla_core.h"
 #include "math/lp/factorization_factory_imp.h"

src/math/lp/nla_grobner.h

@@ -6,6 +6,7 @@
   Lev Nachmanson (levnach)
 
   --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/nla_common.h"

src/math/lp/nla_intervals.h

@@ -6,6 +6,7 @@
     Nikolaj Bjorner (nbjorner)
 
   --*/
+// clang-format off
 #pragma once
 #include "util/dependency.h"
 #include "util/region.h"

src/math/lp/nla_monotone_lemmas.cpp

@@ -5,6 +5,7 @@
     Lev Nachmanson (levnach)
     Nikolaj Bjorner (nbjorner)
   --*/
+// clang-format off
 #include "math/lp/nla_basics_lemmas.h"
 #include "math/lp/nla_core.h"
 namespace nla {

src/math/lp/nla_monotone_lemmas.h

@@ -5,6 +5,7 @@
    Lev Nachmanson (levnach)
    Nikolaj Bjorner (nbjorner)
   --*/
+// clang-format off
 #pragma once
 namespace nla {
 class core;

src/math/lp/nla_order_lemmas.cpp

@@ -5,6 +5,7 @@
     Nikolaj Bjorner (nbjorner)
     Lev Nachmanson (levnach)
   --*/
+// clang-format off
 
 #include "math/lp/nla_order_lemmas.h"
 #include "math/lp/nla_core.h"

src/math/lp/nla_order_lemmas.h

@@ -6,6 +6,7 @@
     Nikolaj Bjorner (nbjorner)
 
   --*/
+// clang-format off
 #pragma once
 #include "math/lp/factorization.h"
 #include "math/lp/nla_common.h"

src/math/lp/nla_powers.cpp

@@ -74,6 +74,7 @@ am().set(yval, am_value(y));
 am().set(rval, am_value(r));
 
 --*/
+// clang-format off
 #include "math/lp/nla_core.h"
 
 namespace nla {

src/math/lp/nla_powers.h

@@ -13,6 +13,7 @@ Description:
   Refines bounds on powers.
   
 --*/
+// clang-format off
 
 #include "math/lp/nla_types.h"
 

src/math/lp/nla_settings.h

@@ -6,6 +6,7 @@ Author:
     Lev Nachmanson (levnach)
 
 --*/
+// clang-format off
 
 #pragma once
 namespace nla {

src/math/lp/nla_solver.cpp

@@ -5,6 +5,7 @@
     Lev Nachmanson (levnach)
     Nikolaj Bjorner (nbjorner)
   --*/
+// clang-format off
 #include "math/lp/nla_solver.h"
 #include <map>
 #include "math/lp/monic.h"

src/math/lp/nla_solver.h

@@ -6,6 +6,7 @@ Author:
     Nikolaj Bjorner (nbjorner)
 
 --*/
+// clang-format off
 #pragma once
 #include "util/vector.h"
 #include "math/lp/lp_settings.h"

src/math/lp/nla_tangent_lemmas.cpp

@@ -6,6 +6,7 @@
     Nikolaj Bjorner (nbjorner)
 
   --*/
+// clang-format off
 #include "math/lp/nla_tangent_lemmas.h"
 #include "math/lp/nla_core.h"
 
@@ -75,8 +76,8 @@ private:
         c().negate_relation(lemma, m_jy, m_y.rat_sign()*pl.y);
 #if Z3DEBUG
         SASSERT(c().val(m_x) == m_xy.x && c().val(m_y) == m_xy.y);
-        int mult_sign = nla::rat_sign(pl.x - m_xy.x)*nla::rat_sign(pl.y - m_xy.y);
-        SASSERT((mult_sign == 1) == m_below);
+        // int mult_sign = nla::rat_sign(pl.x - m_xy.x)*nla::rat_sign(pl.y - m_xy.y);
+        SASSERT((nla::rat_sign(pl.x - m_xy.x)*nla::rat_sign(pl.y - m_xy.y) == 1) == m_below);
         // If "mult_sign is 1"  then (a - x)(b-y) > 0 and ab - bx - ay + xy > 0
         // or -ab + bx + ay < xy or -ay - bx + xy > -ab
         // val(j) stands for xy. So, finally we have  -ay - bx + j > - ab

src/math/lp/nla_tangent_lemmas.h

@@ -5,6 +5,7 @@
     Lev Nachmanson (levnach)
     Nikolaj Bjorner (nbjorner)
   --*/
+// clang-format off
 #pragma once
 #include "util/rational.h"
 #include "math/lp/factorization.h"

src/math/lp/nla_types.h

@@ -13,6 +13,7 @@ Description:
   Types used for nla solver.
   
 --*/
+// clang-format off
 
 #pragma once
 

src/math/lp/numeric_pair.h

@@ -17,6 +17,7 @@
 
 
   --*/
+// clang-format off
 #pragma once
 #define lp_for_z3
 #include <string>

src/math/lp/permutation_matrix.cpp

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include <memory>
 #include "util/vector.h"
 #include "math/lp/permutation_matrix_def.h"

src/math/lp/permutation_matrix.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 #include "util/vector.h"
 #include <algorithm>

src/math/lp/permutation_matrix_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include "util/vector.h"

src/math/lp/random_updater.cpp

@@ -17,5 +17,6 @@ Revision History:
 
 
 --*/
+// clang-format off
 #include "math/lp/random_updater_def.h"
 

src/math/lp/random_updater.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 
 #pragma once
 #include <set>

src/math/lp/random_updater_def.h

@@ -17,6 +17,7 @@ Revision History:
 
 
 --*/
+// clang-format off
 #pragma once
 
 #include "math/lp/random_updater.h"