3
0
Fork 0
mirror of https://github.com/Z3Prover/z3 synced 2025-10-02 22:19:30 +00:00

replace lean to lp

Signed-off-by: Lev Nachmanson <levnach@microsoft.com>
This commit is contained in:
Lev Nachmanson 2017-07-10 11:06:37 -07:00 committed by Lev Nachmanson
parent db0a3f4358
commit d41c65a4f9
72 changed files with 1334 additions and 1213 deletions

View file

@ -25,7 +25,7 @@ Revision History:
#include "util/debug.h"
#include "util/lp/lu.h"
namespace lp {
#ifdef Z3DEBUG
#ifdef LEAN_DEBUG
template <typename T, typename X> // print the nr x nc submatrix at the top left corner
void print_submatrix(sparse_matrix<T, X> & m, unsigned mr, unsigned nc, std::ostream & out) {
vector<vector<std::string>> A;
@ -137,29 +137,29 @@ lu<T, X>::lu(static_matrix<T, X> const & A,
m_failure(false),
m_row_eta_work_vector(A.row_count()),
m_refactor_counter(0) {
SASSERT(!(numeric_traits<T>::precise() && settings.use_tableau()));
#ifdef Z3DEBUG
lp_assert(!(numeric_traits<T>::precise() && settings.use_tableau()));
#ifdef LEAN_DEBUG
debug_test_of_basis(A, basis);
#endif
++m_settings.st().m_num_factorizations;
create_initial_factorization();
#ifdef Z3DEBUG
// SASSERT(check_correctness());
#ifdef LEAN_DEBUG
// lp_assert(check_correctness());
#endif
}
template <typename T, typename X>
void lu<T, X>::debug_test_of_basis(static_matrix<T, X> const & A, vector<unsigned> & basis) {
std::set<unsigned> set;
for (unsigned i = 0; i < A.row_count(); i++) {
SASSERT(basis[i]< A.column_count());
lp_assert(basis[i]< A.column_count());
set.insert(basis[i]);
}
SASSERT(set.size() == A.row_count());
lp_assert(set.size() == A.row_count());
}
template <typename T, typename X>
void lu<T, X>::solve_By(indexed_vector<X> & y) {
SASSERT(false); // not implemented
lp_assert(false); // not implemented
// init_vector_y(y);
// solve_By_when_y_is_ready(y);
}
@ -292,20 +292,20 @@ void lu<T, X>::solve_yB(vector<T>& y) {
template <typename T, typename X>
void lu<T, X>::solve_yB_indexed(indexed_vector<T>& y) {
SASSERT(y.is_OK());
lp_assert(y.is_OK());
// first solve yU = cb*R(-1)
m_R.apply_reverse_from_right_to_T(y); // got y = cb*R(-1)
SASSERT(y.is_OK());
lp_assert(y.is_OK());
m_U.solve_y_U_indexed(y, m_settings); // got y*U=cb*R(-1)
SASSERT(y.is_OK());
lp_assert(y.is_OK());
m_Q.apply_reverse_from_right_to_T(y);
SASSERT(y.is_OK());
lp_assert(y.is_OK());
for (auto e = m_tail.rbegin(); e != m_tail.rend(); ++e) {
#ifdef Z3DEBUG
(*e)->set_number_of_columns(m_dim);
#endif
(*e)->apply_from_right(y);
SASSERT(y.is_OK());
lp_assert(y.is_OK());
}
}
@ -319,8 +319,8 @@ void lu<T, X>::add_delta_to_solution(const vector<T>& yc, vector<T>& y){
template <typename T, typename X>
void lu<T, X>::add_delta_to_solution_indexed(indexed_vector<T>& y) {
// the delta sits in m_y_copy, put result into y
SASSERT(y.is_OK());
SASSERT(m_y_copy.is_OK());
lp_assert(y.is_OK());
lp_assert(m_y_copy.is_OK());
m_ii.clear();
m_ii.resize(y.data_size());
for (unsigned i : y.m_index)
@ -330,7 +330,7 @@ void lu<T, X>::add_delta_to_solution_indexed(indexed_vector<T>& y) {
if (m_ii[i] == 0)
m_ii.set_value(1, i);
}
SASSERT(m_ii.is_OK());
lp_assert(m_ii.is_OK());
y.m_index.clear();
for (unsigned i : m_ii.m_index) {
@ -341,7 +341,7 @@ void lu<T, X>::add_delta_to_solution_indexed(indexed_vector<T>& y) {
v = zero_of_type<T>();
}
SASSERT(y.is_OK());
lp_assert(y.is_OK());
}
template <typename T, typename X>
@ -358,7 +358,7 @@ void lu<T, X>::find_error_of_yB_indexed(const indexed_vector<T>& y, const vector
// it is a non efficient version
indexed_vector<T> yc = m_y_copy;
yc.m_index.clear();
SASSERT(!numeric_traits<T>::precise());
lp_assert(!numeric_traits<T>::precise());
{
vector<unsigned> d_basis(y.m_data.size());
@ -379,10 +379,10 @@ void lu<T, X>::find_error_of_yB_indexed(const indexed_vector<T>& y, const vector
}
}
#endif
SASSERT(m_ii.is_OK());
lp_assert(m_ii.is_OK());
m_ii.clear();
m_ii.resize(y.data_size());
SASSERT(m_y_copy.is_OK());
lp_assert(m_y_copy.is_OK());
// put the error into m_y_copy
for (auto k : y.m_index) {
auto & row = m_A.m_rows[k];
@ -414,7 +414,7 @@ void lu<T, X>::find_error_of_yB_indexed(const indexed_vector<T>& y, const vector
m_y_copy.set_value(v, k);
}
}
SASSERT(m_y_copy.is_OK());
lp_assert(m_y_copy.is_OK());
}
@ -430,31 +430,31 @@ void lu<T, X>::solve_yB_with_error_check_indexed(indexed_vector<T> & y, const ve
solve_yB_indexed(y);
} else {
solve_yB(y.m_data);
y.restore_index_and_clean_from_data();
y.restore_index_and_clp_from_data();
}
return;
}
SASSERT(m_y_copy.is_OK());
SASSERT(y.is_OK());
lp_assert(m_y_copy.is_OK());
lp_assert(y.is_OK());
if (y.m_index.size() * ratio_of_index_size_to_all_size<T>() < m_A.column_count()) {
m_y_copy = y;
solve_yB_indexed(y);
SASSERT(y.is_OK());
lp_assert(y.is_OK());
if (y.m_index.size() * ratio_of_index_size_to_all_size<T>() >= m_A.column_count()) {
find_error_of_yB(m_y_copy.m_data, y.m_data, basis);
solve_yB(m_y_copy.m_data);
add_delta_to_solution(m_y_copy.m_data, y.m_data);
y.restore_index_and_clean_from_data();
y.restore_index_and_clp_from_data();
m_y_copy.clear_all();
} else {
find_error_of_yB_indexed(y, heading, settings); // this works with m_y_copy
solve_yB_indexed(m_y_copy);
add_delta_to_solution_indexed(y);
}
SASSERT(m_y_copy.is_OK());
lp_assert(m_y_copy.is_OK());
} else {
solve_yB_with_error_check(y.m_data, basis);
y.restore_index_and_clean_from_data();
y.restore_index_and_clp_from_data();
}
}
@ -504,7 +504,7 @@ template <typename T, typename X>
void lu<T, X>::perform_transformations_on_w(indexed_vector<T>& w) {
apply_lp_list_to_w(w);
m_Q.apply_reverse_from_left(w);
// TBD does not compile: SASSERT(numeric_traits<T>::precise() || check_vector_for_small_values(w, m_settings));
// TBD does not compile: lp_assert(numeric_traits<T>::precise() || check_vector_for_small_values(w, m_settings));
}
// see Chvatal 24.3
@ -518,7 +518,7 @@ template <typename T, typename X>
void lu<T, X>::apply_lp_list_to_w(indexed_vector<T> & w) {
for (unsigned i = 0; i < m_tail.size(); i++) {
m_tail[i]->apply_from_left_to_T(w, m_settings);
// TBD does not compile: SASSERT(check_vector_for_small_values(w, m_settings));
// TBD does not compile: lp_assert(check_vector_for_small_values(w, m_settings));
}
}
template <typename T, typename X>
@ -610,7 +610,7 @@ void lu<T, X>::check_apply_lp_lists_to_w(T * w) {
permutation_matrix<T, X> qr = m_Q.get_reverse();
apply_to_vector(qr, w);
for (int i = m_dim - 1; i >= 0; i--) {
SASSERT(abs(w[i] - w[i]) < 0.0000001);
lp_assert(abs(w[i] - w[i]) < 0.0000001);
}
}
@ -655,7 +655,7 @@ bool lu<T, X>::is_correct(const vector<unsigned>& basis) {
#ifdef Z3DEBUG
template <typename T, typename X>
dense_matrix<T, X> lu<T, X>::tail_product() {
SASSERT(tail_size() > 0);
lp_assert(tail_size() > 0);
dense_matrix<T, X> left_side = permutation_matrix<T, X>(m_dim);
for (unsigned i = 0; i < tail_size(); i++) {
matrix<T, X>* lp = get_lp_matrix(i);
@ -705,8 +705,8 @@ template <typename T, typename X>
bool lu<T, X>::all_columns_and_rows_are_active() {
unsigned i = m_dim;
while (i--) {
SASSERT(m_U.col_is_active(i));
SASSERT(m_U.row_is_active(i));
lp_assert(m_U.col_is_active(i));
lp_assert(m_U.row_is_active(i));
}
return true;
}
@ -748,9 +748,9 @@ void lu<T, X>::create_initial_factorization(){
}
}
if (j == m_dim) {
// TBD does not compile: SASSERT(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// SASSERT(is_correct());
// SASSERT(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// TBD does not compile: lp_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// lp_assert(is_correct());
// lp_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
return;
}
j++;
@ -763,12 +763,12 @@ void lu<T, X>::create_initial_factorization(){
}
}
m_dense_LU->update_parent_matrix(m_settings);
SASSERT(m_dense_LU->is_L_matrix());
lp_assert(m_dense_LU->is_L_matrix());
m_dense_LU->conjugate_by_permutation(m_Q);
push_matrix_to_tail(m_dense_LU);
m_refactor_counter = 0;
// SASSERT(is_correct());
// SASSERT(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// lp_assert(is_correct());
// lp_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
}
template <typename T, typename X>
@ -795,7 +795,7 @@ void lu<T, X>::scan_last_row_to_work_vector(unsigned lowest_row_of_the_bump) {
vector<indexed_value<T>> & last_row_vec = m_U.get_row_values(m_U.adjust_row(lowest_row_of_the_bump));
for (auto & iv : last_row_vec) {
if (is_zero(iv.m_value)) continue;
SASSERT(!m_settings.abs_val_is_smaller_than_drop_tolerance(iv.m_value));
lp_assert(!m_settings.abs_val_is_smaller_than_drop_tolerance(iv.m_value));
unsigned adjusted_col = m_U.adjust_column_inverse(iv.m_index);
if (adjusted_col < lowest_row_of_the_bump) {
m_row_eta_work_vector.set_value(-iv.m_value, adjusted_col);
@ -816,14 +816,14 @@ void lu<T, X>::pivot_and_solve_the_system(unsigned replaced_column, unsigned low
vector<indexed_value<T>> & row = m_U.get_row_values(aj);
for (auto & iv : row) {
unsigned col = m_U.adjust_column_inverse(iv.m_index);
SASSERT(col >= j || numeric_traits<T>::is_zero(iv.m_value));
lp_assert(col >= j || numeric_traits<T>::is_zero(iv.m_value));
if (col == j) continue;
if (numeric_traits<T>::is_zero(iv.m_value)) {
continue;
}
// the -v is for solving the system ( to zero the last row), and +v is for pivoting
T delta = col < lowest_row_of_the_bump? -v * iv.m_value: v * iv.m_value;
SASSERT(numeric_traits<T>::is_zero(delta) == false);
lp_assert(numeric_traits<T>::is_zero(delta) == false);
@ -900,16 +900,16 @@ void lu<T, X>::replace_column(T pivot_elem_for_checking, indexed_vector<T> & w,
push_matrix_to_tail(row_eta);
}
calculate_Lwave_Pwave_for_bump(replaced_column, lowest_row_of_the_bump);
// SASSERT(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// SASSERT(w.is_OK() && m_row_eta_work_vector.is_OK());
// lp_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// lp_assert(w.is_OK() && m_row_eta_work_vector.is_OK());
}
template <typename T, typename X>
void lu<T, X>::calculate_Lwave_Pwave_for_bump(unsigned replaced_column, unsigned lowest_row_of_the_bump){
T diagonal_elem;
if (replaced_column < lowest_row_of_the_bump) {
diagonal_elem = m_row_eta_work_vector[lowest_row_of_the_bump];
// SASSERT(m_row_eta_work_vector.is_OK());
m_U.set_row_from_work_vector_and_clean_work_vector_not_adjusted(m_U.adjust_row(lowest_row_of_the_bump), m_row_eta_work_vector, m_settings);
// lp_assert(m_row_eta_work_vector.is_OK());
m_U.set_row_from_work_vector_and_clp_work_vector_not_adjusted(m_U.adjust_row(lowest_row_of_the_bump), m_row_eta_work_vector, m_settings);
} else {
diagonal_elem = m_U(lowest_row_of_the_bump, lowest_row_of_the_bump); // todo - get it more efficiently
}
@ -919,7 +919,7 @@ void lu<T, X>::calculate_Lwave_Pwave_for_bump(unsigned replaced_column, unsigned
}
calculate_Lwave_Pwave_for_last_row(lowest_row_of_the_bump, diagonal_elem);
// SASSERT(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
// lp_assert(m_U.is_upper_triangular_and_maximums_are_set_correctly_in_rows(m_settings));
}
template <typename T, typename X>
@ -945,8 +945,8 @@ void init_factorization(lu<T, X>* & factorization, static_matrix<T, X> & m_A, ve
#ifdef Z3DEBUG
template <typename T, typename X>
dense_matrix<T, X> get_B(lu<T, X>& f, const vector<unsigned>& basis) {
SASSERT(basis.size() == f.dimension());
SASSERT(basis.size() == f.m_U.dimension());
lp_assert(basis.size() == f.dimension());
lp_assert(basis.size() == f.m_U.dimension());
dense_matrix<T, X> B(f.dimension(), f.dimension());
for (unsigned i = 0; i < f.dimension(); i++)
for (unsigned j = 0; j < f.dimension(); j++)