3
0
Fork 0
mirror of https://github.com/Z3Prover/z3 synced 2025-04-24 17:45:32 +00:00
Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>
This commit is contained in:
Nikolaj Bjorner 2022-05-10 15:48:06 -07:00
parent 54648f6b50
commit 5ca3bc3212
3 changed files with 54 additions and 76 deletions

View file

@ -136,6 +136,7 @@ namespace simplex {
svector<int> m_var_pos; // temporary map from variables to positions in row
unsigned_vector m_var_pos_idx; // indices in m_var_pos
stats m_stats;
scoped_numeral m_zero;
bool well_formed_row(unsigned row_id) const;
bool well_formed_column(unsigned column_id) const;
@ -144,7 +145,7 @@ namespace simplex {
public:
sparse_matrix(manager& _m): m(_m) {}
sparse_matrix(manager& _m): m(_m), m_zero(m) {}
~sparse_matrix();
void reset();
@ -216,6 +217,7 @@ namespace simplex {
unsigned column_size(var_t v) const { return m_columns[v].size(); }
unsigned num_vars() const { return m_columns.size(); }
unsigned num_rows() const { return m_rows.size(); }
class col_iterator {
friend class sparse_matrix;
@ -306,11 +308,11 @@ namespace simplex {
all_rows get_rows() { return all_rows(*this); }
numeral& get_coeff(row r, unsigned v) {
numeral const& get_coeff(row r, unsigned v) {
for (auto & [coeff, u] : get_row(r))
if (u == v)
return coeff;
throw default_exception("variable not in row");
return m_zero;
}

View file

@ -26,73 +26,44 @@ namespace simplex {
class sparse_matrix_ops {
public:
static void kernel(sparse_matrix<mpq_ext>& M, vector<vector<rational>>& K) {
mpq_ext::numeral coeff;
rational D;
vector<unsigned> d, c;
unsigned m = M.num_vars();
unsigned n = M.num_vars(), m = M.num_rows();
auto& mgr = M.get_manager();
for (unsigned v = 0; v < m; ++v)
c.push_back(0);
c.resize(m, 0u);
d.resize(n, 0u);
unsigned nullity = 0, rank = 0;
for (auto const& row : M.get_rows()) {
// scan for non-zero variable in row
bool found = false;
d.push_back(0);
++nullity;
for (auto& [coeff1, v] : M.get_row(row)) {
if (mpq_manager<false>::is_zero(coeff1))
for (unsigned k = 0; k < n; ++k) {
d[k] = 0;
for (auto [row, row_entry] : M.get_rows(k)) {
if (c[row.id()] != 0)
continue;
auto& m_jk = row_entry->m_coeff;
if (mpq_manager<false>::is_zero(m_jk))
continue;
if (c[v] != 0)
continue;
--nullity;
++rank;
d.back() = v + 1;
c[v] = row.id() + 1;
D = rational(-1) / coeff1;
mgr.set(coeff1, mpq(-1));
// eliminate v from other rows.
for (auto [row2, row_entry2] : M.get_rows(v)) {
if (row.id() >= row2.id())
continue;
mpq & m_js = row_entry2->m_coeff;
mgr.set(m_js, (D * m_js).to_mpq());
}
for (auto& [m_ik, w] : M.get_row(row)) {
if (v == w)
continue;
D = m_ik;
mgr.set(m_ik, mpq(0));
for (auto [row2, row_entry2] : M.get_rows(w)) {
if (row.id() >= row2.id())
continue;
auto& m_js = M.get_coeff(row2, v);
auto & m_is = row_entry2->m_coeff;
mgr.set(m_is, (m_is + D * m_js).to_mpq());
}
}
D = rational(-1) / m_jk;
M.mul(row, D.to_mpq());
for (auto [row_i, row_i_entry] : M.get_rows(k)) {
if (row_i.id() == row.id())
continue;
auto& m_ik = row_i_entry->m_coeff;
// row_i += m_ik * row
M.add(row_i, m_ik, row);
}
c[row.id()] = k + 1;
d[k] = row.id() + 1;
break;
}
}
}
mgr.del(coeff);
std::cout << "nullity " << nullity << "\n";
std::cout << "rank " << rank << "\n";
unsigned_vector pivots(rank, 0u);
unsigned_vector nonpivots(nullity, 0u);
// TODO: extract kernel using d
for (unsigned k = 0; k < d.size(); ++k) {
for (unsigned k = 0; k < n; ++k) {
if (d[k] != 0)
continue;
K.push_back(vector<rational>());
for (unsigned i = 0; i < d.size(); ++i) {
for (unsigned i = 0; i < n; ++i) {
if (d[i] > 0) {
// row r = row(i);
// K.back().push_back(M[d[i]-1][k]);
auto r = sparse_matrix<mpq_ext>::row(d[i]-1);
K.back().push_back(rational(M.get_coeff(r, k)));
}
else if (i == k)
K.back().push_back(rational(1));
@ -102,6 +73,7 @@ namespace simplex {
}
}
};

View file

@ -31,23 +31,23 @@ static vector<R> vec(int i, int j, int k) {
return nv;
}
// static vector<R> vec(int i, int j, int k, int l) {
// vector<R> nv = vec(i, j, k);
// nv.push_back(R(l));
// return nv;
// }
static vector<R> vec(int i, int j, int k, int l) {
vector<R> nv = vec(i, j, k);
nv.push_back(R(l));
return nv;
}
/// static vector<R> vec(int i, int j, int k, int l, int x) {
/// vector<R> nv = vec(i, j, k, l);
/// nv.push_back(R(x));
/// return nv;
/// }
static vector<R> vec(int i, int j, int k, int l, int x) {
vector<R> nv = vec(i, j, k, l);
nv.push_back(R(x));
return nv;
}
// static vector<R> vec(int i, int j, int k, int l, int x, int y) {
// vector<R> nv = vec(i, j, k, l, x);
// nv.push_back(R(y));
// return nv;
// }
static vector<R> vec(int i, int j, int k, int l, int x, int y) {
vector<R> nv = vec(i, j, k, l, x);
nv.push_back(R(y));
return nv;
}
// static vector<R> vec(int i, int j, int k, int l, int x, int y, int z) {
// vector<R> nv = vec(i, j, k, l, x, y);
@ -133,7 +133,7 @@ static void test4() {
}
static void add(qmatrix& m, vector<R> const& v) {
m.ensure_var(v.size());
m.ensure_var(v.size()-1);
auto r = m.mk_row();
for (unsigned u = 0; u < v.size(); ++u)
m.add_var(r, v[u].to_mpq(), u);
@ -142,12 +142,16 @@ static void add(qmatrix& m, vector<R> const& v) {
static void test5() {
unsynch_mpq_manager m;
qmatrix M(m);
add(M, vec(1, 2, 3));
add(M, vec(2, 2, 4));
add(M, vec(1, 0, -3, 0, 2, -8));
add(M, vec(0, 1, 5, 0, -1, 4));
add(M, vec(0, 0, 0, 1, 7, -9));
add(M, vec(0, 0, 0, 0, 0, 0));
M.display(std::cout);
vector<vector<rational>> K;
kernel(M, K);
std::cout << "after\n";
for (auto const& v : K)
std::cout << v << "\n";
M.display(std::cout);
}