From d823349d2e08c373cdbbc95636787c1cc85a2a87 Mon Sep 17 00:00:00 2001
From: Lev Nachmanson <levnach@microsoft.com>
Date: Mon, 9 May 2016 15:24:47 -0700
Subject: [PATCH] dev(lp): integrate with z3

Signed-off-by: Lev Nachmanson <levnach@microsoft.com>
---
 src/tests/util/lp/lp.cpp                     | 74 ++++++++-------
 src/tests/util/lp/mps_reader.h               | 99 ++++++++++++--------
 src/tests/util/lp/smt_reader.h               | 32 +++----
 src/util/lp/binary_heap_priority_queue.cpp   |  2 +-
 src/util/lp/binary_heap_priority_queue.h     |  2 +
 src/util/lp/binary_heap_upair_queue.cpp      |  5 +-
 src/util/lp/binary_heap_upair_queue.h        |  2 +-
 src/util/lp/canonic_left_side.h              |  4 +-
 src/util/lp/column_info.h                    |  2 +-
 src/util/lp/core_solver_pretty_printer.cpp   | 44 ++++-----
 src/util/lp/indexed_value.h                  |  3 -
 src/util/lp/indexed_vector.h                 | 10 +-
 src/util/lp/indexed_vector_instances.cpp     |  3 +-
 src/util/lp/lar_constraints.h                |  8 +-
 src/util/lp/lar_core_solver.cpp              |  3 +-
 src/util/lp/lar_core_solver.h                |  7 +-
 src/util/lp/lar_solver.cpp                   | 38 ++++----
 src/util/lp/lar_solver.h                     | 21 ++---
 src/util/lp/lar_solver_instances.cpp         |  2 +-
 src/util/lp/lp_core_solver_base.cpp          | 17 +++-
 src/util/lp/lp_core_solver_base.h            |  3 +-
 src/util/lp/lp_dual_core_solver.cpp          | 14 ++-
 src/util/lp/lp_dual_simplex.cpp              | 12 ++-
 src/util/lp/lp_dual_simplex.h                |  3 +-
 src/util/lp/lp_primal_core_solver.cpp        | 37 +++++---
 src/util/lp/lp_primal_core_solver.h          |  7 +-
 src/util/lp/lp_primal_simplex.cpp            | 20 ++--
 src/util/lp/lp_primal_simplex.h              |  3 +-
 src/util/lp/lp_settings.cpp                  | 10 +-
 src/util/lp/lp_settings.h                    | 16 +---
 src/util/lp/lp_solver.cpp                    | 28 +++---
 src/util/lp/lp_solver.h                      |  6 +-
 src/util/lp/lu.cpp                           | 19 +++-
 src/util/lp/lu.h                             |  2 -
 src/util/lp/lu_instances.cpp                 |  1 -
 src/util/lp/matrix.h                         |  2 +-
 src/util/lp/numeric_pair.h                   | 93 +++++++++++++-----
 src/util/lp/permutation_matrix.cpp           |  8 +-
 src/util/lp/permutation_matrix.h             | 10 +-
 src/util/lp/permutation_matrix_instances.cpp |  7 +-
 src/util/lp/row_eta_matrix.cpp               |  2 +-
 src/util/lp/row_eta_matrix.h                 |  8 --
 src/util/lp/scaler.h                         |  2 +-
 src/util/lp/sparse_matrix.cpp                | 40 ++++----
 src/util/lp/sparse_matrix.h                  | 11 +--
 src/util/lp/sparse_matrix_instances.cpp      | 11 +--
 src/util/lp/sparse_vector.h                  | 13 +--
 src/util/lp/square_dense_submatrix.h         |  1 -
 src/util/lp/static_matrix.cpp                | 16 ++--
 src/util/lp/static_matrix.h                  | 16 +---
 src/util/lp/static_matrix_instances.cpp      |  3 +-
 src/util/lp/tail_matrix.h                    |  1 +
 52 files changed, 418 insertions(+), 385 deletions(-)

diff --git a/src/tests/util/lp/lp.cpp b/src/tests/util/lp/lp.cpp
index 6c9fd4f71..f940f5aa0 100644
--- a/src/tests/util/lp/lp.cpp
+++ b/src/tests/util/lp/lp.cpp
@@ -87,7 +87,7 @@ void test_matrix(sparse_matrix<T, X> & a) {
 
 void tst1() {
     std::cout << "testing the minimial matrix with 1 row and 1 column" << std::endl;
-    sparse_matrix<float, float> m0(1);
+    sparse_matrix<double, double> m0(1);
     m0.set(0, 0, 1);
     // print_matrix(m0);
     m0.set(0, 0, 0);
@@ -95,17 +95,17 @@ void tst1() {
     test_matrix(m0);
 
     unsigned rows = 2;
-    sparse_matrix<float, float> m(rows);
+    sparse_matrix<double, double> m(rows);
     std::cout << "setting m(0,1)=" << std::endl;
 
     m.set(0, 1,  11);
     m.set(0, 0,  12);
 
-    // print_matrix<float>(m);
+    // print_matrix<double>(m);
 
     test_matrix(m);
 
-    sparse_matrix<float, float> m1(2);
+    sparse_matrix<double, double> m1(2);
     m1.set(0, 0, 2);
     m1.set(1, 0, 3);
     // print_matrix(m1);
@@ -118,7 +118,7 @@ void tst1() {
 
 
     std::cout << "printing zero matrix 3 by 1" << std::endl;
-    sparse_matrix<float, float> m2(3);
+    sparse_matrix<double, double> m2(3);
     // print_matrix(m2);
 
     m2.set(0, 0, 1);
@@ -128,12 +128,12 @@ void tst1() {
 
     test_matrix(m2);
 
-    sparse_matrix<float, float> m10by9(10);
+    sparse_matrix<double, double> m10by9(10);
     m10by9.set(0, 1, 1);
 
     m10by9(0, 1) = 4;
 
-    float test = m10by9(0, 1);
+    double test = m10by9(0, 1);
 
     std::cout << "got " << test << std::endl;
 
@@ -629,7 +629,7 @@ void fill_matrix(sparse_matrix<T, X>& m){
 }
 
 void test_pivot_like_swaps_and_pivot(){
-    sparse_matrix<float, float> m(10);
+    sparse_matrix<double, double> m(10);
     fill_matrix(m);
     // print_matrix(m);
 // pivot at 2,7
@@ -648,12 +648,12 @@ void test_pivot_like_swaps_and_pivot(){
     m.swap_rows(1, 3);
     // print_matrix(m);
 
-    vector<float> row;
-    float alpha = 2.33;
+    vector<double> row;
+    double alpha = 2.33;
     unsigned pivot_row = 1;
     unsigned target_row = 2;
     unsigned pivot_row_0 = 3;
-    float beta = 3.1;
+    double beta = 3.1;
     m(target_row, 3) = 0;
     m(target_row, 5) = 0;
     m(pivot_row, 6) = 0;
@@ -680,7 +680,7 @@ void test_pivot_like_swaps_and_pivot(){
 
 #ifdef LEAN_DEBUG
 void test_swap_rows() {
-    sparse_matrix<float, float> m(10);
+    sparse_matrix<double, double> m(10);
     fill_matrix(m);
     // print_matrix(m);
     test_swap_rows(m, 3, 5);
@@ -701,7 +701,7 @@ void test_swap_rows() {
     test_swap_rows(m, 0, 7);
 
 // go over some corner cases
-    sparse_matrix<float, float> m0(2);
+    sparse_matrix<double, double> m0(2);
     test_swap_rows(m0, 0, 1);
     m0(0, 0) = 3;
     test_swap_rows(m0, 0, 1);
@@ -709,7 +709,7 @@ void test_swap_rows() {
     test_swap_rows(m0, 0, 1);
 
 
-    sparse_matrix<float, float> m1(10);
+    sparse_matrix<double, double> m1(10);
     test_swap_rows(m1, 0, 1);
     m1(0, 0) = 3;
     test_swap_rows(m1, 0, 1);
@@ -721,7 +721,7 @@ void test_swap_rows() {
 
     test_swap_rows(m1, 0, 1);
 
-    sparse_matrix<float, float> m2(3);
+    sparse_matrix<double, double> m2(3);
     test_swap_rows(m2, 0, 1);
     m2(0, 0) = 3;
     test_swap_rows(m2, 0, 1);
@@ -833,7 +833,7 @@ void sparse_matrix_with_permutaions_test() {
 }
 
 void test_swap_columns() {
-    sparse_matrix<float, float> m(10);
+    sparse_matrix<double, double> m(10);
     fill_matrix(m);
     // print_matrix(m);
 
@@ -853,7 +853,7 @@ void test_swap_columns() {
     test_swap_columns(m, 0, 7);
 
 // go over some corner cases
-    sparse_matrix<float, float> m0(2);
+    sparse_matrix<double, double> m0(2);
     test_swap_columns(m0, 0, 1);
     m0(0, 0) = 3;
     test_swap_columns(m0, 0, 1);
@@ -861,7 +861,7 @@ void test_swap_columns() {
     test_swap_columns(m0, 0, 1);
 
 
-    sparse_matrix<float, float> m1(10);
+    sparse_matrix<double, double> m1(10);
     test_swap_columns(m1, 0, 1);
     m1(0, 0) = 3;
     test_swap_columns(m1, 0, 1);
@@ -873,7 +873,7 @@ void test_swap_columns() {
 
     test_swap_columns(m1, 0, 1);
 
-    sparse_matrix<float, float> m2(3);
+    sparse_matrix<double, double> m2(3);
     test_swap_columns(m2, 0, 1);
     m2(0, 0) = 3;
     test_swap_columns(m2, 0, 1);
@@ -1015,6 +1015,13 @@ void update_settings(argument_parser & args_parser, lp_settings& settings) {
         cout << "setting partial pivot constant to " << n << std::endl;
         settings.c_partial_pivoting = n;
     }
+#if LEAN_DEBUG
+    if (get_int_from_args_parser("--ddd", args_parser, n)) {
+        cout << "setting ddd to " << n << std::endl;
+        settings.ddd = n;
+    }
+#endif
+
     if (get_int_from_args_parser("--density", args_parser, n)) {
         double density = static_cast<double>(n) / 100.0;
         cout << "setting density to " << density << std::endl;
@@ -1029,8 +1036,8 @@ void update_settings(argument_parser & args_parser, lp_settings& settings) {
     }
 }
 
-
-void setup_solver(unsigned max_iterations, unsigned time_limit, bool look_for_min, argument_parser & args_parser, lp_solver<double, double> * solver) {
+template <typename T, typename X>
+void setup_solver(unsigned max_iterations, unsigned time_limit, bool look_for_min, argument_parser & args_parser, lp_solver<T, X> * solver) {
     if (max_iterations > 0)
         solver->set_max_iterations_per_stage(max_iterations);
 
@@ -1120,7 +1127,7 @@ void solve_mps_double(std::string file_name, bool look_for_min, unsigned max_ite
     delete solver;
 }
 
-void solve_mps_rational(std::string file_name, bool look_for_min, unsigned max_iterations, unsigned time_limit, bool dual, argument_parser & /*args_parser*/) {
+void solve_mps_rational(std::string file_name, bool look_for_min, unsigned max_iterations, unsigned time_limit, bool dual, argument_parser & args_parser) {
     mps_reader<mpq, mpq> reader(file_name);
     reader.read();
     if (reader.is_ok()) {
@@ -1128,14 +1135,8 @@ void solve_mps_rational(std::string file_name, bool look_for_min, unsigned max_i
         if (look_for_min) {
             solver->flip_costs();
         }
+        setup_solver(max_iterations, time_limit, look_for_min, args_parser, solver);
         int begin = get_millisecond_count();
-        if (max_iterations > 0) {
-            solver->set_max_iterations_per_stage(max_iterations);
-        }
-
-        if (time_limit > 0) {
-            solver->set_time_limit(time_limit);
-        }
         solver->find_maximal_solution();
         std::cout << "Status: " << lp_status_to_string(solver->get_status()) << std::endl;
 
@@ -1756,7 +1757,7 @@ void test_init_U() {
 }
 
 void test_replace_column() {
-    sparse_matrix<float, float> m(10);
+    sparse_matrix<double, double> m(10);
     fill_matrix(m);
     m.swap_columns(0, 7);
     m.swap_columns(6, 3);
@@ -1765,7 +1766,7 @@ void test_replace_column() {
         m(i, 0) = 0;
     }
 
-    indexed_vector<float> w(m.dimension());
+    indexed_vector<double> w(m.dimension());
     for (unsigned i = 0; i < m.dimension(); i++) {
         w.set_value(i % 3, i);
     }
@@ -1783,6 +1784,7 @@ void test_replace_column() {
 
 void setup_args_parser(argument_parser & parser) {
     parser.add_option_with_after_string_with_help("--density", "the percentage of non-zeroes in the matrix below which it is not dense");
+    parser.add_option_with_after_string_with_help("--ddd", "the debug var");
     parser.add_option_with_after_string_with_help("--harris_toler", "harris tolerance");
     parser.add_option_with_help_string("--test_swaps", "test row swaps with a permutation");
     parser.add_option_with_after_string_with_help("--checklu", "the file name for lu checking");
@@ -1862,8 +1864,8 @@ struct mem_cpy_place_holder {
 };
 
 int finalize(unsigned ret) {
-    finalize_util_module();
-    finalize_numerics_module();
+    lean::finalize_util_module();
+    lean::finalize_numerics_module();
     return ret;
 }
 
@@ -2450,8 +2452,8 @@ void test_square_dense_submatrix() {
 }
 
 int main(int argn, char * const * argv) {
-    initialize_util_module();
-    initialize_numerics_module();
+    lean::initialize_util_module();
+    lean::initialize_numerics_module();
     int ret;
     argument_parser args_parser(argn, argv);
     setup_args_parser(args_parser);
@@ -2560,5 +2562,5 @@ int main(int argn, char * const * argv) {
 #endif
     tst1();
     std::cout<< "done with LP tests\n";
-    return finalize(has_violations() ? 1 : 0);
+    return finalize(lean::has_violations() ? 1 : 0);
 }
diff --git a/src/tests/util/lp/mps_reader.h b/src/tests/util/lp/mps_reader.h
index ebd290588..5591fd09b 100644
--- a/src/tests/util/lp/mps_reader.h
+++ b/src/tests/util/lp/mps_reader.h
@@ -8,22 +8,42 @@
 #pragma once
 
 // reads an MPS file reperesenting a Mixed Integer Program
+#include <functional>
+#include <algorithm>
 #include <string>
 #include <vector>
 #include <unordered_map>
+#include <iostream>
+#include <fstream>
+#include <locale>
 #include "util/lp/lp_primal_simplex.h"
 #include "util/lp/lp_dual_simplex.h"
 #include "util/lp/lar_solver.h"
-#include <iostream>
-#include <fstream>
-#include <util/numerics/mpq.h>
-#include <functional>
-#include <algorithm>
+#include "util/lp/lp_utils.h"
+#include "util/lp/lp_solver.h"
 namespace lean {
 using namespace std;
+inline bool my_white_space(const char & a) {
+    return a == ' ' || a == '\t';
+}
+inline size_t number_of_whites(const string & s)  {
+    size_t i = 0;
+    for (; i < s.size(); i++)
+        if (!my_white_space(s[i])) return i;
+    return i;
+}
+inline size_t number_of_whites_from_end(const string & s)  {
+    size_t ret = 0;
+    for (int i = static_cast<int>(s.size()) - 1; i >= 0; i--)
+        if (my_white_space(s[i])) ret++; else break;
+
+    return ret;
+}
+
+
     // trim from start
-inline string &ltrim(std::string &s) {
-    s.erase(s.begin(), find_if(s.begin(), s.end(), std::not1(std::ptr_fun<int, int>(std::isspace))));
+inline std::string &ltrim(std::string &s) {
+    s.erase(0, number_of_whites(s));
     return s;
 }
 
@@ -32,7 +52,8 @@ inline string &ltrim(std::string &s) {
 
     // trim from end
 inline string &rtrim(std::string &s) {
-    s.erase(std::find_if(s.rbegin(), s.rend(), std::not1(std::ptr_fun<int, int>(std::isspace))).base(), s.end());
+    //       s.erase(std::find_if(s.rbegin(), s.rend(), std::not1(std::ptr_fun<int, int>(std::isspace))).base(), s.end());
+    s.erase(s.end() - number_of_whites_from_end(s), s.end());
     return s;
 }
     // trim from both ends
@@ -46,8 +67,8 @@ inline string trim(std::string const &r) {
 }
 
 
-inline vector<string> string_split(const string &source, const char *delimiter, bool keep_empty)  {
-    vector<string> results;
+inline std::vector<string> string_split(const string &source, const char *delimiter, bool keep_empty)  {
+    std::vector<string> results;
     size_t prev = 0;
     size_t next = 0;
     while ((next = source.find_first_of(delimiter, prev)) != std::string::npos) {
@@ -62,9 +83,9 @@ inline vector<string> string_split(const string &source, const char *delimiter,
     return results;
 }
 
-inline vector<string> split_and_trim(string line) {
+inline std::vector<string> split_and_trim(string line) {
     auto split = string_split(line, " \t", false);
-    vector<string> ret;
+    std::vector<string> ret;
     for (auto s : split) {
         ret.push_back(trim(s));
     }
@@ -220,7 +241,7 @@ class mps_reader {
                 goto fail;
             } else {
                 row * row = t->second;
-                row->m_row_columns[column_name] = atof(column_data.substr(8).c_str());
+                row->m_row_columns[column_name] = from_string<T>(column_data.substr(8));
                 if (column_data.size() > 24) {
                     column_data = column_data.substr(25);
                     if (column_data.size() >= 22) {
@@ -248,7 +269,7 @@ class mps_reader {
                 return;
             }
             row *r = t->second;
-            r->m_row_columns[column_name] = atof(tokens[i + 1].c_str());
+            r->m_row_columns[column_name] = from_string<T>(tokens[i + 1]);
         }
     }
 
@@ -307,7 +328,7 @@ class mps_reader {
                 goto fail;
             } else {
                 row * row = t->second;
-                row->m_right_side = atof(rhsides.substr(8).c_str());
+                row->m_right_side = from_string<T>(rhsides.substr(8));
                 if (rhsides.size() > 24) {
                     rhsides = rhsides.substr(25);
                     if (rhsides.size() >= 22) {
@@ -333,7 +354,7 @@ class mps_reader {
             return;
         }
         string rhsides = m_line.substr(14);
-        vector<string> splitted_line = split_and_trim(rhsides);
+        std::vector<string> splitted_line = split_and_trim(rhsides);
 
         for (unsigned i = 0; i < splitted_line.size() - 1; i += 2) {
             auto t = m_rows.find(splitted_line[i]);
@@ -342,7 +363,7 @@ class mps_reader {
                 return;
             }
             row * row = t->second;
-            row->m_right_side = atof(splitted_line[i + 1].c_str());
+            row->m_right_side = from_string<X>(splitted_line[i + 1]);
         }
     }
 
@@ -393,7 +414,7 @@ class mps_reader {
                 cout << "cannot find " << column_name << endl;
                 goto fail;
             } else {
-                vector<string> bound_string;
+                std::vector<string> bound_string;
                 bound_string.push_back(column_name);
                 if (colstr.size() > 14) {
                     bound_string.push_back(colstr.substr(14));
@@ -414,7 +435,7 @@ class mps_reader {
         }
     }
 
-    void update_bound(bound * b, vector<string> bound_string) {
+    void update_bound(bound * b, std::vector<string> bound_string) {
         /*
           UP means an upper bound is applied to the variable. A bound of type LO means a lower bound is applied. A bound type of FX ("fixed") means that the variable has upper and lower bounds equal to a single value. A bound type of FR ("free") means the variable has neither lower nor upper bounds and so can take on negative values. A variation on that is MI for free negative, giving an upper bound of 0 but no lower bound. Bound type PL is for a free positive for zero to plus infinity, but as this is the normal default, it is seldom used. There are also bound types for use in MIP models - BV for binary, being 0 or 1. UI for upper integer and LI for lower integer. SC stands for semi-continuous and indicates that the variable may be zero, but if not must be equal to at least the value given.
         */
@@ -431,19 +452,16 @@ class mps_reader {
                 set_m_ok_to_false();
                 return;
             }
-            double val = atof(bound_string[1].c_str());
-
             b->m_upper_is_set = true;
-            b->m_upper= val;
+            b->m_upper= from_string<T>(bound_string[1]);
         } else if (bound_type == "LO" || bound_type == "LI") {
             if (bound_string.size() <= 1){
                 set_m_ok_to_false();
                 return;
             }
 
-            double val = atof(bound_string[1].c_str());
             b->m_low_is_set = true;
-            b->m_low = val;
+            b->m_low = from_string<T>(bound_string[1]);
         } else if (bound_type == "FR") {
             b->m_free = true;
         } else if (bound_type == "FX") {
@@ -452,9 +470,8 @@ class mps_reader {
                 return;
             }
 
-            double val = atof(bound_string[1].c_str());
             b->m_value_is_fixed = true;
-            b->m_fixed_value = val;
+            b->m_fixed_value =  from_string<T>(bound_string[1]);
         } else if (bound_type == "PL") {
             b->m_low_is_set = true;
             b->m_low = 0;
@@ -471,7 +488,7 @@ class mps_reader {
     void create_or_update_bound() {
         const unsigned name_offset = 14;
         lean_assert(m_line.size() >= 14);
-        vector<string> bound_string = split_and_trim(m_line.substr(name_offset, m_line.size()));
+        std::vector<string> bound_string = split_and_trim(m_line.substr(name_offset, m_line.size()));
 
         if (bound_string.size() == 0) {
             set_m_ok_to_false();
@@ -514,7 +531,7 @@ class mps_reader {
                 goto fail;
             } else {
                 row * row = t->second;
-                row->m_range = atof(rhsides.substr(8).c_str());
+                row->m_range = from_string<T>(rhsides.substr(8));
                 maybe_modify_current_row_and_add_row_for_range(row);
                 if (rhsides.size() > 24) {
                     rhsides = rhsides.substr(25);
@@ -533,7 +550,7 @@ class mps_reader {
     }
 
 
-    void read_range(vector<string> & splitted_line){
+    void read_range(std::vector<string> & splitted_line){
         for (unsigned i = 1; i < splitted_line.size() - 1; i += 2) {
             auto it = m_rows.find(splitted_line[i]);
             if (it == m_rows.end()) {
@@ -541,7 +558,7 @@ class mps_reader {
                 return;
             }
             row * row = it->second;
-            row->m_range = atof(splitted_line[i + 1].c_str());
+            row->m_range = from_string<X>(splitted_line[i + 1]);
             maybe_modify_current_row_and_add_row_for_range(row);
         }
     }
@@ -632,13 +649,13 @@ class mps_reader {
       where |range| is range's absolute value.
     */
 
-    lp_relation get_relation_from_row(mps_reader::row_type row_type) {
-        switch (row_type) {
-        case mps_reader::row_type::Less_or_equal: return lp_relation::Less_or_equal;
-        case mps_reader::row_type::Greater_or_equal: return lp_relation::Greater_or_equal;
-        case mps_reader::row_type::Equal: return lp_relation::Equal;
+    lp_relation get_relation_from_row(row_type rt) {
+        switch (rt) {
+        case mps_reader::Less_or_equal: return lp_relation::Less_or_equal;
+        case mps_reader::Greater_or_equal: return lp_relation::Greater_or_equal;
+        case mps_reader::Equal: return lp_relation::Equal;
         default:
-            std::cout << "Unexpected row_type " << row_type << endl;
+            std::cout << "Unexpected rt " << rt << endl;
             set_m_ok_to_false();
             throw;
         }
@@ -694,8 +711,8 @@ class mps_reader {
     }
 
 public:
-    vector<string> column_names() {
-        vector<string> v;
+    std::vector<string> column_names() {
+        std::vector<string> v;
         for (auto s : m_columns) {
             v.push_back(s.first);
         }
@@ -748,13 +765,13 @@ public:
         return solver;
     }
 
-    lconstraint_kind get_lar_relation_from_row(mps_reader::row_type row_type) {
-        switch (row_type) {
+    lconstraint_kind get_lar_relation_from_row(row_type rt) {
+        switch (rt) {
         case Less_or_equal: return LE;
         case Greater_or_equal: return GE;
         case Equal: return EQ;
         default:
-            std::cout << "Unexpected row_type " << row_type << endl;
+            std::cout << "Unexpected rt " << rt << endl;
             set_m_ok_to_false();
             throw;
         }
diff --git a/src/tests/util/lp/smt_reader.h b/src/tests/util/lp/smt_reader.h
index 87b06b204..9eb691242 100644
--- a/src/tests/util/lp/smt_reader.h
+++ b/src/tests/util/lp/smt_reader.h
@@ -16,29 +16,20 @@
 #include "util/lp/lar_solver.h"
 #include <iostream>
 #include <fstream>
-#include <util/numerics/mpq.h>
 #include <functional>
 #include <algorithm>
 #include "tests/util/lp/mps_reader.h"
+#include "util/lp/canonic_left_side.h"
 #include "util/lp/lar_constraints.h"
 #include <sstream>
 #include <cstdlib>
 namespace lean {
-    using namespace std;
-
-    template<typename T>
-    T from_string(const std::string& str) {
-        std::istringstream ss(str);
-        T ret;
-        ss >> ret;
-        return ret;
-    }
 
     class smt_reader {
     public:
         struct lisp_elem {
             string m_head;
-            vector<lisp_elem> m_elems;
+            std::vector<lisp_elem> m_elems;
             void print() {
                 if (m_elems.size()) {
                     cout << '(';
@@ -51,12 +42,12 @@ namespace lean {
                     cout << " " << m_head;
                 }
             }
-            unsigned size() const { return m_elems.size(); }
+            unsigned size() const { return static_cast<unsigned>(m_elems.size()); }
             bool is_simple() const { return size() == 0; }
         };
         struct formula_constraint {
             lconstraint_kind m_kind;
-            vector<pair<mpq, string>> m_coeffs;
+            std::vector<pair<mpq, string>> m_coeffs;
             mpq m_right_side = numeric_traits<mpq>::zero();
             void add_pair(mpq c, string name) {
                 m_coeffs.push_back(make_pair(c, name));
@@ -65,13 +56,12 @@ namespace lean {
 
         lisp_elem m_formula_lisp_elem;
 
-        vector<formula_constraint> m_constraints;
+        std::vector<formula_constraint> m_constraints;
         string m_file_name;
         ifstream m_file_stream;
         string m_line;
         bool m_is_OK = true;
         unsigned m_line_number = 0;
-
         smt_reader(string file_name):
             m_file_name(file_name), m_file_stream(file_name) {
         }
@@ -90,9 +80,9 @@ namespace lean {
         }
 
         int first_separator() {
-            unsigned blank_pos = m_line.find(' ');
-            unsigned br_pos = m_line.find('(');
-            unsigned reverse_br_pos = m_line.find(')');
+            unsigned blank_pos = static_cast<unsigned>(m_line.find(' '));
+            unsigned br_pos = static_cast<unsigned>(m_line.find('('));
+            unsigned reverse_br_pos = static_cast<unsigned>(m_line.find(')'));
             return min(blank_pos, min(br_pos, reverse_br_pos));
         }
 
@@ -150,7 +140,7 @@ namespace lean {
 
         void parse_line() {
             if (m_line.find(":formula") == 0) {
-                int first_br = m_line.find('(');
+                int first_br = static_cast<int>(m_line.find('('));
                 if (first_br == -1) {
                     cout << "empty formula" << endl;
                     return;
@@ -214,7 +204,7 @@ namespace lean {
         }
 
 
-        void add_mult_elem(formula_constraint & c, vector<lisp_elem> & els) {
+        void add_mult_elem(formula_constraint & c, std::vector<lisp_elem> & els) {
             lean_assert(els.size() == 2);
             mpq coeff = get_coeff(els[0]);
             string col_name = get_name(els[1]);
@@ -249,7 +239,7 @@ namespace lean {
             }
         }
 
-        void add_sum(formula_constraint & c, vector<lisp_elem> & sum_els) {
+        void add_sum(formula_constraint & c, std::vector<lisp_elem> & sum_els) {
             for (auto & el : sum_els)
                 add_sum_elem(c, el);
         }
diff --git a/src/util/lp/binary_heap_priority_queue.cpp b/src/util/lp/binary_heap_priority_queue.cpp
index 9d18c1198..37e26e202 100644
--- a/src/util/lp/binary_heap_priority_queue.cpp
+++ b/src/util/lp/binary_heap_priority_queue.cpp
@@ -169,7 +169,7 @@ template <typename T> void binary_heap_priority_queue<T>::put_the_last_at_the_to
 }
 /// return the first element of the queue and removes it from the queue
 template <typename T> unsigned binary_heap_priority_queue<T>::dequeue() {
-    lean_assert(m_heap_size);
+    lean_assert(m_heap_size > 0);
     int ret = m_heap[1];
     put_the_last_at_the_top_and_fix_the_heap();
     m_heap_inverse[ret] = -1;
diff --git a/src/util/lp/binary_heap_priority_queue.h b/src/util/lp/binary_heap_priority_queue.h
index 0ff18fcaf..47350bd73 100644
--- a/src/util/lp/binary_heap_priority_queue.h
+++ b/src/util/lp/binary_heap_priority_queue.h
@@ -1,3 +1,4 @@
+
 /*
   Copyright (c) 2013 Microsoft Corporation. All rights reserved.
   Released under Apache 2.0 license as described in the file LICENSE.
@@ -7,6 +8,7 @@
 #pragma once
 #include <vector>
 #include "util/debug.h"
+#include "util/lp/lp_utils.h"
 namespace lean {
 // the elements with the smallest priority are dequeued first
 template <typename T>
diff --git a/src/util/lp/binary_heap_upair_queue.cpp b/src/util/lp/binary_heap_upair_queue.cpp
index a2d22fa04..b0b3d6e81 100644
--- a/src/util/lp/binary_heap_upair_queue.cpp
+++ b/src/util/lp/binary_heap_upair_queue.cpp
@@ -6,10 +6,11 @@
 */
 
 #include <set>
+#include "util/lp/lp_utils.h"
 #include "util/lp/binary_heap_upair_queue.h"
 namespace lean {
 template <typename T> binary_heap_upair_queue<T>::binary_heap_upair_queue(unsigned size) : m_q(size), m_pairs(size) {
-    lean_assert(size);
+    lean_assert(size > 0);
     for (unsigned i = 0; i < size; i++)
         m_available_spots.push_back(i);
 }
@@ -49,7 +50,7 @@ template <typename T> void binary_heap_upair_queue<T>::enqueue(unsigned i, unsig
         // it is a new pair, let us find a spot for it
         if (m_available_spots.empty()) {
             // we ran out of empty spots
-            unsigned size_was = m_pairs.size();
+            unsigned size_was = static_cast<unsigned>(m_pairs.size());
             unsigned new_size = size_was << 1;
             for (unsigned i = size_was; i < new_size; i++)
                 m_available_spots.push_back(i);
diff --git a/src/util/lp/binary_heap_upair_queue.h b/src/util/lp/binary_heap_upair_queue.h
index 0571944eb..ddb94d638 100644
--- a/src/util/lp/binary_heap_upair_queue.h
+++ b/src/util/lp/binary_heap_upair_queue.h
@@ -13,7 +13,7 @@
 #include <set>
 #include <utility>
 #include "util/lp/binary_heap_priority_queue.h"
-#include "util/lp/hash_helper.h"
+
 
 typedef std::pair<unsigned, unsigned> upair;
 
diff --git a/src/util/lp/canonic_left_side.h b/src/util/lp/canonic_left_side.h
index ac5a4a559..64b5c3cfe 100644
--- a/src/util/lp/canonic_left_side.h
+++ b/src/util/lp/canonic_left_side.h
@@ -10,9 +10,7 @@
 #include <string>
 #include <algorithm>
 #include <utility>
-#include "util/numerics/mpq.h"
 #include "util/lp/column_info.h"
-#include "util/lp/hash_helper.h"
 namespace lean {
 typedef unsigned var_index;
 typedef unsigned constraint_index;
@@ -48,7 +46,7 @@ public:
         m_column_info.set_name(name);
     }
 
-    unsigned size() const { return m_coeffs.size(); }
+    unsigned size() const { return static_cast<unsigned>(m_coeffs.size()); }
 
     void normalize() {
         if (m_coeffs.size() == 0) return;
diff --git a/src/util/lp/column_info.h b/src/util/lp/column_info.h
index 3e31d6c6f..b7435284d 100644
--- a/src/util/lp/column_info.h
+++ b/src/util/lp/column_info.h
@@ -100,7 +100,7 @@ public:
     }
 
     T get_fixed_value() const {
-        lean_assert(m_is_fixed)
+        lean_assert(m_is_fixed);
         return m_fixed_value;
     }
 
diff --git a/src/util/lp/core_solver_pretty_printer.cpp b/src/util/lp/core_solver_pretty_printer.cpp
index 0d535d673..a40e53ffa 100644
--- a/src/util/lp/core_solver_pretty_printer.cpp
+++ b/src/util/lp/core_solver_pretty_printer.cpp
@@ -7,8 +7,10 @@
 #include <limits>
 #include <string>
 #include <algorithm>
+#include "util/lp/lp_utils.h"
 #include "util/lp/lp_core_solver_base.h"
 #include "util/lp/core_solver_pretty_printer.h"
+#include "util/lp/numeric_pair.h"
 namespace lean {
 
 
@@ -32,7 +34,7 @@ core_solver_pretty_printer<T, X>::core_solver_pretty_printer(lp_core_solver_base
     m_cost_title = "costs";
     m_basis_heading_title = "heading";
     m_x_title = "x*";
-    m_title_width = std::max(std::max(m_cost_title.size(), std::max(m_basis_heading_title.size(), m_x_title.size())), m_approx_norm_title.size());
+    m_title_width = static_cast<unsigned>(std::max(std::max(m_cost_title.size(), std::max(m_basis_heading_title.size(), m_x_title.size())), m_approx_norm_title.size()));
 }
 
 template <typename T, typename X> void core_solver_pretty_printer<T, X>::init_costs() {
@@ -52,9 +54,9 @@ template <typename T, typename X> core_solver_pretty_printer<T, X>::~core_solver
     delete [] m_ed_buff;
 }
 template <typename T, typename X> void core_solver_pretty_printer<T, X>::init_rs_width() {
-    m_rs_width = T_to_string(m_core_solver.get_cost()).size();
+    m_rs_width = static_cast<unsigned>(T_to_string(m_core_solver.get_cost()).size());
     for (unsigned i = 0; i < nrows(); i++) {
-        auto wt = T_to_string(m_rs[i]).size();
+        unsigned wt = static_cast<unsigned>(T_to_string(m_rs[i]).size());
         if (wt > m_rs_width) {
             m_rs_width = wt;
         }
@@ -81,7 +83,7 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::init_m_
                       name);
             m_rs[row] += m_core_solver.m_ed[row] * m_core_solver.m_x[column];
         }
-        m_exact_column_norms.push_back(current_column_norm() + 1);
+        m_exact_column_norms.push_back(current_column_norm() + T(1)); // a conversion missing 1 -> T
     }
 }
 
@@ -122,11 +124,11 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::adjust_
 
 
 template <typename T, typename X> unsigned core_solver_pretty_printer<T, X>:: get_column_width(unsigned column) {
-    unsigned w = std::max(m_costs[column].size(), T_to_string(m_core_solver.m_x[column]).size());
+    unsigned w = static_cast<unsigned>(std::max(m_costs[column].size(), T_to_string(m_core_solver.m_x[column]).size()));
     adjust_width_with_bounds(column, w);
     adjust_width_with_basis_heading(column, w);
     for (unsigned i = 0; i < nrows(); i++) {
-        unsigned cellw =  m_A[i][column].size();
+        unsigned cellw = static_cast<unsigned>(m_A[i][column].size());
         if (cellw > w) {
             w = cellw;
         }
@@ -171,14 +173,14 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_x
         return;
     }
 
-    int blanks = m_title_width + 1 - m_x_title.size();
+    int blanks = m_title_width + 1 - static_cast<int>(m_x_title.size());
     m_out << m_x_title;
     print_blanks(blanks, m_out);
 
     auto bh = m_core_solver.m_x;
     for (unsigned i = 0; i < ncols(); i++) {
         string s = T_to_string(bh[i]);
-        int blanks = m_column_widths[i] - s.size();
+        int blanks = m_column_widths[i] - static_cast<int>(s.size());
         print_blanks(blanks, m_out);
         m_out << s << "   "; // the column interval
     }
@@ -217,13 +219,13 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_l
     if (ncols() == 0) {
         return;
     }
-    int blanks = m_title_width + 1 - m_low_bounds_title.size();
+    int blanks = m_title_width + 1 - static_cast<unsigned>(m_low_bounds_title.size());
     m_out << m_low_bounds_title;
     print_blanks(blanks, m_out);
 
     for (unsigned i = 0; i < ncols(); i++) {
         string s = get_low_bound_string(i);
-        int blanks = m_column_widths[i] - s.size();
+        int blanks = m_column_widths[i] - static_cast<unsigned>(s.size());
         print_blanks(blanks, m_out);
         m_out << s << "   "; // the column interval
     }
@@ -234,13 +236,13 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_u
     if (ncols() == 0) {
         return;
     }
-    int blanks = m_title_width + 1 - m_upp_bounds_title.size();
+    int blanks = m_title_width + 1 - static_cast<unsigned>(m_upp_bounds_title.size());
     m_out << m_upp_bounds_title;
     print_blanks(blanks, m_out);
 
     for (unsigned i = 0; i < ncols(); i++) {
         string s = get_upp_bound_string(i);
-        int blanks = m_column_widths[i] - s.size();
+        int blanks = m_column_widths[i] - static_cast<unsigned>(s.size());
         print_blanks(blanks, m_out);
         m_out << s << "   "; // the column interval
     }
@@ -248,12 +250,12 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_u
 }
 
 template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_exact_norms() {
-    int blanks = m_title_width + 1 - m_exact_norm_title.size();
+    int blanks = m_title_width + 1 - static_cast<int>(m_exact_norm_title.size());
     m_out << m_exact_norm_title;
     print_blanks(blanks, m_out);
     for (unsigned i = 0; i < ncols(); i++) {
         string s = get_exact_column_norm_string(i);
-        int blanks = m_column_widths[i] - s.size();
+        int blanks = m_column_widths[i] - static_cast<int>(s.size());
         print_blanks(blanks, m_out);
         m_out << s << "   ";
     }
@@ -261,12 +263,12 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_e
 }
 
 template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_approx_norms() {
-    int blanks = m_title_width + 1 - m_approx_norm_title.size();
+    int blanks = m_title_width + 1 - static_cast<int>(m_approx_norm_title.size());
     m_out << m_approx_norm_title;
     print_blanks(blanks, m_out);
     for (unsigned i = 0; i < ncols(); i++) {
         string s = T_to_string(m_core_solver.m_column_norms[i]);
-        int blanks = m_column_widths[i] - s.size();
+        int blanks = m_column_widths[i] - static_cast<int>(s.size());
         print_blanks(blanks, m_out);
         m_out << s << "   ";
     }
@@ -289,7 +291,7 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print()
 }
 
 template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_basis_heading() {
-    int blanks = m_title_width + 1 - m_basis_heading_title.size();
+    int blanks = m_title_width + 1 - static_cast<int>(m_basis_heading_title.size());
     m_out << m_basis_heading_title;
     print_blanks(blanks, m_out);
 
@@ -299,7 +301,7 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_b
     auto bh = m_core_solver.m_basis_heading;
     for (unsigned i = 0; i < ncols(); i++) {
         string s = T_to_string(bh[i]);
-        int blanks = m_column_widths[i] - s.size();
+        int blanks = m_column_widths[i] - static_cast<unsigned>(s.size());
         print_blanks(blanks, m_out);
         m_out << s << "   "; // the column interval
     }
@@ -307,7 +309,7 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_b
 }
 
 template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_cost() {
-    int blanks = m_title_width + 1 - m_cost_title.size();
+    int blanks = m_title_width + 1 - static_cast<int>(m_cost_title.size());
     m_out << m_cost_title;
     print_blanks(blanks, m_out);
     print_given_rows(m_costs, m_cost_signs, m_core_solver.get_cost());
@@ -317,7 +319,7 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_g
     for (unsigned col = 0; col < row.size(); col++) {
         unsigned width = m_column_widths[col];
         string s = row[col];
-        int number_of_blanks = width - s.size();
+        int number_of_blanks = width - static_cast<unsigned>(s.size());
         lean_assert(number_of_blanks >= 0);
         print_blanks(number_of_blanks, m_out);
         m_out << s << ' ';
@@ -328,7 +330,7 @@ template <typename T, typename X> void core_solver_pretty_printer<T, X>::print_g
     m_out << '=';
 
     string rs = T_to_string(rst);
-    int nb = m_rs_width - rs.size();
+    int nb = m_rs_width - static_cast<int>(rs.size());
     lean_assert(nb >= 0);
     print_blanks(nb + 1, m_out);
     m_out << rs << std::endl;
diff --git a/src/util/lp/indexed_value.h b/src/util/lp/indexed_value.h
index 0d30cf512..b88e9fbd5 100644
--- a/src/util/lp/indexed_value.h
+++ b/src/util/lp/indexed_value.h
@@ -5,9 +5,6 @@
   Author: Lev Nachmanson
 */
 #pragma once
-#include "util/numerics/float.h"
-#include "util/numerics/double.h"
-#include "util/numerics/mpq.h"
 
 namespace lean {
 template <typename T>
diff --git a/src/util/lp/indexed_vector.h b/src/util/lp/indexed_vector.h
index 28398b10a..4217ad2e8 100644
--- a/src/util/lp/indexed_vector.h
+++ b/src/util/lp/indexed_vector.h
@@ -8,17 +8,9 @@
 #pragma once
 #include <vector>
 #include "util/debug.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
-#include "util/numerics/mpq.h"
-#include "util/numerics/mpz.h"
-#include "util/numerics/mpbq.h"
-#include "util/numerics/double.h"
-#include "util/numerics/float.h"
-#include "util/numerics/mpfp.h"
 #include <string>
-#include "util/lp/sparse_vector.h"
 #include <iomanip>
+#include "util/lp/lp_utils.h"
 namespace lean {
 
 template <typename T> void print_vector(const std::vector<T> & t, std::ostream & out);
diff --git a/src/util/lp/indexed_vector_instances.cpp b/src/util/lp/indexed_vector_instances.cpp
index 20808d9ed..7ddc82de3 100644
--- a/src/util/lp/indexed_vector_instances.cpp
+++ b/src/util/lp/indexed_vector_instances.cpp
@@ -11,7 +11,6 @@ template void indexed_vector<double>::clear();
 template void indexed_vector<double>::clear_all();
 template void indexed_vector<double>::erase_from_index(unsigned int);
 template void indexed_vector<double>::set_value(double, unsigned int);
-template void indexed_vector<float>::set_value(float, unsigned int);
 template void indexed_vector<mpq>::clear();
 template void indexed_vector<mpq>::clear_all();
 template void indexed_vector<mpq>::erase_from_index(unsigned int);
@@ -19,6 +18,6 @@ template void indexed_vector<mpq>::resize(unsigned int);
 template void indexed_vector<mpq>::set_value(mpq, unsigned int);
 #ifdef LEAN_DEBUG
 template bool indexed_vector<double>::is_OK() const;
-template bool lean::indexed_vector<mpq>::is_OK() const;
+template bool indexed_vector<mpq>::is_OK() const;
 #endif
 }
diff --git a/src/util/lp/lar_constraints.h b/src/util/lp/lar_constraints.h
index 1b0323dee..0852bb6bb 100644
--- a/src/util/lp/lar_constraints.h
+++ b/src/util/lp/lar_constraints.h
@@ -8,13 +8,10 @@
 #pragma once
 #include <vector>
 #include <utility>
-#include "util/debug.h"
-#include "util/buffer.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
 #include <unordered_map>
 #include <string>
 #include <algorithm>
+#include "util/lp/lp_utils.h"
 #include "util/lp/canonic_left_side.h"
 namespace lean {
 inline lconstraint_kind flip_kind(lconstraint_kind t) {
@@ -30,6 +27,7 @@ inline std::string lconstraint_kind_string(lconstraint_kind t) {
     case EQ: return std::string("=");
     }
     lean_unreachable();
+    return std::string(); // it is unreachable
 }
 
 class lar_base_constraint {
@@ -54,7 +52,7 @@ public:
     lar_constraint(const lar_base_constraint & c);
 
     unsigned size() const {
-        return m_left_side.size();
+        return static_cast<unsigned>(m_left_side.size());
     }
 
     buffer<std::pair<mpq, var_index>> get_left_side_coefficients() const;
diff --git a/src/util/lp/lar_core_solver.cpp b/src/util/lp/lar_core_solver.cpp
index af1a52b10..820f0dec3 100644
--- a/src/util/lp/lar_core_solver.cpp
+++ b/src/util/lp/lar_core_solver.cpp
@@ -647,6 +647,7 @@ template <typename T, typename X>    bool lar_core_solver<T, X>::improves_pivot_
     default:
         lean_unreachable();
     }
+    return false; // it is unreachable
 }
 
 template <typename T, typename X>    int lar_core_solver<T, X>::choose_entering_column_for_row_inf_strategy(int inf_sign) {
@@ -744,7 +745,7 @@ template <typename T, typename X>    void lar_core_solver<T, X>::advance_on_infe
     int entering_delta_sign;
     std::vector<unsigned> leaving_candidates;
     X delta = find_initial_delta_and_its_sign(inf_row, entering, inf_sign, entering_delta_sign, leaving_candidates);
-    lean_assert(leaving_candidates.size());
+    lean_assert(leaving_candidates.size() > 0);
     lean_assert(delta > zero_of_type<X>());
     unsigned row = my_random() % this->m_m;
     unsigned initial_row = row;
diff --git a/src/util/lp/lar_core_solver.h b/src/util/lp/lar_core_solver.h
index da5e8b76c..6222f7710 100644
--- a/src/util/lp/lar_core_solver.h
+++ b/src/util/lp/lar_core_solver.h
@@ -7,6 +7,7 @@
 #pragma once
 #include <vector>
 #include <string>
+#include <utility>
 #include "util/lp/lp_core_solver_base.h"
 #include <algorithm>
 #include "util/lp/indexed_vector.h"
@@ -23,7 +24,7 @@ class lar_core_solver : public lp_core_solver_base<T, X> {
     std::vector<unsigned> m_tight_basic_columns;
     std::vector<breakpoint<X>> m_breakpoints;
     binary_heap_priority_queue<X> m_breakpoint_indices_queue;
-    std::vector<pair<mpq, unsigned>> m_infeasible_row;
+    std::vector<std::pair<mpq, unsigned>> m_infeasible_row;
     int m_infeasible_row_sign = 0;
     // with a breakpoint at this delta
 public:
@@ -38,7 +39,7 @@ public:
 
     int get_infeasible_row_sign() const { return m_infeasible_row_sign;   }
 
-    const std::vector<pair<mpq, unsigned>> & get_infeasibility_info(int & inf_sign) const {
+    const std::vector<std::pair<mpq, unsigned>> & get_infeasibility_info(int & inf_sign) const {
         inf_sign = m_infeasible_row_sign;
         return m_infeasible_row;
     }
@@ -155,7 +156,7 @@ public:
                                   std::vector<unsigned> & leaving_candidates);
 
     unsigned find_leaving_for_inf_row_strategy(std::vector<unsigned> & leaving_candidates) {
-        lean_assert(leaving_candidates.size());
+        lean_assert(leaving_candidates.size() > 0);
         return leaving_candidates[my_random() % leaving_candidates.size()]; // more randomness
     }
 
diff --git a/src/util/lp/lar_solver.cpp b/src/util/lp/lar_solver.cpp
index 1e8f50f01..cec41b6aa 100644
--- a/src/util/lp/lar_solver.cpp
+++ b/src/util/lp/lar_solver.cpp
@@ -7,6 +7,7 @@
 #include <string>
 #include <algorithm>
 #include <vector>
+#include <utility>
 #include "util/lp/lar_solver.h"
 namespace lean {
 double conversion_helper <double>::get_low_bound(const column_info<mpq> & ci) {
@@ -29,7 +30,7 @@ double conversion_helper <double>::get_upper_bound(const column_info<mpq> & ci)
     return ci.get_upper_bound().get_double() - eps;
 }
 
-canonic_left_side * lar_solver::create_or_fetch_existing_left_side(const buffer<pair<mpq, var_index>>& left_side_par) {
+canonic_left_side * lar_solver::create_or_fetch_existing_left_side(const buffer<std::pair<mpq, var_index>>& left_side_par) {
     auto left_side = new canonic_left_side(left_side_par);
     lean_assert(left_side->size() > 0);
     auto it = m_canonic_left_sides.find(left_side);
@@ -60,7 +61,7 @@ mpq lar_solver::find_ratio(canonic_left_side * ls, const lar_constraint & constr
 }
 
 void lar_solver::add_canonic_left_side_for_var(var_index i, std::string var_name) {
-    buffer<pair<mpq, var_index>> b;
+    buffer<std::pair<mpq, var_index>> b;
     b.push_back(std::make_pair(numeric_traits<mpq>::one(), i));
     auto can_ls = new canonic_left_side(b);
     can_ls->set_name(var_name);
@@ -74,7 +75,7 @@ void lar_solver::map_left_side_to_column_of_A(canonic_left_side*  left_side, uns
     lean_assert(m_column_indices_to_canonic_left_sides.find(j) == m_column_indices_to_canonic_left_sides.end());
     left_side->m_column_index = j; // assigning this index does not change the hash of canonic_left_side
     if (left_side->size() > 1) { // if size is one we will not create a row for this left side
-        left_side->m_row_index = m_basis.size();
+        left_side->m_row_index = static_cast<int>(m_basis.size());
         m_basis.push_back(j); // j will be a basis column, so we put it into the basis as well
     }
     m_column_indices_to_canonic_left_sides[j++] = left_side; // pointing from the column to the left side
@@ -105,8 +106,8 @@ void lar_solver::add_row_to_A(static_matrix<U, V> & A, unsigned i, canonic_left_
 
 template <typename U, typename V>
 void lar_solver::create_matrix_A(static_matrix<U, V> & A) {
-    unsigned n = m_column_indices_to_canonic_left_sides.size();
-    unsigned m = m_basis.size();
+    unsigned n = static_cast<unsigned>(m_column_indices_to_canonic_left_sides.size());
+    unsigned m = static_cast<unsigned>(m_basis.size());
     A.init_empty_matrix(m, n);
     unsigned i = 0;
     for (auto & t : m_column_indices_to_canonic_left_sides)
@@ -244,7 +245,7 @@ void lar_solver::fill_column_types() {
 
 template <typename V>
 void lar_solver::fill_bounds_for_core_solver(std::vector<V> & lb, std::vector<V> & ub) {
-    unsigned n = m_canonic_left_sides.size(); // this is the number of columns
+    unsigned n = static_cast<unsigned>(m_canonic_left_sides.size()); // this is the number of columns
     lb.resize(n);
     ub.resize(n);
     for (auto t : m_canonic_left_sides) {
@@ -284,6 +285,7 @@ template <typename V> V lar_solver::get_column_val(std::vector<V> & low_bound, s
     default:
         lean_unreachable();
     }
+    return zero_of_type<V>(); // it is unreachable
 }
 
 lar_solver::~lar_solver() {
@@ -307,8 +309,8 @@ var_index lar_solver::add_var(std::string s) {
     return i;
 }
 
-constraint_index lar_solver::add_constraint(const buffer<pair<mpq, var_index>>& left_side, lconstraint_kind kind_par, mpq right_side_par) {
-    lean_assert(left_side.size());
+constraint_index lar_solver::add_constraint(const buffer<std::pair<mpq, var_index>>& left_side, lconstraint_kind kind_par, mpq right_side_par) {
+    lean_assert(left_side.size() > 0);
     constraint_index i = m_available_constr_index++;
     lar_constraint original_constr(left_side, kind_par, right_side_par, i);
     canonic_left_side * ls = create_or_fetch_existing_left_side(left_side);
@@ -341,6 +343,7 @@ bool lar_solver::constraint_holds(const lar_constraint & constr, std::unordered_
     default:
         lean_unreachable();
     }
+    return false; // it is unreachable
 }
 
 void lar_solver::solve_with_core_solver() {
@@ -352,7 +355,7 @@ void lar_solver::solve_with_core_solver() {
 #endif
 }
 
-bool lar_solver::the_relations_are_of_same_type(const buffer<pair<mpq, unsigned>> & evidence, lconstraint_kind & the_kind_of_sum) {
+bool lar_solver::the_relations_are_of_same_type(const buffer<std::pair<mpq, unsigned>> & evidence, lconstraint_kind & the_kind_of_sum) {
     unsigned n_of_G = 0, n_of_L = 0;
     bool strict = false;
     for (auto & it : evidence) {
@@ -381,7 +384,7 @@ void lar_solver::register_in_map(std::unordered_map<var_index, mpq> & coeffs, la
         else p->second += it.second * a;
     }
 }
-bool lar_solver::the_left_sides_sum_to_zero(const buffer<pair<mpq, unsigned>> & evidence) {
+bool lar_solver::the_left_sides_sum_to_zero(const buffer<std::pair<mpq, unsigned>> & evidence) {
     std::unordered_map<var_index, mpq> coeff_map;
     for (auto & it : evidence) {
         mpq coeff = it.first;
@@ -395,7 +398,7 @@ bool lar_solver::the_left_sides_sum_to_zero(const buffer<pair<mpq, unsigned>> &
     return true;
 }
 
-bool lar_solver::the_righ_sides_do_not_sum_to_zero(const buffer<pair<mpq, unsigned>> & evidence) {
+bool lar_solver::the_righ_sides_do_not_sum_to_zero(const buffer<std::pair<mpq, unsigned>> & evidence) {
     mpq ret = numeric_traits<mpq>::zero();
     for (auto & it : evidence) {
         mpq coeff = it.first;
@@ -407,7 +410,7 @@ bool lar_solver::the_righ_sides_do_not_sum_to_zero(const buffer<pair<mpq, unsign
 }
 #ifdef LEAN_DEBUG
 bool lar_solver::the_evidence_is_correct() {
-    buffer<pair<mpq, unsigned>> evidence;
+    buffer<std::pair<mpq, unsigned>> evidence;
     get_infeasibility_evidence(evidence);
     lconstraint_kind kind;
     lean_assert(the_relations_are_of_same_type(evidence, kind));
@@ -442,7 +445,7 @@ void lar_solver::init_right_sides_with_zeros(std::vector<V> & rs) {
     rs.resize(m_basis.size(), zero_of_type<V>());
 }
 
-mpq lar_solver::sum_of_right_sides_of_evidence(const buffer<pair<mpq, unsigned>> & evidence) {
+mpq lar_solver::sum_of_right_sides_of_evidence(const buffer<std::pair<mpq, unsigned>> & evidence) {
     mpq ret = numeric_traits<mpq>::zero();
     for (auto & it : evidence) {
         mpq coeff = it.first;
@@ -499,7 +502,7 @@ void lar_solver::find_solution_signature_with_doubles(lar_solution_signature & s
         A.clear();
         create_matrix_A(A);
         for (auto & s : column_scale_vector)
-            s = one_of_type<double>();
+            s = 1.0;
     }
     std::vector<double> costs(A.column_count());
     auto core_solver = lp_primal_core_solver<double, double>(A,
@@ -545,7 +548,7 @@ lp_status lar_solver::check() {
     solve();
     return m_status;
 }
-void lar_solver::get_infeasibility_evidence(buffer<pair<mpq, constraint_index>> & evidence){
+void lar_solver::get_infeasibility_evidence(buffer<std::pair<mpq, constraint_index>> & evidence){
     if (!m_mpq_core_solver.get_infeasible_row_sign()) {
         return;
     }
@@ -556,8 +559,8 @@ void lar_solver::get_infeasibility_evidence(buffer<pair<mpq, constraint_index>>
     get_infeasibility_evidence_for_inf_sign(evidence, inf_row, inf_sign);
 }
 
-void lar_solver::get_infeasibility_evidence_for_inf_sign(buffer<pair<mpq, constraint_index>> & evidence,
-                                                         const std::vector<pair<mpq, unsigned>> & inf_row,
+void lar_solver::get_infeasibility_evidence_for_inf_sign(buffer<std::pair<mpq, constraint_index>> & evidence,
+                                                         const std::vector<std::pair<mpq, unsigned>> & inf_row,
                                                          int inf_sign) {
     for (auto & it : inf_row) {
         mpq coeff = it.first;
@@ -727,6 +730,7 @@ mpq lar_solver::get_infeasibility_of_constraint(const lar_normalized_constraint
     default:
         lean_unreachable();
     }
+    return mpq(0); // it is unreachable
 }
 
 mpq lar_solver::get_canonic_left_side_val(canonic_left_side * ls, std::unordered_map<std::string, mpq> & solution) {
diff --git a/src/util/lp/lar_solver.h b/src/util/lp/lar_solver.h
index 4020441e1..279583bd7 100644
--- a/src/util/lp/lar_solver.h
+++ b/src/util/lp/lar_solver.h
@@ -6,10 +6,9 @@
 */
 #pragma once
 #include <vector>
+#include <utility>
 #include "util/debug.h"
 #include "util/buffer.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
 #include <unordered_map>
 #include <unordered_set>
 #include <string>
@@ -61,7 +60,7 @@ class lar_solver {
     std::vector<numeric_pair<mpq>> m_right_side_vector; // this vector will be all zeroes, it might change when the optimization with fixed variables will used
     std::vector<mpq> m_costs;
     canonic_left_side * m_infeasible_canonic_left_side = nullptr; // such can be found at the initialization step
-    canonic_left_side * create_or_fetch_existing_left_side(const buffer<pair<mpq, var_index>>& left_side_par);
+    canonic_left_side * create_or_fetch_existing_left_side(const buffer<std::pair<mpq, var_index>>& left_side_par);
 
 
     bool var_is_fixed(unsigned j);
@@ -148,7 +147,7 @@ public:
 
     var_index add_var(std::string s);
 
-    constraint_index add_constraint(const buffer<pair<mpq, var_index>>& left_side, lconstraint_kind kind_par, mpq right_side_par);
+    constraint_index add_constraint(const buffer<std::pair<mpq, var_index>>& left_side, lconstraint_kind kind_par, mpq right_side_par);
 
     bool is_infeasible(const column_info<mpq> & ci) {
         return ci.low_bound_is_set() && ci.upper_bound_is_set() && ci.get_low_bound() > ci.get_upper_bound();
@@ -162,12 +161,12 @@ public:
 
     void solve_with_core_solver();
 
-    bool the_relations_are_of_same_type(const buffer<pair<mpq, unsigned>> & evidence, lconstraint_kind & the_kind_of_sum);
+    bool the_relations_are_of_same_type(const buffer<std::pair<mpq, unsigned>> & evidence, lconstraint_kind & the_kind_of_sum);
 
     void register_in_map(std::unordered_map<var_index, mpq> & coeffs, lar_constraint & cn, const mpq & a);
-    bool the_left_sides_sum_to_zero(const buffer<pair<mpq, unsigned>> & evidence);
+    bool the_left_sides_sum_to_zero(const buffer<std::pair<mpq, unsigned>> & evidence);
 
-    bool the_righ_sides_do_not_sum_to_zero(const buffer<pair<mpq, unsigned>> & evidence);
+    bool the_righ_sides_do_not_sum_to_zero(const buffer<std::pair<mpq, unsigned>> & evidence);
 
     bool the_evidence_is_correct();
 
@@ -176,7 +175,7 @@ public:
     template <typename V>
     void init_right_sides_with_zeros(std::vector<V> & rs);
 
-    mpq sum_of_right_sides_of_evidence(const buffer<pair<mpq, unsigned>> & evidence);
+    mpq sum_of_right_sides_of_evidence(const buffer<std::pair<mpq, unsigned>> & evidence);
     void prepare_independently_of_numeric_type();
 
     template <typename U, typename V>
@@ -201,10 +200,10 @@ public:
     void solve();
 
     lp_status check();
-    void get_infeasibility_evidence(buffer<pair<mpq, constraint_index>> & evidence);
+    void get_infeasibility_evidence(buffer<std::pair<mpq, constraint_index>> & evidence);
 
-    void get_infeasibility_evidence_for_inf_sign(buffer<pair<mpq, constraint_index>> & evidence,
-                                                 const std::vector<pair<mpq, unsigned>> & inf_row,
+    void get_infeasibility_evidence_for_inf_sign(buffer<std::pair<mpq, constraint_index>> & evidence,
+                                                 const std::vector<std::pair<mpq, unsigned>> & inf_row,
                                                  int inf_sign);
 
 
diff --git a/src/util/lp/lar_solver_instances.cpp b/src/util/lp/lar_solver_instances.cpp
index 84f7f4e06..8d980296c 100644
--- a/src/util/lp/lar_solver_instances.cpp
+++ b/src/util/lp/lar_solver_instances.cpp
@@ -4,5 +4,5 @@
 
   Author: Lev Nachmanson
 */
-#
+
 #include "util/lp/lar_solver.cpp"
diff --git a/src/util/lp/lp_core_solver_base.cpp b/src/util/lp/lp_core_solver_base.cpp
index aee6821b6..078c6cf20 100644
--- a/src/util/lp/lp_core_solver_base.cpp
+++ b/src/util/lp/lp_core_solver_base.cpp
@@ -7,6 +7,7 @@
 #include <set>
 #include <string>
 #include <vector>
+#include "util/lp/lp_utils.h"
 #include "util/lp/lp_core_solver_base.h"
 namespace lean {
 void init_basic_part_of_basis_heading(std::vector<unsigned> & basis, unsigned m, std::vector<int> & basis_heading) {
@@ -21,7 +22,7 @@ void init_non_basic_part_of_basis_heading(std::vector<int> & basis_heading, std:
         if (basis_heading[j] < 0) {
             non_basic_columns.push_back(j);
             // the index of column j in m_non_basic_columns is (- basis_heading[j] - 1)
-            basis_heading[j] = - non_basic_columns.size();
+            basis_heading[j] = - static_cast<int>(non_basic_columns.size());
         }
     }
 }
@@ -116,7 +117,7 @@ solve_yB(std::vector<T> & y) {
 template <typename T, typename X> void lp_core_solver_base<T, X>::
 update_index_of_ed() {
     m_index_of_ed.clear();
-    unsigned i = m_ed.size();
+    unsigned i = static_cast<unsigned>(m_ed.size());
     while (i--) {
         if (!is_zero(m_ed[i]))
             m_index_of_ed.push_back(i);
@@ -248,7 +249,7 @@ calculate_pivot_row_when_pivot_row_of_B1_is_ready() {
         }
     }
 
-    unsigned j = m_pivot_row.size();
+    unsigned j = static_cast<unsigned>(m_pivot_row.size());
     while (j--) {
         if (!is_zero(m_pivot_row[j]))
             m_pivot_row_index.push_back(j);
@@ -274,7 +275,7 @@ print_statistics(X cost) {
 }
 template <typename T, typename X> bool lp_core_solver_base<T, X>::
 print_statistics_with_iterations_and_check_that_the_time_is_over(unsigned total_iterations) {
-    if (total_iterations % m_settings.report_frequency == 0) {
+    if (m_settings.print_statistics && total_iterations % m_settings.report_frequency == 0) {
         std::cout << "iterations = " << total_iterations  <<  ", nonzeros = " << m_factorization->get_number_of_nonzeroes() << std::endl;
         if (time_is_over()) {
             return true;
@@ -356,6 +357,8 @@ column_is_dual_feasible(unsigned j) const {
         std::cout << "unexpected column type = " << column_type_to_string(m_column_type[j]) << std::endl;
         lean_unreachable();
     }
+    lean_unreachable();
+    return false;
 }
 template <typename T, typename X> bool lp_core_solver_base<T, X>::
 d_is_not_negative(unsigned j) const {
@@ -419,7 +422,9 @@ update_basis_and_x(int entering, int leaving, X const & tt) {
                 m_refactor_counter = 0;
                 m_iters_with_no_cost_growing++;
                 if (m_factorization->get_status() != LU_status::OK) {
-                    throw exception(sstream() << "failing refactor on off_result for entering = " << entering << ", leaving = " << leaving << " total_iterations = " << m_total_iterations);
+                    std::stringstream s;
+                    s << "failing refactor on off_result for entering = " << entering << ", leaving = " << leaving << " total_iterations = " << m_total_iterations;
+                    throw_exception(s.str());
                 }
                 return false;
             }
@@ -691,6 +696,8 @@ get_non_basic_column_value_position(unsigned j) {
     default:
         lean_unreachable();
     }
+    lean_unreachable();
+    return at_low_bound;
 }
 
 template <typename T, typename X> void lp_core_solver_base<T, X>::init_lu() {
diff --git a/src/util/lp/lp_core_solver_base.h b/src/util/lp/lp_core_solver_base.h
index 697a797b2..c31c901b0 100644
--- a/src/util/lp/lp_core_solver_base.h
+++ b/src/util/lp/lp_core_solver_base.h
@@ -8,8 +8,7 @@
 #include <set>
 #include <vector>
 #include <string>
-#include "util/sstream.h"
-#include "util/exception.h"
+#include "util/lp/lp_utils.h"
 #include "util/lp/core_solver_pretty_printer.h"
 #include "util/lp/numeric_pair.h"
 #include "util/lp/static_matrix.h"
diff --git a/src/util/lp/lp_dual_core_solver.cpp b/src/util/lp/lp_dual_core_solver.cpp
index cf93ce995..059c62892 100644
--- a/src/util/lp/lp_dual_core_solver.cpp
+++ b/src/util/lp/lp_dual_core_solver.cpp
@@ -55,7 +55,7 @@ template <typename T, typename X> void lp_dual_core_solver<T, X>::fill_non_basis
     while (j--) {
         if (this->m_basis_heading[j] >= 0 || !m_can_enter_basis[j]) continue;
         nb.push_back(j);
-        this->m_basis_heading[j] = - nb.size();
+        this->m_basis_heading[j] = - static_cast<int>(nb.size());
     }
 }
 
@@ -68,7 +68,7 @@ template <typename T, typename X> void lp_dual_core_solver<T, X>::restore_non_ba
         if (m_can_enter_basis[j]) {
             lean_assert(std::find(nb.begin(), nb.end(), j) == nb.end());
             nb.push_back(j);
-            this->m_basis_heading[j] = - nb.size();
+            this->m_basis_heading[j] = - static_cast<int>(nb.size());
         }
     }
 }
@@ -183,6 +183,8 @@ template <typename T, typename X> T lp_dual_core_solver<T, X>::pricing_for_row(u
     default:
         lean_unreachable();
     }
+    lean_unreachable();
+    return numeric_traits<T>::zero();
 }
 
 template <typename T, typename X> void lp_dual_core_solver<T, X>::pricing_loop(unsigned number_of_rows_to_try, unsigned offset_in_rows) {
@@ -283,6 +285,8 @@ template <typename T, typename X> int lp_dual_core_solver<T, X>::define_sign_of_
     default:
         lean_unreachable();
     }
+    lean_unreachable();
+    return 0;
 }
 
 template <typename T, typename X> bool lp_dual_core_solver<T, X>::can_be_breakpoint(unsigned j) {
@@ -329,7 +333,7 @@ template <typename T, typename X> T lp_dual_core_solver<T, X>::get_delta() {
             return this->m_x[m_p] - this->m_low_bound_values[m_p];
         }
         if (this->x_above_upper_bound(m_p)) {
-            return this->m_x[m_p] - this->m_upper_bound_values[m_p];;
+            return this->m_x[m_p] - this->m_upper_bound_values[m_p];
         }
         lean_unreachable();
     case low_bound:
@@ -343,10 +347,12 @@ template <typename T, typename X> T lp_dual_core_solver<T, X>::get_delta() {
         }
         lean_unreachable();
     case fixed:
-        return this->m_x[m_p] - this->m_upper_bound_values[m_p];;
+        return this->m_x[m_p] - this->m_upper_bound_values[m_p];
     default:
         lean_unreachable();
     }
+    lean_unreachable();
+    return zero_of_type<T>();
 }
 
 template <typename T, typename X> void lp_dual_core_solver<T, X>::restore_d() {
diff --git a/src/util/lp/lp_dual_simplex.cpp b/src/util/lp/lp_dual_simplex.cpp
index c16a59a2c..6638bd679 100644
--- a/src/util/lp/lp_dual_simplex.cpp
+++ b/src/util/lp/lp_dual_simplex.cpp
@@ -5,7 +5,7 @@
   Author: Lev Nachmanson
 */
 #include "util/lp/lp_dual_simplex.h"
-namespace lean {
+namespace lean{
 
 template <typename T, typename X> void lp_dual_simplex<T, X>::decide_on_status_after_stage1() {
     switch (m_core_solver->get_status()) {
@@ -193,9 +193,13 @@ template <typename T, typename X> void lp_dual_simplex<T, X>::fill_costs_bounds_
     lean_assert(this->m_columns.find(jj) != this->m_columns.end());
     column_info<T> * ci = this->m_columns[jj];
     switch (ci->get_column_type()) {
-    case upper_bound:
-        throw exception(sstream() << "unexpected bound type " << j << " "
-                        << column_type_to_string(get_column_type(j)));
+    case upper_bound: {
+        std::stringstream s;
+        s << "unexpected bound type " << j << " "
+          << column_type_to_string(get_column_type(j));
+        throw_exception(s.str());
+        break;
+    }
     case low_bound: {
         m_can_enter_basis[j] = true;
         this->set_scaled_cost(j);
diff --git a/src/util/lp/lp_dual_simplex.h b/src/util/lp/lp_dual_simplex.h
index dbec4c301..77ac6e32a 100644
--- a/src/util/lp/lp_dual_simplex.h
+++ b/src/util/lp/lp_dual_simplex.h
@@ -6,8 +6,7 @@
 */
 #pragma once
 #include <vector>
-#include "util/sstream.h"
-#include "util/exception.h"
+#include "util/lp/lp_utils.h"
 #include "util/lp/lp_solver.h"
 #include "util/lp/lp_dual_core_solver.h"
 namespace lean {
diff --git a/src/util/lp/lp_primal_core_solver.cpp b/src/util/lp/lp_primal_core_solver.cpp
index 6702c51bf..42fdd2171 100644
--- a/src/util/lp/lp_primal_core_solver.cpp
+++ b/src/util/lp/lp_primal_core_solver.cpp
@@ -30,7 +30,7 @@ void lp_primal_core_solver<T, X>::sort_non_basis() {
     // reinit m_basis_heading
     for (unsigned j = 0; j < this->m_non_basic_columns.size(); j++) {
         unsigned col = this->m_non_basic_columns[j];
-        this->m_basis_heading[col] = - j - 1;
+        this->m_basis_heading[col] = - static_cast<int>(j) - 1;
         m_non_basis_list.push_back(col);
     }
 }
@@ -38,7 +38,7 @@ void lp_primal_core_solver<T, X>::sort_non_basis() {
 
 template <typename T, typename X>
 int lp_primal_core_solver<T, X>::choose_entering_column(unsigned number_of_benefitial_columns_to_go_over) { // at this moment m_y = cB * B(-1)
-    lean_assert(number_of_benefitial_columns_to_go_over);
+    lean_assert(number_of_benefitial_columns_to_go_over > 0);
     if (this->m_sort_counter == 0) {
         sort_non_basis();
         this->m_sort_counter = 20;
@@ -162,7 +162,8 @@ template <typename T, typename X>    int lp_primal_core_solver<T, X>::find_leavi
         }
         if (++i == this->m_m) i = 0;
     } while ( i != initial_i);
-    restore_harris_eps();
+    if (!precise<T>())
+        restore_harris_eps();
     return leaving;
 }
 
@@ -189,7 +190,7 @@ template <typename T, typename X>    bool lp_primal_core_solver<T, X>::try_jump_
 }
 
 template <typename T, typename X>    int lp_primal_core_solver<T, X>::find_leaving_and_t(unsigned entering, X & t){
-    if (numeric_traits<T>::precise()) return find_leaving_and_t_precisely(entering, t);
+    // if (numeric_traits<T>::precise()) return find_leaving_and_t_precisely(entering, t);
     X theta = get_harris_theta();
     lean_assert(theta >= zero_of_type<X>());
     if (try_jump_to_another_bound_on_entering(entering, theta, t)) return entering;
@@ -251,7 +252,7 @@ template <typename T, typename X> lp_primal_core_solver<T, X>::lp_primal_core_so
                                                                                      std::vector<X> & upper_bound_values,
                                                                                      lp_settings & settings,
                                                                                      std::unordered_map<unsigned, std::string> const & column_names):
-    lp_core_solver_base<T, X>(A, b,
+lp_core_solver_base<T, X>(A, b,
                               basis,
                               x,
                               costs,
@@ -260,8 +261,12 @@ template <typename T, typename X> lp_primal_core_solver<T, X>::lp_primal_core_so
                               column_type_array,
                               low_bound_values,
                               upper_bound_values),
-    m_beta(A.row_count()),
-    m_converted_harris_eps(convert_struct<T, double>::convert(this->m_settings.harris_feasibility_tolerance)) {
+    m_beta(A.row_count()) {
+    if (!(numeric_traits<T>::precise())) {
+        m_converted_harris_eps = convert_struct<T, double>::convert(this->m_settings.harris_feasibility_tolerance);
+    } else {
+        m_converted_harris_eps = zero_of_type<T>();
+    }
     this->m_status = UNKNOWN;
     this->m_column_norm_update_counter = settings.column_norms_update_frequency;
 }
@@ -504,7 +509,7 @@ template <typename T, typename X>void lp_primal_core_solver<T, X>::advance_on_en
 //     this->init_reduced_costs_for_one_iteration();
 // }
 
-template <typename T, typename X>    void lp_primal_core_solver<T, X>::advance_on_entering(int entering) {
+template <typename T, typename X> void lp_primal_core_solver<T, X>::advance_on_entering(int entering) {
     lean_assert(entering > -1);
     this->solve_Bd(entering);
     int refresh_result = refresh_reduced_cost_at_entering_and_check_that_it_is_off(entering);
@@ -545,7 +550,7 @@ template <typename T, typename X>    void lp_primal_core_solver<T, X>::push_forw
 }
 
 template <typename T, typename X>  unsigned lp_primal_core_solver<T, X>::get_number_of_non_basic_column_to_try_for_enter() {
-    unsigned ret = this->m_non_basic_columns.size();
+    unsigned ret = static_cast<unsigned>(this->m_non_basic_columns.size());
     if (this->m_status == TENTATIVE_UNBOUNDED)
         return ret; // we really need to find entering with a large reduced cost
     if (ret > 300) {
@@ -663,7 +668,11 @@ template <typename T, typename X>    void lp_primal_core_solver<T, X>::update_or
 // following Swietanowski - A new steepest ...
 template <typename T, typename X>    void lp_primal_core_solver<T, X>::update_column_norms(unsigned entering, unsigned leaving) {
     T pivot = this->m_pivot_row[entering];
-    T g_ent = std::max(calculate_norm_of_entering_exactly() / pivot /pivot, T(0.000001));
+    T g_ent = calculate_norm_of_entering_exactly() / pivot / pivot;
+    if (!numeric_traits<T>::precise()) {
+        if (g_ent < T(0.000001))
+            g_ent = T(0.000001);
+    }
     this->m_column_norms[leaving] = g_ent;
 
     for (unsigned j : this->m_pivot_row_index) {
@@ -733,7 +742,7 @@ template <typename T, typename X>    void lp_primal_core_solver<T, X>::add_colum
     }
 }
 
-template <typename T, typename X>    void lp_primal_core_solver<T, X>::one_iteration() {
+template <typename T, typename X> void lp_primal_core_solver<T, X>::one_iteration() {
     this->m_total_iterations++;
     unsigned number_of_benefitial_columns_to_go_over = get_number_of_non_basic_column_to_try_for_enter();
     int entering = choose_entering_column(number_of_benefitial_columns_to_go_over);
@@ -743,7 +752,7 @@ template <typename T, typename X>    void lp_primal_core_solver<T, X>::one_itera
         advance_on_entering(entering);
 }
 
-template <typename T, typename X>    void lp_primal_core_solver<T, X>::fill_breakpoints_array(unsigned entering) {
+template <typename T, typename X>  void lp_primal_core_solver<T, X>::fill_breakpoints_array(unsigned entering) {
     clear_breakpoints();
     for (unsigned i = this->m_m; i--;)
         try_add_breakpoint_in_row(i);
@@ -846,7 +855,8 @@ template <typename T, typename X>    void lp_primal_core_solver<T, X>::update_ba
     else
         this->update_x(entering, delta);
 }
-template <typename T, typename X>    bool lp_primal_core_solver<T, X>::column_is_feasible(unsigned j) const {
+
+template <typename T, typename X> bool lp_primal_core_solver<T, X>::column_is_feasible(unsigned j) const {
     const X& x = this->m_x[j];
     switch (this->m_column_type[j]) {
     case fixed:
@@ -879,6 +889,7 @@ template <typename T, typename X>    bool lp_primal_core_solver<T, X>::column_is
     default:
         lean_unreachable();
     }
+    return false; // it is unreachable
 }
 
 template <typename T, typename X> bool lp_primal_core_solver<T, X>::calc_current_x_is_feasible() const {
diff --git a/src/util/lp/lp_primal_core_solver.h b/src/util/lp/lp_primal_core_solver.h
index 27e26b174..68c85e4e1 100644
--- a/src/util/lp/lp_primal_core_solver.h
+++ b/src/util/lp/lp_primal_core_solver.h
@@ -16,7 +16,6 @@
 #include <math.h>
 #include <cstdlib>
 #include <algorithm>
-#include "util/numerics/double.h"
 #include "util/lp/lu.h"
 #include "util/lp/lp_solver.h"
 #include "util/lp/static_matrix.h"
@@ -43,7 +42,7 @@ public:
     bool m_recalc_reduced_costs = false;
     std::set<unsigned> m_forbidden_enterings;
     std::vector<T> m_beta; // see Swietanowski working vector beta for column norms
-    T m_epsilon_of_reduced_cost = T(1e-7);
+    T m_epsilon_of_reduced_cost = T(1)/T(10000000);
     bool m_exit_on_feasible_solution = false;
     std::vector<T> m_costs_backup;
     bool m_current_x_is_feasible;
@@ -56,7 +55,7 @@ public:
     int find_leaving_and_t_precisely(unsigned entering, X & t);
 
     static X positive_infinity() {
-        return convert_struct<X, double>::convert(std::numeric_limits<unsigned>::max());
+        return convert_struct<X, unsigned>::convert(std::numeric_limits<unsigned>::max());
     }
 
     X get_harris_theta();
@@ -91,7 +90,7 @@ public:
         limit_inf_on_bound_m_pos(m, this->m_x[j], this->m_upper_bound_values[j], theta);
     };
 
-    X harris_eps_for_bound(const X & bound) const { return ( convert_struct<X, double>::convert(1) +  abs(bound)/10) * m_converted_harris_eps/3;
+    X harris_eps_for_bound(const X & bound) const { return ( convert_struct<X, int>::convert(1) +  abs(bound)/10) * m_converted_harris_eps/3;
     }
 
 
diff --git a/src/util/lp/lp_primal_simplex.cpp b/src/util/lp/lp_primal_simplex.cpp
index ab907c240..cdd117ed5 100644
--- a/src/util/lp/lp_primal_simplex.cpp
+++ b/src/util/lp/lp_primal_simplex.cpp
@@ -31,7 +31,7 @@ template <typename T, typename X> void lp_primal_simplex<T, X>::init_buffer(unsi
 template <typename T, typename X> void lp_primal_simplex<T, X>::refactor() {
     m_core_solver->init_lu();
     if (m_core_solver->factorization()->get_status() != LU_status::OK) {
-        throw exception("cannot refactor");
+        throw_exception("cannot refactor");
     }
 }
 
@@ -226,6 +226,10 @@ template <typename T, typename X> void lp_primal_simplex<T, X>::fill_A_x_and_bas
 
 template <typename T, typename X> void lp_primal_simplex<T, X>::solve_with_total_inf() {
     int total_vars = this->m_A->column_count() + this->row_count();
+    if (total_vars == 0) {
+        this->m_status = OPTIMAL;
+        return;
+    }
     m_low_bounds.clear();
     m_low_bounds.resize(total_vars, zero_of_type<X>());  // low bounds are shifted ot zero
     this->m_x.resize(total_vars, numeric_traits<T>::zero());
@@ -237,10 +241,9 @@ template <typename T, typename X> void lp_primal_simplex<T, X>::solve_with_total
     this->fill_column_names_for_core_solver();
     unsigned j = this->m_A->column_count() - 1;
     unsigned core_solver_cols = this->number_of_core_structurals();
-
-    while (j >= core_solver_cols)
+    while (j >= core_solver_cols) {
         this->m_costs[j--] = numeric_traits<T>::zero();
-
+    }
     set_scaled_costs();
     m_core_solver = new lp_primal_core_solver<T, X>(*this->m_A,
                                                     this->m_b,
@@ -267,7 +270,9 @@ template <typename T, typename X> bool lp_primal_simplex<T, X>::bounds_hold(std:
     for (auto it : this->m_columns) {
         auto sol_it = solution.find(it.second->get_name());
         if (sol_it == solution.end()) {
-            throw exception(sstream() << "cannot find column " << it.first << " in solutio");
+            std::stringstream s;
+            s << "cannot find column " << it.first << " in solution";
+            throw_exception(s.str() );
         }
 
         if (!it.second->bounds_hold(sol_it->second)) {
@@ -282,7 +287,9 @@ template <typename T, typename X> bool lp_primal_simplex<T, X>::bounds_hold(std:
 template <typename T, typename X> T lp_primal_simplex<T, X>::get_row_value(unsigned i, std::unordered_map<std::string, T> const & solution, bool print) {
     auto it = this->m_A_values.find(i);
     if (it == this->m_A_values.end()) {
-        throw exception(sstream() << "cannot find row " << i);
+        std::stringstream s;
+        s << "cannot find row " << i;
+        throw_exception(s.str() );
     }
     T ret = numeric_traits<T>::zero();
     for (auto & pair : it->second) {
@@ -340,6 +347,7 @@ template <typename T, typename X> bool lp_primal_simplex<T, X>::row_constraint_h
         return true;;
     }
     lean_unreachable();
+    return false; // it is unreachable
 }
 
 template <typename T, typename X> bool lp_primal_simplex<T, X>::row_constraints_hold(std::unordered_map<std::string, T> const & solution) {
diff --git a/src/util/lp/lp_primal_simplex.h b/src/util/lp/lp_primal_simplex.h
index 5706c3012..11fa0d28d 100644
--- a/src/util/lp/lp_primal_simplex.h
+++ b/src/util/lp/lp_primal_simplex.h
@@ -9,8 +9,7 @@
 #include <unordered_map>
 #include <string>
 #include <algorithm>
-#include "util/exception.h"
-#include "util/sstream.h"
+#include "util/lp/lp_utils.h"
 #include "util/lp/column_info.h"
 #include "util/lp/lp_primal_core_solver.h"
 #include "util/lp/lp_solver.h"
diff --git a/src/util/lp/lp_settings.cpp b/src/util/lp/lp_settings.cpp
index 4e0210d25..402e96731 100644
--- a/src/util/lp/lp_settings.cpp
+++ b/src/util/lp/lp_settings.cpp
@@ -7,7 +7,6 @@
 #include <cmath>
 #include <string>
 #include <vector>
-#include "util/numerics/double.h"
 #include "util/lp/lp_settings.h"
 namespace lean {
 std::string column_type_to_string(column_type t) {
@@ -25,6 +24,7 @@ std::string column_type_to_string(column_type t) {
     default:
         lean_unreachable();
     }
+    return "unknown"; // it is unreachable
 }
 
 const char* lp_status_to_string(lp_status status) {
@@ -44,6 +44,7 @@ const char* lp_status_to_string(lp_status status) {
     default:
         lean_unreachable();
     }
+    return "UNKNOWN";  // it is unreachable
 }
 
 lp_status lp_status_from_string(std::string status) {
@@ -57,6 +58,7 @@ lp_status lp_status_from_string(std::string status) {
     if (status == "ITERATIONS_EXHAUSTED") return lp_status::ITERATIONS_EXHAUSTED;
     if (status == "EMPTY") return lp_status::EMPTY;
     lean_unreachable();
+    return lp_status::UNKNOWN; // it is unreachable
 }
 int get_millisecond_count() {
     timeb tb;
@@ -71,14 +73,14 @@ int get_millisecond_span(int start_time) {
     return span;
 }
 void my_random_init(unsigned * seed) {
-#ifdef LEAN_WINDOWS
+#if defined(LEAN_WINDOWS) || defined(WIN32)
     srand(*seed);
 #else
     rand_r(seed);
 #endif
 }
 unsigned my_random() {
-#ifdef LEAN_WINDOWS
+#if defined(LEAN_WINDOWS) || defined(WIN32)
     return rand(); // NOLINT
 #else
     return lrand48();
@@ -129,7 +131,7 @@ bool vectors_are_equal(const std::vector<T> & a, const buffer<T>  &b) {
 
 template <typename T>
 bool vectors_are_equal(const std::vector<T> & a, const std::vector<T>  &b) {
-    unsigned n = a.size();
+    unsigned n = static_cast<unsigned>(a.size());
     if (n != b.size()) return false;
     if (numeric_traits<T>::precise()) {
         for (unsigned i = 0; i < n; i ++){
diff --git a/src/util/lp/lp_settings.h b/src/util/lp/lp_settings.h
index dc27a8f33..7ec3c82d0 100644
--- a/src/util/lp/lp_settings.h
+++ b/src/util/lp/lp_settings.h
@@ -11,10 +11,7 @@
 #include <algorithm>
 #include <limits>
 #include <sys/timeb.h>
-#include "util/debug.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/sstream.h"
-#include "util/numerics/mpq.h"
+#include "util/lp/lp_utils.h"
 
 namespace lean {
 
@@ -51,11 +48,6 @@ lp_status lp_status_from_string(std::string status);
 enum non_basic_column_value_position { at_low_bound, at_upper_bound, at_fixed, free_of_bounds };
 
 template <typename X> bool is_epsilon_small(const X & v, const double& eps);    // forward definition
-template <typename X> bool precise() { return numeric_traits<X>::precise();}
-template <typename X> X zero_of_type() { return numeric_traits<X>::zero(); }
-template <typename X> X one_of_type() { return numeric_traits<X>::one(); }
-template <typename X> bool is_zero(const X & v) { return numeric_traits<X>::is_zero(v); }
-template <typename X> double  get_double(const X & v) { return numeric_traits<X>::get_double(v); }
 
 struct lp_settings {
     // when the absolute value of an element is less than pivot_epsilon
@@ -66,7 +58,7 @@ struct lp_settings {
     // a quatation "if some choice of the entering vairable leads to an eta matrix
     // whose diagonal element in the eta column is less than e2 (entering_diag_epsilon) in magnitude, the this choice is rejected ...
     double entering_diag_epsilon = 1e-8;
-    double c_partial_pivoting = 10; // this is the constant c from page 410
+    int c_partial_pivoting = 10; // this is the constant c from page 410
     unsigned depth_of_rook_search = 4;
     bool using_partial_pivoting = true;
     // dissertation of Achim Koberstein
@@ -151,7 +143,8 @@ struct lp_settings {
     // the method of lar solver to use
     bool row_feasibility = true;
     bool use_double_solver_for_lar = true;
-    int report_frequency = 10000000;
+    int report_frequency = 1000;
+    bool print_statistics = false;
     unsigned column_norms_update_frequency = 1000;
     bool scale_with_ratio = true;
     double density_threshold = 0.7; // need to tune it up, todo
@@ -210,4 +203,3 @@ inline void print_blanks(int n, std::ostream & out) {
     while (n--) {out << ' '; }
 }
 }
-
diff --git a/src/util/lp/lp_solver.cpp b/src/util/lp/lp_solver.cpp
index a7783bf22..77e19ec75 100644
--- a/src/util/lp/lp_solver.cpp
+++ b/src/util/lp/lp_solver.cpp
@@ -43,7 +43,9 @@ template <typename T, typename X> void lp_solver<T, X>::give_symbolic_name_to_co
 template <typename T, typename X>  T lp_solver<T, X>::get_column_value_by_name(std::string name) const {
     auto it = m_names_to_columns.find(name);
     if (it == m_names_to_columns.end()) {
-        throw exception(sstream() << "get_column_value_by_name " << name);
+        std::stringstream s;
+        s << "get_column_value_by_name " << name;
+        throw_exception(s.str());
     }
     return get_column_value(it -> second);
 }
@@ -249,7 +251,7 @@ template <typename T, typename X>    bool lp_solver<T, X>::row_e_is_obsolete(std
     return false;
 }
 
-template <typename T, typename X>    int lp_solver<T, X>::row_ge_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
+template <typename T, typename X>  bool lp_solver<T, X>::row_ge_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
     T rs = m_constraints[row_index].m_rs;
     if (row_is_zero(row)) {
         if (rs > zero_of_type<X>())
@@ -274,7 +276,7 @@ template <typename T, typename X>    int lp_solver<T, X>::row_ge_is_obsolete(std
     return false;
 }
 
-template <typename T, typename X>    bool lp_solver<T, X>::row_le_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
+template <typename T, typename X> bool lp_solver<T, X>::row_le_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index) {
     T low_bound;
     T rs = m_constraints[row_index].m_rs;
     if (row_is_zero(row)) {
@@ -318,6 +320,7 @@ template <typename T, typename X>    bool lp_solver<T, X>::row_is_obsolete(std::
         return row_le_is_obsolete(row, row_index);
     }
     lean_unreachable();
+    return false; // it is unreachable
 }
 
 template <typename T, typename X> void lp_solver<T, X>::remove_fixed_or_zero_columns() {
@@ -349,7 +352,7 @@ template <typename T, typename X> void lp_solver<T, X>::remove_fixed_or_zero_col
     }
 }
 
-template <typename T, typename X>    unsigned lp_solver<T, X>::try_to_remove_some_rows() {
+template <typename T, typename X> unsigned lp_solver<T, X>::try_to_remove_some_rows() {
     std::vector<unsigned> rows_to_delete;
     for (auto & t : m_A_values) {
         if (row_is_obsolete(t.second, t.first)) {
@@ -366,13 +369,13 @@ template <typename T, typename X>    unsigned lp_solver<T, X>::try_to_remove_som
         }
     }
     remove_fixed_or_zero_columns();
-    return rows_to_delete.size();
+    return static_cast<unsigned>(rows_to_delete.size());
 }
 
 template <typename T, typename X> void lp_solver<T, X>::cleanup() {
     int n = 0; // number of deleted rows
     int d;
-    while ((d = try_to_remove_some_rows()))
+    while ((d = try_to_remove_some_rows() > 0))
         n += d;
 
     // if (n == 1)
@@ -395,7 +398,7 @@ template <typename T, typename X> void lp_solver<T, X>::map_external_columns_to_
     for (auto & row : m_A_values) {
         for (auto & col : row.second) {
             if (col.second == numeric_traits<T>::zero() || m_columns[col.first]->is_fixed()) {
-                throw exception("found fixed column");
+                throw_exception("found fixed column");
             }
             unsigned j = col.first;
             auto j_place = m_external_columns_to_core_solver_columns.find(j);
@@ -426,7 +429,7 @@ template <typename T, typename X> void lp_solver<T, X>::unscale() {
 }
 
 template <typename T, typename X> void lp_solver<T, X>::fill_A_from_A_values() {
-    m_A = new static_matrix<T, X>(m_A_values.size(), number_of_core_structurals());
+    m_A = new static_matrix<T, X>(static_cast<unsigned>(m_A_values.size()), number_of_core_structurals());
     for (auto & t : m_A_values) {
         lean_assert(m_external_rows_to_core_solver_rows.find(t.first) != m_external_rows_to_core_solver_rows.end());
         unsigned row =  m_external_rows_to_core_solver_rows[t.first];
@@ -461,22 +464,19 @@ template <typename T, typename X> void lp_solver<T, X>::count_slacks_and_artific
 template <typename T, typename X> void lp_solver<T, X>::count_slacks_and_artificials_for_row(unsigned i) {
     lean_assert(this->m_constraints.find(this->m_core_solver_rows_to_external_rows[i]) != this->m_constraints.end());
     auto & constraint = this->m_constraints[this->m_core_solver_rows_to_external_rows[i]];
-    T rs;
     switch (constraint.m_relation) {
     case Equal:
         m_artificials++;
         break;
     case Greater_or_equal:
         m_slacks++;
-        rs = this->m_b[i];
-        if (rs > 0) {
+        if (this->m_b[i] > 0) {
             m_artificials++;
         }
         break;
     case Less_or_equal:
         m_slacks++;
-        rs = this->m_b[i];
-        if (rs < 0) {
+        if (this->m_b[i] < 0) {
             m_artificials++;
         }
         break;
@@ -488,7 +488,7 @@ template <typename T, typename X>    T lp_solver<T, X>::low_bound_shift_for_row(
 
     auto row = this->m_A_values.find(i);
     if (row == this->m_A_values.end()) {
-        throw exception("cannot find row");
+        throw_exception("cannot find row");
     }
     for (auto col : row->second) {
         ret += col.second * this->m_columns[col.first]->get_shift();
diff --git a/src/util/lp/lp_solver.h b/src/util/lp/lp_solver.h
index 7e52f69a0..63a1e8020 100644
--- a/src/util/lp/lp_solver.h
+++ b/src/util/lp/lp_solver.h
@@ -10,8 +10,6 @@
 #include <unordered_map>
 #include <algorithm>
 #include <vector>
-#include "util/exception.h"
-#include "util/sstream.h"
 #include "util/lp/lp_settings.h"
 #include "util/lp/column_info.h"
 #include "util/lp/static_matrix.h"
@@ -173,7 +171,7 @@ protected:
 
     bool row_e_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
 
-    int row_ge_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
+    bool row_ge_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
 
     bool row_le_is_obsolete(std::unordered_map<unsigned, T> & row, unsigned row_index);
 
@@ -199,7 +197,7 @@ protected:
 
     void fill_column_names_for_core_solver();
 
-    unsigned number_of_core_structurals() { return m_external_columns_to_core_solver_columns.size(); }
+    unsigned number_of_core_structurals() { return static_cast<unsigned>(m_external_columns_to_core_solver_columns.size()); }
 
 
     void restore_column_scales_to_one() {
diff --git a/src/util/lp/lu.cpp b/src/util/lp/lu.cpp
index a6cc52765..d07d2dfbc 100644
--- a/src/util/lp/lu.cpp
+++ b/src/util/lp/lu.cpp
@@ -8,6 +8,7 @@
 #include <algorithm>
 #include <set>
 #include <vector>
+#include <utility>
 #include "util/lp/lu.h"
 namespace lean {
 #ifdef LEAN_DEBUG
@@ -33,7 +34,7 @@ template<typename T, typename X>
 void print_matrix(static_matrix<T, X> &m, std::ostream & out) {
     std::vector<std::vector<std::string>> A;
     std::vector<unsigned> widths;
-    std::set<pair<unsigned, unsigned>> domain = m.get_domain();
+    std::set<std::pair<unsigned, unsigned>> domain = m.get_domain();
     for (unsigned i = 0; i < m.row_count(); i++) {
         A.push_back(std::vector<std::string>());
         for (unsigned j = 0; j < m.column_count(); j++) {
@@ -169,9 +170,13 @@ void lu<T, X>::solve_Bd_when_w_is_ready(std::vector<T> & d, indexed_vector<T>& w
 template <typename T, typename X>
 template <typename L>
 void lu<T, X>::solve_By_when_y_is_ready(std::vector<L> & y) {
+    if (numeric_traits<T>::precise()) {
+        m_U.solve_U_y(y);
+        m_R.apply_reverse_from_left(y); // see 24.3 from Chvatal
+        return;
+    }
     m_U.double_solve_U_y(y);
     m_R.apply_reverse_from_left(y); // see 24.3 from Chvatal
-    if (precise<X>()) return;
     unsigned i = m_dim;
     while (i--) {
         if (is_zero(y[i])) continue;
@@ -232,7 +237,7 @@ void lu<T, X>::solve_yB_internal(std::vector<T>& y) {
 }
 template <typename T, typename X>
 void lu<T, X>::add_delta_to_solution(std::vector<T>& yc, std::vector<T>& y){
-    unsigned i = y.size();
+    unsigned i = static_cast<unsigned>(y.size());
     while (i--)
         y[i]+=yc[i];
 }
@@ -249,6 +254,10 @@ void lu<T, X>::find_error_of_yB(std::vector<T>& yc, const std::vector<T>& y) {
 // y is the input
 template <typename T, typename X>
 void lu<T, X>::solve_yB(std::vector<T> & y) {
+    if (numeric_traits<T>::precise()) {
+        solve_yB_internal(y);
+        return;
+    }
     std::vector<T> yc(y); // copy y aside
     solve_yB_internal(y);
     find_error_of_yB(yc, y);
@@ -427,7 +436,7 @@ void lu<T, X>::check_apply_lp_lists_to_w(T * w) {
 template <typename T, typename X>
 void lu<T, X>::process_column(int j) {
     unsigned pi, pj;
-    m_U.get_pivot_for_column(pi, pj, T(m_settings.c_partial_pivoting), j);
+    m_U.get_pivot_for_column(pi, pj, m_settings.c_partial_pivoting, j);
     if (static_cast<int>(pi) == -1) {
         std::cout << "cannot find the pivot for column " << j << std::endl;
         m_failure = true;
@@ -693,7 +702,7 @@ void lu<T, X>::replace_column(unsigned leaving, T pivot_elem, indexed_vector<T>
         push_matrix_to_tail(row_eta);
     }
     calculate_Lwave_Pwave_for_bump(replaced_column, lowest_row_of_the_bump);
-    lean_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
+     lean_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
     lean_assert(w.is_OK() && m_row_eta_work_vector.is_OK());
 }
 template <typename T, typename X>
diff --git a/src/util/lp/lu.h b/src/util/lp/lu.h
index 664509501..d8dfd5862 100644
--- a/src/util/lp/lu.h
+++ b/src/util/lp/lu.h
@@ -9,8 +9,6 @@
 
 #include <vector>
 #include "util/debug.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
 #include <algorithm>
 #include <set>
 #include "util/lp/sparse_matrix.h"
diff --git a/src/util/lp/lu_instances.cpp b/src/util/lp/lu_instances.cpp
index 2a4d07ac4..c137e8c52 100644
--- a/src/util/lp/lu_instances.cpp
+++ b/src/util/lp/lu_instances.cpp
@@ -33,7 +33,6 @@ template void lean::init_factorization<lean::mpq, lean::mpq>(lean::lu<lean::mpq,
 template void lean::init_factorization<lean::mpq, lean::numeric_pair<lean::mpq> >(lean::lu<lean::mpq, lean::numeric_pair<lean::mpq> >*&, lean::static_matrix<lean::mpq, lean::numeric_pair<lean::mpq> >&, std::vector<unsigned int, std::allocator<unsigned int> >&, std::vector<int, std::allocator<int> >&, lean::lp_settings&, std::vector<unsigned int, std::allocator<unsigned int> >&);
 #ifdef LEAN_DEBUG
 template void lean::print_matrix<double, double>(lean::sparse_matrix<double, double>&, std::ostream & out);
-template void lean::print_matrix<float, float>(lean::sparse_matrix<float, float>&, std::ostream & out);
 template void lean::print_matrix<double, double>(lean::static_matrix<double, double>&, std::ostream & out);
 template bool lean::lu<double, double>::is_correct();
 template lean::dense_matrix<double, double> lean::get_B<double, double>(lean::lu<double, double>&);
diff --git a/src/util/lp/matrix.h b/src/util/lp/matrix.h
index ed0dc7e72..c4dd06bab 100644
--- a/src/util/lp/matrix.h
+++ b/src/util/lp/matrix.h
@@ -7,7 +7,7 @@
 #ifdef LEAN_DEBUG
 #pragma once
 #include "util/numerics/numeric_traits.h"
-#include "util/numerics/double.h"
+// #include "util/numerics/double.h"
 #include <vector>
 #include <string>
 #include "util/lp/lp_settings.h"
diff --git a/src/util/lp/numeric_pair.h b/src/util/lp/numeric_pair.h
index 432ec664d..01d75d377 100644
--- a/src/util/lp/numeric_pair.h
+++ b/src/util/lp/numeric_pair.h
@@ -3,25 +3,67 @@
   Released under Apache 2.0 license as described in the file LICENSE.
 
   Author: Lev Nachmanson
+  The idea is that it is only one different file in Lean and z3 source inside of LP
 */
 #pragma once
+// #define lp_for_z3
 #include <string>
 #include <cmath>
 #include <algorithm>
-#include "util/numerics/mpq.h"
-#include "util/numerics/double.h"
-#include "util/lp/lp_settings.h"
+#ifdef lp_for_z3
+#include "../rational.h"
+#include "../sstream.h"
+#include "../z3_exception.h"
 
+#else
+#include "util/numerics/mpq.h"
+#include "util/numerics/numeric_traits.h"
+#endif
 namespace lean {
+#ifdef lp_for_z3 // rename rationals
+    typedef rational mpq;
+#else
+    typedef lean::mpq mpq;
+#endif
+
+
 template <typename T>
 std::string T_to_string(const T & t); // forward definition
+#ifdef lp_for_z3
+template <typename T> class numeric_traits {};
+
+template <>  class numeric_traits<double> {
+    public:
+        static bool precise() { return false; }
+        static double g_zero;
+        static double const &zero() { return g_zero;  }
+        static double g_one;
+        static double const &one() { return g_one; }
+        static bool is_zero(double v) { return v == 0.0; }
+        static double const & get_double(double const & d) { return d;}
+        static double log(double const & d) { NOT_IMPLEMENTED_YET(); return d;}
+        static double from_string(std::string const & str) { return atof(str.c_str()); }
+    };
+
+    template<>
+    class numeric_traits<rational> {
+    public:
+        static bool precise() { return true; }
+        static rational const & zero() { return rational::zero(); }
+        static rational const & one() { return rational::one(); }
+        static bool is_zero(const rational & v) { return v.is_zero(); }
+        static double const  get_double(const rational  & d) { return d.get_double();}
+        static rational log(rational const& r) { UNREACHABLE(); return r; }
+        static rational from_string(std::string const & str) { return rational(str.c_str()); }
+    };
+#endif
 
 template <typename X, typename Y>
 struct convert_struct {
     static X convert(const Y & y){ return X(y);}
-    static bool is_epsilon_small(const X & x,  const double & y) { return std::abs(get_double(x)) < y; }
-    static bool below_bound_numeric(const X &, const X &, const Y &) { lean_unreachable(); }
-    static bool above_bound_numeric(const X &, const X &, const Y &) { lean_unreachable(); }
+    static bool is_epsilon_small(const X & x,  const double & y) { return std::abs(numeric_traits<X>::get_double(x)) < y; }
+    static bool below_bound_numeric(const X &, const X &, const Y &) { /*lean_unreachable();*/ return false;}
+    static bool above_bound_numeric(const X &, const X &, const Y &) { /*lean_unreachable();*/ return false; }
 };
 
 
@@ -75,7 +117,7 @@ struct numeric_pair {
     }
 
     numeric_pair operator/(const numeric_pair &) const {
-        lean_unreachable();
+        //        lean_unreachable();
     }
 
 
@@ -84,7 +126,7 @@ struct numeric_pair {
     }
 
     numeric_pair operator*(const numeric_pair & /*a*/) const  {
-        lean_unreachable();
+        //        lean_unreachable();
     }
 
     numeric_pair&  operator+=(const numeric_pair & a) {
@@ -119,7 +161,7 @@ struct numeric_pair {
         return numeric_pair(-x, -y);
     }
 
-    static bool precize() { return numeric_traits<T>::precize();}
+    static bool precize() { return lean::numeric_traits<T>::precize();}
 
     std::string to_string() const { return std::string("(") + T_to_string(x) + ", "  + T_to_string(y) + ")"; }
 };
@@ -141,20 +183,23 @@ numeric_pair<T> operator*(const numeric_pair<T> & r, const X & a) {
     return numeric_pair<T>(a * r.x, a * r.y);
 }
 
+} // close namespace lean
 
-
-
+namespace lean {
 // template <numeric_pair, typename T>  bool precise() { return numeric_traits<T>::precise();}
-template <typename T> double get_double(const numeric_pair<T> & ) { lean_unreachable(); }
+template <typename T> double get_double(const lean::numeric_pair<T> & ) { /* lean_unreachable(); */ return 0;}
 template <typename T>
-class numeric_traits<numeric_pair<T>> {
+class numeric_traits<lean::numeric_pair<T>> {
   public:
     static bool precise() { return numeric_traits<T>::precise();}
-    static numeric_pair<T> zero() { return numeric_pair<T>(numeric_traits<T>::zero(), numeric_traits<T>::zero()); }
-    static bool is_zero(const numeric_pair<T> & v) { return numeric_traits<T>::is_zero(v.x) && numeric_traits<T>::is_zero(v.y); }
-    static double get_double(const numeric_pair<T> & v){ return numeric_traits<T>::get_double(v.x); } // just return the double of the first coordinate
-    static double one() { lean_unreachable(); }
+    static lean::numeric_pair<T> zero() { return lean::numeric_pair<T>(numeric_traits<T>::zero(), numeric_traits<T>::zero()); }
+    static bool is_zero(const lean::numeric_pair<T> & v) { return numeric_traits<T>::is_zero(v.x) && numeric_traits<T>::is_zero(v.y); }
+    static double get_double(const lean::numeric_pair<T> & v){ return numeric_traits<T>::get_double(v.x); } // just return the double of the first coordinate
+    static double one() { /*lean_unreachable();*/ return 0;}
 };
+} // close namespace lean
+
+namespace lean {
 template <>
 struct convert_struct<double, numeric_pair<double>> {
     static double convert(const numeric_pair<double> & q) {return q.x;}
@@ -165,19 +210,20 @@ template <typename X> bool is_epsilon_small(const X & v, const double& eps);   /
 template <typename T>
 struct convert_struct<numeric_pair<T>, double> {
     static numeric_pair<T> convert(const double & q) {
-        return numeric_pair<T>(convert_struct<T, double>::convert(q), zero_of_type<T>());
+        return numeric_pair<T>(convert_struct<T, double>::convert(q), numeric_traits<T>::zero());
     }
     static bool is_epsilon_small(const numeric_pair<T> & p, const double & eps) {
         return convert_struct<T, double>::is_epsilon_small(p.x, eps) && convert_struct<T, double>::is_epsilon_small(p.y, eps);
     }
     static bool below_bound_numeric(const numeric_pair<T> &, const numeric_pair<T> &, const double &) {
-        lean_unreachable();
+        // lean_unreachable();
+        return false;
     }
     static bool above_bound_numeric(const numeric_pair<T> &, const numeric_pair<T> &, const double &) {
-        lean_unreachable();
+        // lean_unreachable();
+        return false;
     }
 };
-
 template <>
 struct convert_struct<numeric_pair<double>, double> {
     static numeric_pair<double> convert(const double & q) {
@@ -188,7 +234,6 @@ struct convert_struct<numeric_pair<double>, double> {
     }
 
     static int compare_on_coord(const double & x, const double & bound, const double eps) {
-        lean_assert(eps > 0);
         if (bound == 0) return (x < - eps)? -1: (x > eps? 1 : 0); // it is an important special case
         double relative = (bound > 0)? - eps: eps;
         return (x < bound * (1.0 + relative) - eps)? -1 : ((x > bound * (1.0 - relative) + eps)? 1 : 0);
@@ -199,7 +244,6 @@ struct convert_struct<numeric_pair<double>, double> {
         if (r == 1) return false;
         if (r == -1) return true;
         // the first coordinates are almost the same
-        lean_assert(r == 0);
         return compare_on_coord(x.y, bound.y, eps) == -1;
     }
 
@@ -208,7 +252,6 @@ struct convert_struct<numeric_pair<double>, double> {
         if (r == -1) return false;
         if (r ==  1) return true;
         // the first coordinates are almost the same
-        lean_assert(r == 0);
         return compare_on_coord(x.y, bound.y, eps) == 1;
     }
 };
@@ -220,13 +263,11 @@ struct convert_struct<double, double> {
     }
     static double convert(const double & y){ return y;}
     static bool below_bound_numeric(const double & x, const double & bound, const double & eps) {
-        lean_assert(eps > 0);
         if (bound == 0) return x < - eps;
         double relative = (bound > 0)? - eps: eps;
         return x < bound * (1.0 + relative) - eps;
     }
     static bool above_bound_numeric(const double & x, const double & bound, const double & eps) {
-        lean_assert(eps > 0);
         if (bound == 0) return x > eps;
         double relative = (bound > 0)?  eps: - eps;
         return x > bound * (1.0 + relative) + eps;
diff --git a/src/util/lp/permutation_matrix.cpp b/src/util/lp/permutation_matrix.cpp
index 8e0bd2f21..95c9fee65 100644
--- a/src/util/lp/permutation_matrix.cpp
+++ b/src/util/lp/permutation_matrix.cpp
@@ -77,7 +77,7 @@ void permutation_matrix<T, X>::apply_from_left_perm(indexed_vector<L> & w, lp_se
     copy_aside(t, tmp_index, w);
     clear_data(w);
     // set the new values
-    for (unsigned i = t.size(); i > 0;) {
+    for (unsigned i = static_cast<unsigned>(t.size()); i > 0;) {
         i--;
         unsigned j = m_rev[tmp_index[i]];
         w[j] = t[i];
@@ -112,7 +112,7 @@ template <typename T, typename X> void permutation_matrix<T, X>::apply_from_righ
 
 template <typename T, typename X> template <typename L>
 void permutation_matrix<T, X>::copy_aside(std::vector<L> & t, std::vector<unsigned> & tmp_index, indexed_vector<L> & w) {
-    for (unsigned i = t.size(); i > 0;) {
+    for (unsigned i = static_cast<unsigned>(t.size()); i > 0;) {
         i--;
         unsigned j = w.m_index[i];
         t[i] = w[j]; // copy aside all non-zeroes
@@ -123,7 +123,7 @@ void permutation_matrix<T, X>::copy_aside(std::vector<L> & t, std::vector<unsign
 template <typename T, typename X> template <typename L>
 void permutation_matrix<T, X>::clear_data(indexed_vector<L> & w) {
     // clear old non-zeroes
-    for (unsigned i = w.m_index.size(); i > 0;) {
+    for (unsigned i = static_cast<unsigned>(w.m_index.size()); i > 0;) {
         i--;
         unsigned j = w.m_index[i];
         w[j] = zero_of_type<L>();
@@ -145,7 +145,7 @@ void permutation_matrix<T, X>::apply_reverse_from_left(indexed_vector<L> & w) {
     clear_data(w);
 
     // set the new values
-    for (unsigned i = t.size(); i > 0;) {
+    for (unsigned i = static_cast<unsigned>(t.size()); i > 0;) {
         i--;
         unsigned j = m_permutation[tmp_index[i]];
         w[j] = t[i];
diff --git a/src/util/lp/permutation_matrix.h b/src/util/lp/permutation_matrix.h
index 8e10e1d0f..27d971415 100644
--- a/src/util/lp/permutation_matrix.h
+++ b/src/util/lp/permutation_matrix.h
@@ -8,14 +8,6 @@
 #include <vector>
 #include <algorithm>
 #include "util/debug.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
-#include "util/numerics/mpq.h"
-#include "util/numerics/mpz.h"
-#include "util/numerics/mpbq.h"
-#include "util/numerics/double.h"
-#include "util/numerics/float.h"
-#include "util/numerics/mpfp.h"
 #include <string>
 #include "util/lp/sparse_vector.h"
 #include "util/lp/indexed_vector.h"
@@ -129,7 +121,7 @@ namespace lean {
 
         bool is_identity() const;
 
-        unsigned size() const { return m_rev.size(); }
+        unsigned size() const { return static_cast<unsigned>(m_rev.size()); }
 
         unsigned * values() const { return m_permutation; }
     }; // end of the permutation class
diff --git a/src/util/lp/permutation_matrix_instances.cpp b/src/util/lp/permutation_matrix_instances.cpp
index f8a0d9bac..952cf1da4 100644
--- a/src/util/lp/permutation_matrix_instances.cpp
+++ b/src/util/lp/permutation_matrix_instances.cpp
@@ -16,10 +16,7 @@ template void lean::permutation_matrix<double, double>::multiply_by_permutation_
 template void lean::permutation_matrix<double, double>::multiply_by_reverse_from_right(lean::permutation_matrix<double, double>&);
 template lean::permutation_matrix<double, double>::permutation_matrix(unsigned int, std::vector<unsigned int, std::allocator<unsigned int> > const&);
 template void lean::permutation_matrix<double, double>::transpose_from_left(unsigned int, unsigned int);
-template void lean::permutation_matrix<float, float>::apply_from_right(std::vector<float, std::allocator<float> >&);
-template lean::permutation_matrix<float, float>::permutation_matrix(unsigned int);
-template void lean::permutation_matrix<float, float>::transpose_from_left(unsigned int, unsigned int);
-template void lean::permutation_matrix<float, float>::transpose_from_right(unsigned int, unsigned int);
+
 template void lean::permutation_matrix<lean::mpq, lean::mpq>::apply_from_right(std::vector<lean::mpq, std::allocator<lean::mpq> >&);
 template bool lean::permutation_matrix<lean::mpq, lean::mpq>::is_identity() const;
 template void lean::permutation_matrix<lean::mpq, lean::mpq>::multiply_by_permutation_from_left(lean::permutation_matrix<lean::mpq, lean::mpq>&);
@@ -44,8 +41,6 @@ template void lean::permutation_matrix<double, double>::apply_reverse_from_left<
 template void lean::permutation_matrix<double, double>::apply_reverse_from_left<double>(std::vector<double, std::allocator<double> >&);
 template void lean::permutation_matrix<double, double>::apply_reverse_from_right<double>(std::vector<double, std::allocator<double> >&);
 template void lean::permutation_matrix<double, double>::transpose_from_right(unsigned int, unsigned int);
-template void lean::permutation_matrix<float, float>::apply_from_left_perm<float>(lean::indexed_vector<float>&, lean::lp_settings&);
-template void lean::permutation_matrix<float, float>::apply_from_left_perm<float>(std::vector<float, std::allocator<float> >&);
 template void lean::permutation_matrix<lean::mpq, lean::mpq>::apply_from_left_perm<lean::mpq>(lean::indexed_vector<lean::mpq>&, lean::lp_settings&);
 template void lean::permutation_matrix<lean::mpq, lean::mpq>::apply_from_left_perm<lean::mpq>(std::vector<lean::mpq, std::allocator<lean::mpq> >&);
 template void lean::permutation_matrix<lean::mpq, lean::mpq>::apply_reverse_from_left<lean::mpq>(lean::indexed_vector<lean::mpq>&);
diff --git a/src/util/lp/row_eta_matrix.cpp b/src/util/lp/row_eta_matrix.cpp
index f00884627..bb7bfc28d 100644
--- a/src/util/lp/row_eta_matrix.cpp
+++ b/src/util/lp/row_eta_matrix.cpp
@@ -131,7 +131,7 @@ void row_eta_matrix<T, X>::conjugate_by_permutation(permutation_matrix<T, X> & p
     std::vector<unsigned> columns;
     for (auto & it : m_row_vector.m_data)
         columns.push_back(it.first);
-    for (unsigned i = columns.size(); i-- > 0;)
+    for (unsigned i = static_cast<unsigned>(columns.size()); i-- > 0;)
         m_row_vector.m_data[i].first = p.get_rev(columns[i]);
 #ifdef LEAN_DEBUG
     // lean_assert(deb == *this);
diff --git a/src/util/lp/row_eta_matrix.h b/src/util/lp/row_eta_matrix.h
index ec6d819e7..96a8d5faf 100644
--- a/src/util/lp/row_eta_matrix.h
+++ b/src/util/lp/row_eta_matrix.h
@@ -8,14 +8,6 @@
 #pragma once
 #include <vector>
 #include "util/debug.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
-#include "util/numerics/mpq.h"
-#include "util/numerics/mpz.h"
-#include "util/numerics/mpbq.h"
-#include "util/numerics/double.h"
-#include "util/numerics/float.h"
-#include "util/numerics/mpfp.h"
 #include <string>
 #include "util/lp/sparse_vector.h"
 #include "util/lp/indexed_vector.h"
diff --git a/src/util/lp/scaler.h b/src/util/lp/scaler.h
index bfb7ccf3a..bc536cdd1 100644
--- a/src/util/lp/scaler.h
+++ b/src/util/lp/scaler.h
@@ -11,7 +11,7 @@
 #include <algorithm>
 #include <stdio.h>      /* printf, fopen */
 #include <stdlib.h>     /* exit, EXIT_FAILURE */
-#include "util/numerics/double.h"
+#include "util/lp/lp_utils.h"
 #include "util/lp/static_matrix.h"
 namespace lean {
 // for scaling an LP
diff --git a/src/util/lp/sparse_matrix.cpp b/src/util/lp/sparse_matrix.cpp
index da73e6842..59e4dab7b 100644
--- a/src/util/lp/sparse_matrix.cpp
+++ b/src/util/lp/sparse_matrix.cpp
@@ -11,13 +11,13 @@ namespace lean {
 template <typename T, typename X>
 void sparse_matrix<T, X>::copy_column_from_static_matrix(unsigned col, static_matrix<T, X> const &A, unsigned col_index_in_the_new_matrix) {
     std::vector<column_cell<T>> const & A_col_vector = A.m_columns[col];
-    unsigned size = A_col_vector.size();
+    unsigned size = static_cast<unsigned>(A_col_vector.size());
     std::vector<indexed_value<T>> & new_column_vector = m_columns[col_index_in_the_new_matrix].m_values;
     for (unsigned l = 0; l < size; l++) {
         column_cell<T> const & col_cell = A_col_vector[l];
-        unsigned col_offset = new_column_vector.size();
+        unsigned col_offset = static_cast<unsigned>(new_column_vector.size());
         std::vector<indexed_value<T>> & row_vector = m_rows[col_cell.m_i];
-        unsigned row_offset = row_vector.size();
+        unsigned row_offset = static_cast<unsigned>(row_vector.size());
         new_column_vector.push_back(indexed_value<T>(col_cell.m_value, col_cell.m_i, row_offset));
         row_vector.push_back(indexed_value<T>(col_cell.m_value, col_index_in_the_new_matrix, col_offset));
     }
@@ -193,7 +193,7 @@ void sparse_matrix<T, X>::set_max_in_row(std::vector<indexed_value<T>> & row_val
         return;
     T max_val = numeric_traits<T>::zero();
     int max_index = 0;
-    for (unsigned i = row_vals.size(); i-- > 0;) {
+    for (unsigned i = static_cast<unsigned>(row_vals.size()); i-- > 0;) {
         T iabs = abs(row_vals[i].m_value);
         if (iabs > max_val) {
             max_val = iabs;
@@ -324,7 +324,7 @@ bool sparse_matrix<T, X>::set_row_from_work_vector_and_clean_work_vector(unsigne
 template <typename T, typename X>
 void sparse_matrix<T, X>::remove_zero_elements_and_set_data_on_existing_elements(unsigned row) {
     auto & row_vals = m_rows[row];
-    for (unsigned k = row_vals.size(); k-- > 0;) { // we cannot simply run the iterator since we are removing
+    for (unsigned k = static_cast<unsigned>(row_vals.size()); k-- > 0;) { // we cannot simply run the iterator since we are removing
         // elements from row_vals
         auto & row_el_iv = row_vals[k];
         unsigned j = row_el_iv.m_index;
@@ -342,7 +342,7 @@ void sparse_matrix<T, X>::remove_zero_elements_and_set_data_on_existing_elements
 template <typename T, typename X>
 void sparse_matrix<T, X>::remove_zero_elements_and_set_data_on_existing_elements_not_adjusted(unsigned row, indexed_vector<T> & work_vec, lp_settings & settings) {
     auto & row_vals = m_rows[row];
-    for (unsigned k = row_vals.size(); k-- > 0;) { // we cannot simply run the iterator since we are removing
+    for (unsigned k = static_cast<unsigned>(row_vals.size()); k-- > 0;) { // we cannot simply run the iterator since we are removing
         // elements from row_vals
         auto & row_el_iv = row_vals[k];
         unsigned j = row_el_iv.m_index;
@@ -619,7 +619,7 @@ void sparse_matrix<T, X>::remove_elements_that_are_not_in_w_and_update_common_el
     // --------------------------------
     // column_vals represents the old column
     auto & column_vals = m_columns[column_to_replace].m_values;
-    for (unsigned k = column_vals.size(); k-- > 0;) {
+    for (unsigned k = static_cast<unsigned>(column_vals.size()); k-- > 0;) {
         indexed_value<T> & col_el_iv = column_vals[k];
         unsigned i = col_el_iv.m_index;
         T w_data_at_i = w[adjust_row_inverse(i)];
@@ -647,8 +647,8 @@ template <typename T, typename X>
 void sparse_matrix<T, X>::add_new_element(unsigned row, unsigned col, T val) {
     auto & row_vals = m_rows[row];
     auto & col_vals = m_columns[col].m_values;
-    unsigned row_el_offs = row_vals.size();
-    unsigned col_el_offs = col_vals.size();
+    unsigned row_el_offs = static_cast<unsigned>(row_vals.size());
+    unsigned col_el_offs = static_cast<unsigned>(col_vals.size());
     row_vals.push_back(indexed_value<T>(val, col, col_el_offs));
     col_vals.push_back(indexed_value<T>(val, row, row_el_offs));
 }
@@ -666,7 +666,7 @@ void sparse_matrix<T, X>::add_new_elements_of_w_and_clear_w(unsigned column_to_r
             auto & row_chunk = m_rows[ai];
             lean_assert(row_chunk.size() > 0);
             if (abs(w_at_i) > abs(row_chunk[0].m_value))
-                put_max_index_to_0(row_chunk, row_chunk.size() - 1);
+                put_max_index_to_0(row_chunk, static_cast<unsigned>(row_chunk.size()) - 1);
         }
         w[i] = numeric_traits<T>::zero();
     }
@@ -690,7 +690,7 @@ unsigned sparse_matrix<T, X>::pivot_score(unsigned i, unsigned j) {
 
     auto col_header = m_columns[j];
 
-    return get_row_values(i).size() * (col_header.m_values.size() - col_header.m_shortened_markovitz - 1);
+    return static_cast<unsigned>(get_row_values(i).size() * (col_header.m_values.size() - col_header.m_shortened_markovitz - 1));
 }
 
 template <typename T, typename X>
@@ -698,10 +698,10 @@ void sparse_matrix<T, X>::enqueue_domain_into_pivot_queue() {
     lean_assert(m_pivot_queue.size() == 0);
     for (unsigned i = 0; i < dimension(); i++) {
         auto & rh = m_rows[i];
-        unsigned rnz = rh.size();
+        unsigned rnz = static_cast<unsigned>(rh.size());
         for (auto iv : rh) {
             unsigned j = iv.m_index;
-            m_pivot_queue.enqueue(i, j, rnz * (m_columns[j].m_values.size() - 1));
+            m_pivot_queue.enqueue(i, j, rnz * (static_cast<unsigned>(m_columns[j].m_values.size()) - 1));
         }
     }
 }
@@ -735,7 +735,7 @@ void sparse_matrix<T, X>::recover_pivot_queue(std::vector<upair> & rejected_pivo
 }
 
 template <typename T, typename X>
-int sparse_matrix<T, X>::elem_is_too_small(unsigned i, unsigned j, const T & c_partial_pivoting)  {
+int sparse_matrix<T, X>::elem_is_too_small(unsigned i, unsigned j, int c_partial_pivoting)  {
     auto & row_chunk = m_rows[i];
 
     if (j == row_chunk[0].m_index) {
@@ -779,12 +779,12 @@ void sparse_matrix<T, X>::update_active_pivots(unsigned row) {
     unsigned arow = adjust_row(row);
     for (auto & iv : m_rows[arow]) {
         col_header & ch = m_columns[iv.m_index];
-        int cols = ch.m_values.size() - ch.m_shortened_markovitz - 1;
+        int cols = static_cast<int>(ch.m_values.size()) - ch.m_shortened_markovitz - 1;
         lean_assert(cols >= 0);
         for (auto &ivc : ch.m_values) {
             unsigned i = ivc.m_index;
             if (adjust_row_inverse(i) <= row) continue; // the i is not an active row
-            m_pivot_queue.enqueue(i, iv.m_index, m_rows[i].size()*cols);
+            m_pivot_queue.enqueue(i, iv.m_index, static_cast<unsigned>(m_rows[i].size())*cols);
         }
     }
 }
@@ -802,10 +802,10 @@ bool sparse_matrix<T, X>::shorten_active_matrix(unsigned row, eta_matrix<T, X> *
     for (auto & it : eta_matrix->m_column_vector.m_data) {
         unsigned row = adjust_row(it.first);
         auto & row_values = m_rows[row];
-        unsigned rnz = row_values.size();
+        unsigned rnz = static_cast<unsigned>(row_values.size());
         for (auto & iv : row_values) {
             col_header& ch = m_columns[iv.m_index];
-            int cnz = ch.m_values.size() - ch.m_shortened_markovitz - 1;
+            int cnz = static_cast<int>(ch.m_values.size()) - ch.m_shortened_markovitz - 1;
             lean_assert(cnz >= 0);
             m_pivot_queue.enqueue(row, iv.m_index, rnz * cnz);
         }
@@ -903,7 +903,7 @@ bool sparse_matrix<T, X>::pivot_queue_is_correct_after_pivoting(int k) {
 #endif
 
 template <typename T, typename X>
-bool sparse_matrix<T, X>::get_pivot_for_column(unsigned &i, unsigned &j, T const & c_partial_pivoting, unsigned k) {
+bool sparse_matrix<T, X>::get_pivot_for_column(unsigned &i, unsigned &j, int c_partial_pivoting, unsigned k) {
     std::vector<upair> pivots_candidates_that_are_too_small;
     while (!m_pivot_queue.is_empty()) {
         m_pivot_queue.dequeue(i, j);
@@ -933,7 +933,7 @@ bool sparse_matrix<T, X>::get_pivot_for_column(unsigned &i, unsigned &j, T const
 }
 
 template <typename T, typename X>
-bool sparse_matrix<T, X>::elem_is_too_small(std::vector<indexed_value<T>> & row_chunk, indexed_value<T> & iv, T const & c_partial_pivoting) {
+bool sparse_matrix<T, X>::elem_is_too_small(std::vector<indexed_value<T>> & row_chunk, indexed_value<T> & iv, int c_partial_pivoting) {
     if (&iv == &row_chunk[0]) {
         return false; // the max element is at the head
     }
diff --git a/src/util/lp/sparse_matrix.h b/src/util/lp/sparse_matrix.h
index 5aa72122b..61ceb6145 100644
--- a/src/util/lp/sparse_matrix.h
+++ b/src/util/lp/sparse_matrix.h
@@ -7,7 +7,6 @@
 
 #pragma once
 #include <vector>
-#include "util/numerics/float.h"
 #include "util/lp/permutation_matrix.h"
 #include <unordered_map>
 #include "util/lp/static_matrix.h"
@@ -140,7 +139,7 @@ public:
 
 
 
-    unsigned dimension() const {return m_row_permutation.size();}
+    unsigned dimension() const {return static_cast<unsigned>(m_row_permutation.size());}
 
 #ifdef LEAN_DEBUG
     unsigned row_count() const {return dimension();}
@@ -316,7 +315,7 @@ public:
 
     void recover_pivot_queue(std::vector<upair> & rejected_pivots);
 
-    int elem_is_too_small(unsigned i, unsigned j, const T & c_partial_pivoting);
+    int elem_is_too_small(unsigned i, unsigned j, int c_partial_pivoting);
 
     bool remove_row_from_active_pivots_and_shorten_columns(unsigned row);
 
@@ -336,12 +335,12 @@ public:
 
     bool pivot_queue_is_correct_after_pivoting(int k);
 
-    bool get_pivot_for_column(unsigned &i, unsigned &j, T const & c_partial_pivoting, unsigned k);
+    bool get_pivot_for_column(unsigned &i, unsigned &j, int c_partial_pivoting, unsigned k);
 
-    bool elem_is_too_small(std::vector<indexed_value<T>> & row_chunk, indexed_value<T> & iv, T const & c_partial_pivoting);
+    bool elem_is_too_small(std::vector<indexed_value<T>> & row_chunk, indexed_value<T> & iv, int c_partial_pivoting);
 
     unsigned number_of_non_zeroes_in_row(unsigned row) const {
-        return m_rows[row].size();
+        return static_cast<unsigned>(m_rows[row].size());
     }
 
     unsigned number_of_non_zeroes_in_column(unsigned col) const {
diff --git a/src/util/lp/sparse_matrix_instances.cpp b/src/util/lp/sparse_matrix_instances.cpp
index b8cf51a35..5a57e2648 100644
--- a/src/util/lp/sparse_matrix_instances.cpp
+++ b/src/util/lp/sparse_matrix_instances.cpp
@@ -16,7 +16,7 @@ template void sparse_matrix<double, double>::fill_eta_matrix(unsigned int, eta_m
 template const double & sparse_matrix<double, double>::get(unsigned int, unsigned int) const;
 template unsigned sparse_matrix<double, double>::get_number_of_nonzeroes() const;
 template unsigned sparse_matrix<double, double>::get_number_of_nonzeroes_below_row(unsigned int) const;
-template bool sparse_matrix<double, double>::get_pivot_for_column(unsigned int&, unsigned int&, double const&, unsigned int);
+template bool sparse_matrix<double, double>::get_pivot_for_column(unsigned int&, unsigned int&, int, unsigned int);
 template unsigned sparse_matrix<double, double>::lowest_row_in_column(unsigned int);
 template bool sparse_matrix<double, double>::pivot_row_to_row(unsigned int, double, unsigned int, lp_settings&);
 template bool sparse_matrix<double, double>::pivot_with_eta(unsigned int, eta_matrix<double, double>*, lp_settings&);
@@ -30,18 +30,13 @@ template bool sparse_matrix<double, double>::shorten_active_matrix(unsigned int,
 template void sparse_matrix<double, double>::solve_y_U(std::vector<double>&) const;
 template sparse_matrix<double, double>::sparse_matrix(static_matrix<double, double> const&, std::vector<unsigned int>&);
 template sparse_matrix<double, double>::sparse_matrix(unsigned int);
-template float const & sparse_matrix<float, float>::get(unsigned int, unsigned int) const;
-template bool sparse_matrix<float, float>::pivot_row_to_row(unsigned int, float, unsigned int, lp_settings&);
-template void  sparse_matrix<float, float>::replace_column(unsigned int, indexed_vector<float>&, lp_settings&);
-template void     sparse_matrix<float, float>::set(unsigned int, unsigned int, float);
-template sparse_matrix<float, float>::sparse_matrix(unsigned int);
 template void     sparse_matrix<mpq, mpq>::add_new_element(unsigned int, unsigned int, mpq);
 template void     sparse_matrix<mpq, mpq>::divide_row_by_constant(unsigned int, mpq&, lp_settings&);
 template void     sparse_matrix<mpq, mpq>::fill_eta_matrix(unsigned int, eta_matrix<mpq, mpq>**);
 template  mpq const & sparse_matrix<mpq, mpq>::get(unsigned int, unsigned int) const;
 template unsigned sparse_matrix<mpq, mpq>::get_number_of_nonzeroes() const;
 template unsigned sparse_matrix<mpq, mpq>::get_number_of_nonzeroes_below_row(unsigned int) const;
-template bool sparse_matrix<mpq, mpq>::get_pivot_for_column(unsigned int&, unsigned int&, mpq const&, unsigned int);
+template bool sparse_matrix<mpq, mpq>::get_pivot_for_column(unsigned int&, unsigned int&, int, unsigned int);
 template unsigned sparse_matrix<mpq, mpq>::lowest_row_in_column(unsigned int);
 template bool sparse_matrix<mpq, mpq>::pivot_with_eta(unsigned int, eta_matrix<mpq, mpq>*, lp_settings&);
 template void  sparse_matrix<mpq, mpq>::prepare_for_factorization();
@@ -58,7 +53,7 @@ template void     sparse_matrix<mpq, numeric_pair<mpq>>::fill_eta_matrix(unsigne
 template const mpq & sparse_matrix<mpq, numeric_pair<mpq>>::get(unsigned int, unsigned int) const;
 template unsigned sparse_matrix<mpq, numeric_pair<mpq>>::get_number_of_nonzeroes() const;
 template unsigned sparse_matrix<mpq, numeric_pair<mpq>>::get_number_of_nonzeroes_below_row(unsigned int) const;
-template bool sparse_matrix<mpq, numeric_pair<mpq>>::get_pivot_for_column(unsigned int&, unsigned int&, mpq const&, unsigned int);
+template bool sparse_matrix<mpq, numeric_pair<mpq>>::get_pivot_for_column(unsigned int&, unsigned int&, int, unsigned int);
 template unsigned sparse_matrix<mpq, numeric_pair<mpq>>::lowest_row_in_column(unsigned int);
 template bool sparse_matrix<mpq, numeric_pair<mpq>>::pivot_with_eta(unsigned int, eta_matrix<mpq, numeric_pair<mpq> >*, lp_settings&);
 template void     sparse_matrix<mpq, numeric_pair<mpq>>::prepare_for_factorization();
diff --git a/src/util/lp/sparse_vector.h b/src/util/lp/sparse_vector.h
index 421f5ccd4..6e69653ea 100644
--- a/src/util/lp/sparse_vector.h
+++ b/src/util/lp/sparse_vector.h
@@ -7,23 +7,16 @@
 
 #pragma once
 #include <vector>
+#include <utility>
 #include "util/debug.h"
-#include "util/pair.h"
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
-#include "util/numerics/mpq.h"
-#include "util/numerics/mpz.h"
-#include "util/numerics/mpbq.h"
-#include "util/numerics/double.h"
-#include "util/numerics/float.h"
-#include "util/numerics/mpfp.h"
+#include "util/lp/lp_utils.h"
 #include "util/lp/lp_settings.h"
 namespace lean {
 
 template <typename T>
 class sparse_vector {
 public:
-    std::vector<pair<unsigned, T>>  m_data;
+    std::vector<std::pair<unsigned, T>>  m_data;
     void push_back(unsigned index, T val) {
         m_data.emplace_back(index, val);
     }
diff --git a/src/util/lp/square_dense_submatrix.h b/src/util/lp/square_dense_submatrix.h
index 9ffdc96c3..509485801 100644
--- a/src/util/lp/square_dense_submatrix.h
+++ b/src/util/lp/square_dense_submatrix.h
@@ -7,7 +7,6 @@
 
 #pragma once
 #include <vector>
-#include "util/numerics/float.h"
 #include "util/lp/permutation_matrix.h"
 #include <unordered_map>
 #include "util/lp/static_matrix.h"
diff --git a/src/util/lp/static_matrix.cpp b/src/util/lp/static_matrix.cpp
index eece03688..ed64302df 100644
--- a/src/util/lp/static_matrix.cpp
+++ b/src/util/lp/static_matrix.cpp
@@ -153,19 +153,19 @@ template <typename T, typename X>    void static_matrix<T, X>::set(unsigned row,
     if (numeric_traits<T>::is_zero(val)) return;
     lean_assert(row < row_count() && col < column_count());
 #ifdef LEAN_DEBUG
-    pair<unsigned, unsigned> p(row, col);
+    std::pair<unsigned, unsigned> p(row, col);
     lean_assert(m_domain.find(p) == m_domain.end());
     m_domain.insert(p);
 #endif
     auto & r = m_rows[row];
-    unsigned offs_in_cols = m_columns[col].size();
-    m_columns[col].push_back(make_column_cell(row, r.size(), val));
+    unsigned offs_in_cols = static_cast<unsigned>(m_columns[col].size());
+    m_columns[col].push_back(make_column_cell(row, static_cast<unsigned>(r.size()), val));
     r.push_back(make_row_cell(col, offs_in_cols, val));
 }
 
 template <typename T, typename X>
-std::set<pair<unsigned, unsigned>>  static_matrix<T, X>::get_domain() {
-    std::set<pair<unsigned, unsigned>> ret;
+std::set<std::pair<unsigned, unsigned>>  static_matrix<T, X>::get_domain() {
+    std::set<std::pair<unsigned, unsigned>> ret;
     for (unsigned i = 0; i < m_rows.size(); i++) {
         for (auto it : m_rows[i]) {
             ret.insert(std::make_pair(i, it.m_j));
@@ -255,15 +255,15 @@ template <typename T, typename X>    void static_matrix<T, X>::check_consistency
     std::unordered_map<std::pair<unsigned, unsigned>, T> by_rows;
     for (int i = 0; i < m_rows.size(); i++){
         for (auto & t : m_rows[i]) {
-            pair<unsigned, unsigned> p(i, t.m_j);
+            std::pair<unsigned, unsigned> p(i, t.m_j);
             lean_assert(by_rows.find(p) == by_rows.end());
             by_rows[p] = t.get_val();
         }
     }
-    std::unordered_map<pair<unsigned, unsigned>, T> by_cols;
+    std::unordered_map<std::pair<unsigned, unsigned>, T> by_cols;
     for (int i = 0; i < m_columns.size(); i++){
         for (auto & t : m_columns[i]) {
-            pair<unsigned, unsigned> p(t.m_i, i);
+            std::pair<unsigned, unsigned> p(t.m_i, i);
             lean_assert(by_cols.find(p) == by_cols.end());
             by_cols[p] = get_value_of_column_cell(t);
         }
diff --git a/src/util/lp/static_matrix.h b/src/util/lp/static_matrix.h
index e7ceae07e..9d44c4ebf 100644
--- a/src/util/lp/static_matrix.h
+++ b/src/util/lp/static_matrix.h
@@ -10,14 +10,6 @@
 #include <set>
 #include <unordered_map>
 #include <utility>
-#include "util/numerics/numeric_traits.h"
-#include "util/numerics/xnumeral.h"
-#include "util/numerics/mpq.h"
-#include "util/numerics/mpz.h"
-#include "util/numerics/mpbq.h"
-#include "util/numerics/double.h"
-#include "util/numerics/float.h"
-#include "util/numerics/mpfp.h"
 #include "util/lp/sparse_vector.h"
 #include "util/lp/indexed_vector.h"
 #include "util/lp/permutation_matrix.h"
@@ -53,7 +45,7 @@ class static_matrix
 #endif
 {
 #ifdef LEAN_DEBUG
-    std::set<pair<unsigned, unsigned>> m_domain;
+    std::set<std::pair<unsigned, unsigned>> m_domain;
 #endif
 public:
     typedef std::vector<row_cell<T>> row_strip;
@@ -103,9 +95,9 @@ public:
 
     void init_empty_matrix(unsigned m, unsigned n);
 
-    unsigned row_count() const { return m_rows.size(); }
+    unsigned row_count() const { return static_cast<unsigned>(m_rows.size()); }
 
-    unsigned column_count() const { return m_columns.size(); }
+    unsigned column_count() const { return static_cast<unsigned>(m_columns.size()); }
     template <typename L>
     L dot_product_with_row(unsigned row, const std::vector<L> & w);;
 
@@ -144,7 +136,7 @@ public:
 
     ref operator()(unsigned row, unsigned col) { return ref(*this, row, col); }
 
-    std::set<pair<unsigned, unsigned>>  get_domain();
+    std::set<std::pair<unsigned, unsigned>>  get_domain();
 
     void copy_column_to_vector (unsigned j, indexed_vector<T> & v) const;
 
diff --git a/src/util/lp/static_matrix_instances.cpp b/src/util/lp/static_matrix_instances.cpp
index ba6d3030e..51dc564d5 100644
--- a/src/util/lp/static_matrix_instances.cpp
+++ b/src/util/lp/static_matrix_instances.cpp
@@ -7,6 +7,7 @@
 #include <vector>
 #include <memory>
 #include <set>
+#include <utility>
 #include "util/lp/static_matrix.cpp"
 #include "util/lp/lp_core_solver_base.h"
 #include "util/lp/lp_dual_core_solver.h"
@@ -29,7 +30,7 @@ template void static_matrix<double, double>::copy_column_to_vector(unsigned int,
 template void static_matrix<double, double>::divide_row_by_constant(unsigned int, double const&);
 template double static_matrix<double, double>::dot_product_with_column(std::vector<double, std::allocator<double> > const&, unsigned int) const;
 template double static_matrix<double, double>::get_balance() const;
-template std::set<pair<unsigned, unsigned>> static_matrix<double, double>::get_domain();
+    template std::set<std::pair<unsigned, unsigned>> static_matrix<double, double>::get_domain();
 template double static_matrix<double, double>::get_elem(unsigned int, unsigned int) const;
 template double static_matrix<double, double>::get_max_abs_in_column(unsigned int) const;
 template double static_matrix<double, double>::get_min_abs_in_column(unsigned int) const;
diff --git a/src/util/lp/tail_matrix.h b/src/util/lp/tail_matrix.h
index 480919561..4db1c93fe 100644
--- a/src/util/lp/tail_matrix.h
+++ b/src/util/lp/tail_matrix.h
@@ -9,6 +9,7 @@
 #include <vector>
 #include "util/lp/indexed_vector.h"
 #include "util/lp/matrix.h"
+#include "util/lp/lp_settings.h"
 // These matrices appear at the end of the list
 
 namespace lean {