mirrors/z3
3
0
Fork 0
mirror of https://github.com/Z3Prover/z3 synced 2026-02-08 10:07:59 +00:00
Code Activity

t

Browse source 
Signed-off-by: Lev Nachmanson <levnach@hotmail.com>
This commit is contained in:
Lev Nachmanson 2026-01-28 06:19:45 -10:00
parent 7cfafc133f
commit 5e1c104667
6 changed files with 601 additions and 367 deletions
Download patch file Download diff file Expand all files Collapse all files

620
src/nlsat/levelwise.cpp
View file

@ -105,9 +105,9 @@ namespace nlsat {
}
struct root_function {
scoped_anum val;
scoped_anum val;
indexed_root_expr ire;
unsigned ps_idx; // index in m_level_ps
unsigned ps_idx; // index in m_level_ps
root_function(anum_manager& am, poly* p, unsigned idx, anum const& v, unsigned ps_idx)
: val(am), ire{ p, idx }, ps_idx(ps_idx) { am.set(val, v); }
root_function(root_function&& other) noexcept : val(other.val.m()), ire(other.ire), ps_idx(other.ps_idx) { val = other.val; }
@ -345,10 +345,10 @@ namespace nlsat {
void request_factorized(polynomial_ref const& poly, inv_req req) {
for_each_distinct_factor(poly, [&](polynomial_ref const& f) {
TRACE(lws,
tout << " request_factorized: factor=";
m_pm.display(tout, f);
tout << " at level " << m_pm.max_var(f) << "\n";
);
tout << " request_factorized: factor=";
m_pm.display(tout, f);
tout << " at level " << m_pm.max_var(f) << "\n";
);
request(f.get(), req); // inherit tag across factorization (SMT-RAT appendOnCorrectLevel)
});
}
@ -421,12 +421,12 @@ namespace nlsat {
return best;
}
void add_projections_for(polynomial_ref const& p, unsigned x, polynomial_ref const& nonzero_coeff, bool add_leading_coeff, bool add_discriminant) {
void add_projection_for_poly(polynomial_ref const& p, unsigned x, polynomial_ref const& nonzero_coeff, bool add_leading_coeff, bool add_discriminant) {
TRACE(lws,
tout << " add_projections_for: p=";
m_pm.display(tout, p);
tout << " x=" << x << " add_lc=" << add_leading_coeff << " add_disc=" << add_discriminant << "\n";
);
tout << " add_projection_for_poly: p=";
m_pm.display(tout, p);
tout << " x=" << x << " add_lc=" << add_leading_coeff << " add_disc=" << add_discriminant << "\n";
);
// Line 11 (non-null witness coefficient)
if (nonzero_coeff && !is_const(nonzero_coeff))
request_factorized(nonzero_coeff, inv_req::sign);
@ -439,10 +439,10 @@ namespace nlsat {
polynomial_ref lc(m_pm);
lc = m_pm.coeff(p, x, deg);
TRACE(lws,
tout << " adding lc: ";
m_pm.display(tout, lc);
tout << "\n";
);
tout << " adding lc: ";
m_pm.display(tout, lc);
tout << "\n";
);
request_factorized(lc, inv_req::sign);
}
}
@ -455,14 +455,15 @@ namespace nlsat {
// Compute side_mask: track which side(s) each polynomial appears on
// bit0 = lower (<= sample), bit1 = upper (> sample), 3 = both sides
void compute_side_mask() {
m_side_mask.clear();
if (!is_set(m_l_rf))
return;
m_side_mask.resize(m_level_ps.size(), 0);
anum const& v = sample().value(m_level);
for (auto const& rf : m_rel.m_rfunc) {
if (m_am.compare(rf.val, v) <= 0)
m_side_mask[rf.ps_idx] |= 1;
for (unsigned i = 0; i < m_rel.m_rfunc.size(); i ++) {
unsigned ps_idx = m_rel.m_rfunc[i].ps_idx;
if (i <= m_l_rf)
m_side_mask[ps_idx] |= 1;
else
m_side_mask[rf.ps_idx] |= 2;
m_side_mask[ps_idx] |= 2;
}
}
@ -545,7 +546,6 @@ namespace nlsat {
if (m_rel.m_rfunc.empty() || m_rel.m_pairs.empty())
return;
compute_side_mask();
if (require_leaf)
compute_resultant_graph_degree();
@ -564,28 +564,22 @@ namespace nlsat {
}
}
// Compute noLdcf for sector case
void compute_omit_lc_sector() {
if (!is_set(m_l_rf) || !is_set(m_u_rf)) return;
if (m_rel_mode == spanning_tree)
compute_omit_lc_both_sides(true);
else
compute_omit_lc_both_sides(false);
// m_rel == spanning_tree, only leaves can drop ldcfs, but in the biggest cell all non-bounds are leaves,
compute_omit_lc_both_sides(m_rel_mode == spanning_tree);
}
// Compute noLdcf for section case
void compute_omit_lc_section() {
if (m_rel_mode == biggest_cell)
m_omit_lc.clear(); // no omit_lc for biggest_cell
else // spanning_tree
if (m_rel_mode == spanning_tree)
compute_omit_lc_both_sides(true);
}
// Relation construction heuristics (same intent as previous implementation).
void fill_relation_with_biggest_cell_heuristic() {
void fill_relation_sector_biggest_cell() {
TRACE(lws,
tout << " fill_biggest_cell: m_l_rf=" << m_l_rf << ", m_u_rf=" << m_u_rf << ", rfunc.size=" << m_rel.m_rfunc.size() << "\n";
);
tout << " fill_biggest_cell: m_l_rf=" << m_l_rf << ", m_u_rf=" << m_u_rf << ", rfunc.size=" << m_rel.m_rfunc.size() << "\n";
);
if (is_set(m_l_rf))
for (unsigned j = 0; j < m_l_rf; ++j) {
TRACE(lws, tout << " add_pair(" << j << ", " << m_l_rf << ")\n";);
@ -604,8 +598,8 @@ namespace nlsat {
m_rel.add_pair(m_l_rf, m_u_rf);
}
TRACE(lws,
tout << " fill_biggest_cell done: pairs.size=" << m_rel.m_pairs.size() << "\n";
);
tout << " fill_biggest_cell done: pairs.size=" << m_rel.m_pairs.size() << "\n";
);
}
void fill_relation_pairs_for_section_biggest_cell() {
@ -630,229 +624,146 @@ namespace nlsat {
// Build spanning tree on both-side polynomials using the lemma construction.
// Adds pairs directly to m_rel.m_pairs. Returns true if tree was built.
// K = lower rfunc positions, f = upper rfunc positions
bool build_both_side_spanning_tree() {
void build_both_side_spanning_tree() {
auto const& rfs = m_rel.m_rfunc;
unsigned n = rfs.size();
if (n == 0) return false;
SASSERT(n > 0 && is_set(m_l_rf) && is_set(m_u_rf));
SASSERT(!is_section());
SASSERT(!same_boundary_poly());
// Map ps_idx -> rfunc index on lower/upper side
std_vector<unsigned> lower_pos(m_level_ps.size(), UINT_MAX);
std_vector<unsigned> upper_pos(m_level_ps.size(), UINT_MAX);
// Collect both-side polynomials with their rfunc indices on each side
struct both_info {
unsigned ps_idx;
unsigned lower_rf; // rfunc index on lower side
unsigned upper_rf; // rfunc index on upper side
};
std_vector<both_info> both;
if (is_set(m_l_rf)) {
for (unsigned i = 0; i <= m_l_rf; ++i)
lower_pos[rfs[i].ps_idx] = i;
}
// For sector case: upper side is [m_u_rf, n)
// For section case: upper side is (m_l_rf, n) since m_u_rf is not set
unsigned upper_start = is_set(m_u_rf) ? m_u_rf : (is_set(m_l_rf) ? m_l_rf + 1 : 0);
for (unsigned i = upper_start; i < n; ++i)
upper_pos[rfs[i].ps_idx] = i;
// For sector: lower side is [0, m_l_rf], upper side is [m_u_rf, n)
// Collect both-side polynomial ps_idx's (must have valid positions on both sides)
std_vector<unsigned> both;
// Build map from ps_idx to rfunc indices
std_vector<unsigned> lower_rf(m_level_ps.size(), UINT_MAX);
std_vector<unsigned> upper_rf(m_level_ps.size(), UINT_MAX);
for (unsigned i = 0; i <= m_l_rf; ++i)
lower_rf[rfs[i].ps_idx] = i;
for (unsigned i = m_u_rf; i < n; ++i)
upper_rf[rfs[i].ps_idx] = i;
// Collect both-side polynomials
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
if (m_side_mask[i] == 3 && lower_pos[i] != UINT_MAX && upper_pos[i] != UINT_MAX)
both.push_back(i);
if (m_side_mask[i] == 3) {
SASSERT(lower_rf[i] != UINT_MAX && upper_rf[i] != UINT_MAX);
both.push_back({i, lower_rf[i], upper_rf[i]});
}
}
// Need at least threshold polynomials on both sides
if (both.size() < m_spanning_tree_threshold)
return false;
SASSERT(both.size() >= m_spanning_tree_threshold);
// Sort both by lower_pos (this is K in the lemma)
// Sort by lower_rf (root position on lower side)
std::sort(both.begin(), both.end(),
[&](unsigned a, unsigned b) { return lower_pos[a] < lower_pos[b]; });
[](both_info const& a, both_info const& b) { return a.lower_rf < b.lower_rf; });
// Lemma construction: iteratively connect min(remaining) to element with min f-value in tree
// Build spanning tree using Reaching Orders Lemma:
// Start with element that has minimum upper_rf, root of out-arborescence on the right
std_vector<bool> in_tree(both.size(), false);
// Start with element that has minimum upper_pos (root of out-arborescence on f(K))
unsigned root_idx = 0;
for (unsigned i = 1; i < both.size(); ++i) {
if (upper_pos[both[i]] < upper_pos[both[root_idx]])
if (both[i].upper_rf < both[root_idx].upper_rf)
root_idx = i;
}
in_tree[root_idx] = true;
// Add remaining elements
// Iteratively add remaining elements
for (unsigned added = 1; added < both.size(); ++added) {
// Find minimum lower_pos element not yet in tree (element 'a' in lemma)
// Find min lower_rf element not in tree (first in sorted order)
unsigned a_idx = UINT_MAX;
for (unsigned i = 0; i < both.size(); ++i) {
if (!in_tree[i]) {
a_idx = i;
break; // both is sorted by lower_pos, so first not-in-tree is min
}
}
for (unsigned i = 0; i < both.size(); ++i)
if (!in_tree[i]) { a_idx = i; break; }
SASSERT(a_idx != UINT_MAX);
// Find element in tree with minimum upper_pos (element 'c' in lemma)
// Find element in tree with min upper_rf
unsigned c_idx = UINT_MAX;
unsigned min_upper = UINT_MAX;
for (unsigned i = 0; i < both.size(); ++i) {
if (in_tree[i] && upper_pos[both[i]] < min_upper) {
min_upper = upper_pos[both[i]];
if (in_tree[i] && (c_idx == UINT_MAX || both[i].upper_rf < both[c_idx].upper_rf))
c_idx = i;
}
}
SASSERT(c_idx != UINT_MAX);
// Add edge {a, c} as ps_idx pair directly to m_rel.m_pairs
unsigned ps_a = both[a_idx];
unsigned ps_c = both[c_idx];
// Add edge {a, c} as ps_idx pair
unsigned ps_a = both[a_idx].ps_idx;
unsigned ps_c = both[c_idx].ps_idx;
m_rel.m_pairs.insert({std::min(ps_a, ps_c), std::max(ps_a, ps_c)});
in_tree[a_idx] = true;
}
return true;
}
// Sector spanning tree heuristic:
// 1. Build spanning tree on both-side polynomials
// 2. Extend with lowest_degree for single-side polynomials
void fill_relation_with_spanning_tree_heuristic() {
// Helper: Connect non-tree (single-side) polynomials to their respective boundaries
void connect_non_tree_to_bounds() {
auto const& rfs = m_rel.m_rfunc;
unsigned n = rfs.size();
if (n == 0) return;
// Need side_mask computed
compute_side_mask();
// Build spanning tree on both-side polynomials (adds to m_rel.m_pairs)
if (!build_both_side_spanning_tree()) {
// Not enough valid both-side polys, fall back
fill_relation_with_biggest_cell_heuristic();
return;
}
// Both-side polynomials exist, so both boundaries must be set
SASSERT(is_set(m_l_rf) && is_set(m_u_rf));
// Extend for single-side polynomials using biggest_cell style
// (connect all unconnected to the boundary)
// Lower side: connect all unconnected polynomials to the boundary m_l_rf
// Lower side: connect single-side polys to lower boundary
for (unsigned j = 0; j < m_l_rf; ++j) {
if (m_side_mask[rfs[j].ps_idx] != 3) // single-side poly
if (m_side_mask[rfs[j].ps_idx] != 3)
m_rel.add_pair(j, m_l_rf);
}
// If boundary is not a both-side poly, connect spanning tree to boundary
// using the both-side poly with maximum lower position
if (m_side_mask[rfs[m_l_rf].ps_idx] != 3) {
unsigned max_lower_idx = UINT_MAX;
for (unsigned j = 0; j < m_l_rf; ++j) {
if (m_side_mask[rfs[j].ps_idx] == 3)
max_lower_idx = j; // last one found is maximum
}
if (max_lower_idx != UINT_MAX)
m_rel.add_pair(max_lower_idx, m_l_rf);
}
// Upper side: connect all unconnected polynomials to the boundary m_u_rf
// Upper side: connect single-side polys to upper boundary
for (unsigned j = m_u_rf + 1; j < n; ++j) {
if (m_side_mask[rfs[j].ps_idx] != 3) // single-side poly
if (m_side_mask[rfs[j].ps_idx] != 3)
m_rel.add_pair(m_u_rf, j);
}
// If boundary is not a both-side poly, connect spanning tree to boundary
// using the both-side poly with minimum upper position
if (m_side_mask[rfs[m_u_rf].ps_idx] != 3) {
unsigned min_upper_idx = UINT_MAX;
for (unsigned j = m_u_rf + 1; j < n; ++j) {
if (m_side_mask[rfs[j].ps_idx] == 3) {
min_upper_idx = j; // first one found is minimum
break;
}
}
if (min_upper_idx != UINT_MAX)
m_rel.add_pair(m_u_rf, min_upper_idx);
}
// Helper: Connect spanning tree extremes to boundaries (when boundaries are different polys)
void connect_tree_extremes_to_bounds() {
auto const& rfs = m_rel.m_rfunc;
unsigned n = rfs.size();
// Find max lower both-side poly (closest to lower boundary from below)
unsigned max_lower_both = UINT_MAX;
for (unsigned j = 0; j < m_l_rf; ++j) {
if (m_side_mask[rfs[j].ps_idx] == 3)
max_lower_both = j;
}
// Find min upper both-side poly (closest to upper boundary from above)
unsigned min_upper_both = UINT_MAX;
for (unsigned j = m_u_rf + 1; j < n; ++j) {
if (m_side_mask[rfs[j].ps_idx] == 3) {
min_upper_both = j;
break;
}
}
// Connect tree extremes to boundaries
if (max_lower_both != UINT_MAX)
m_rel.add_pair(max_lower_both, m_l_rf);
if (min_upper_both != UINT_MAX)
m_rel.add_pair(m_u_rf, min_upper_both);
}
void fill_relation_sector_spanning_tree() {
build_both_side_spanning_tree();
connect_non_tree_to_bounds();
connect_tree_extremes_to_bounds(); // otherwise the trees on both sides are connected to bounds already
// Connect lower and upper boundaries
SASSERT(m_l_rf + 1 == m_u_rf);
m_rel.add_pair(m_l_rf, m_u_rf);
}
// Section spanning tree heuristic
void fill_relation_pairs_for_section_spanning_tree() {
auto const& rfs = m_rel.m_rfunc;
unsigned n = rfs.size();
if (n == 0) return;
SASSERT(is_set(m_l_rf));
SASSERT(m_l_rf < n);
// Need side_mask computed
compute_side_mask();
// Build spanning tree on both-side polynomials (adds to m_rel.m_pairs)
if (!build_both_side_spanning_tree()) {
// Not enough valid both-side polys, fall back
fill_relation_pairs_for_section_biggest_cell();
return;
}
// Extend for single-side polynomials using biggest_cell style
// (connect all unconnected to the section m_l_rf)
// Below section: connect all unconnected to section
for (unsigned j = 0; j < m_l_rf; ++j) {
if (m_side_mask[rfs[j].ps_idx] != 3) // single-side poly
m_rel.add_pair(j, m_l_rf);
}
// If section poly is not a both-side poly, connect spanning tree to section
// using the both-side poly with maximum lower position
if (m_side_mask[rfs[m_l_rf].ps_idx] != 3) {
unsigned max_lower_idx = UINT_MAX;
for (unsigned j = 0; j < m_l_rf; ++j) {
if (m_side_mask[rfs[j].ps_idx] == 3)
max_lower_idx = j; // last one found is maximum
}
if (max_lower_idx != UINT_MAX)
m_rel.add_pair(max_lower_idx, m_l_rf);
}
// Above section: connect all unconnected to section
for (unsigned j = m_l_rf + 1; j < n; ++j) {
if (m_side_mask[rfs[j].ps_idx] != 3) // single-side poly
m_rel.add_pair(m_l_rf, j);
}
// If section poly is not a both-side poly, connect spanning tree to section
// using the both-side poly with minimum upper position
if (m_side_mask[rfs[m_l_rf].ps_idx] != 3) {
unsigned min_upper_idx = UINT_MAX;
for (unsigned j = m_l_rf + 1; j < n; ++j) {
if (m_side_mask[rfs[j].ps_idx] == 3) {
min_upper_idx = j; // first one found is minimum
break;
}
}
if (min_upper_idx != UINT_MAX)
m_rel.add_pair(m_l_rf, min_upper_idx);
}
}
void fill_relation_pairs_for_section() {
SASSERT(m_I[m_level].section);
if (m_rel_mode == spanning_tree)
fill_relation_pairs_for_section_spanning_tree();
else
fill_relation_pairs_for_section_biggest_cell();
}
void fill_relation_pairs_for_sector() {
SASSERT(!m_I[m_level].section);
if (m_rel_mode == spanning_tree)
fill_relation_with_spanning_tree_heuristic();
else
fill_relation_with_biggest_cell_heuristic();
}
// Extract roots of polynomial p around sample value v at the given level,
// partitioning them into lhalf (roots <= v) and uhalf (roots > v).
// ps_idx is the index of p in m_level_ps.
// Returns whether the polynomial has any roots.
bool isolate_roots_around_sample(poly* p, unsigned ps_idx,
anum const& v, std_vector<root_function>& lhalf,
std_vector<root_function>& uhalf) {
anum const& v, std_vector<root_function>& lhalf,
std_vector<root_function>& uhalf) {
scoped_anum_vector roots(m_am);
@ -910,14 +821,15 @@ namespace nlsat {
}
bool collect_partitioned_root_functions_around_sample(anum const& v,
std_vector<root_function>& lhalf, std_vector<root_function>& uhalf) {
std_vector<root_function>& lhalf, std_vector<root_function>& uhalf) {
for (unsigned i = 0; i < m_level_ps.size(); ++i)
m_poly_has_roots[i] = isolate_roots_around_sample(m_level_ps.get(i), i, v, lhalf, uhalf);
return !lhalf.empty() || !uhalf.empty();
}
std_vector<root_function>::iterator set_relation_root_functions_from_partitions(std_vector<root_function>& lhalf,
std_vector<root_function>& uhalf) {
std_vector<root_function>::iterator set_relation_root_functions_from_partitions(
std_vector<root_function>& lhalf,
std_vector<root_function>& uhalf) {
auto& rfs = m_rel.m_rfunc;
size_t mid_pos = lhalf.size();
rfs.clear();
@ -945,9 +857,9 @@ namespace nlsat {
void sort_root_function_partitions(std_vector<root_function>::iterator mid) {
auto& rfs = m_rel.m_rfunc;
std::sort(rfs.begin(), mid,
[&](root_function const& a, root_function const& b) { return root_function_lt(a, b, true); });
[&](root_function const& a, root_function const& b) { return root_function_lt(a, b, true); });
std::sort(mid, rfs.end(),
[&](root_function const& a, root_function const& b) { return root_function_lt(a, b, false); });
[&](root_function const& a, root_function const& b) { return root_function_lt(a, b, false); });
}
// Populate Θ (root functions) around the sample, partitioned at `mid`, and sort each partition.
@ -981,7 +893,7 @@ namespace nlsat {
m_I[m_level].section = true;
m_I[m_level].l = r.ire.p;
m_I[m_level].l_index = r.ire.i;
m_u_rf = UINT_MAX;
m_u_rf = m_l_rf;
}
else {
m_I[m_level].l = r.ire.p;
@ -1018,35 +930,34 @@ namespace nlsat {
return false;
set_interval_from_root_partition(v, mid);
compute_side_mask();
return true;
}
void add_relation_resultants() {
TRACE(lws,
tout << " add_relation_resultants: " << m_rel.m_pairs.size() << " pairs\n";
);
TRACE(lws, tout << " add_relation_resultants: " << m_rel.m_pairs.size() << " pairs\n";);
for (auto const& pr : m_rel.m_pairs) {
poly* p1 = m_level_ps.get(pr.first);
poly* p2 = m_level_ps.get(pr.second);
TRACE(lws,
tout << " resultant(" << pr.first << ", " << pr.second << "): ";
m_pm.display(tout, p1);
tout << " and ";
m_pm.display(tout, p2);
tout << "\n";
);
tout << " resultant(" << pr.first << ", " << pr.second << "): ";
m_pm.display(tout, p1);
tout << " and ";
m_pm.display(tout, p2);
tout << "\n";
);
polynomial_ref res = psc_resultant(p1, p2, m_level);
TRACE(lws,
tout << " resultant poly: ";
if (res) {
m_pm.display(tout, res);
tout << "\n resultant sign at sample: " << m_am.eval_sign_at(res, sample());
} else {
tout << "(null)";
}
tout << "\n";
);
tout << " resultant poly: ";
if (res) {
m_pm.display(tout, res);
tout << "\n resultant sign at sample: " << m_am.eval_sign_at(res, sample());
} else {
tout << "(null)";
}
tout << "\n";
);
request_factorized(res, inv_req::ord);
}
}
@ -1063,6 +974,14 @@ namespace nlsat {
scoped_anum_vector roots(m_am);
m_am.isolate_roots(polynomial_ref(p, m_pm), undef_var_assignment(sample(), m_n), roots);
m_poly_has_roots[i] = !roots.empty();
TRACE(lws,
tout << " poly[" << i << "] has " << roots.size() << " roots: ";
for (unsigned k = 0; k < roots.size(); ++k) {
if (k > 0) tout << ", ";
m_am.display_decimal(tout, roots[k], 5);
}
tout << "\n";
);
for (unsigned k = 0; k < roots.size(); ++k) {
scoped_anum root_v(m_am);
m_am.set(root_v, roots[k]);
@ -1075,6 +994,18 @@ namespace nlsat {
std::sort(root_vals.begin(), root_vals.end(), [&](auto const& a, auto const& b) {
return m_am.lt(a.first, b.first);
});
TRACE(lws,
tout << " Sorted roots:\n";
for (unsigned j = 0; j < root_vals.size(); ++j) {
tout << " [" << j << "] val=";
m_am.display_decimal(tout, root_vals[j].first, 5);
tout << " poly=";
m_pm.display(tout, root_vals[j].second);
tout << "\n";
}
);
std::set<std::pair<poly*, poly*>> added_pairs;
for (unsigned j = 0; j + 1 < root_vals.size(); ++j) {
poly* p1 = root_vals[j].second;
@ -1084,6 +1015,13 @@ namespace nlsat {
if (p1 > p2) std::swap(p1, p2);
if (!added_pairs.insert({p1, p2}).second)
continue;
TRACE(lws,
tout << " Adjacent resultant pair: ";
m_pm.display(tout, p1);
tout << " and ";
m_pm.display(tout, p2);
tout << "\n";
);
request_factorized(psc_resultant(p1, p2, m_n), inv_req::ord);
}
}
@ -1100,19 +1038,19 @@ namespace nlsat {
void add_level_projections_sector() {
TRACE(lws,
tout << "\n add_level_projections_sector at level " << m_level << "\n";
tout << " Lower bound rf=" << m_l_rf << ", Upper bound rf=" << m_u_rf << "\n";
if (is_set(m_l_rf)) {
tout << " lower poly idx=" << m_rel.m_rfunc[m_l_rf].ps_idx << ": ";
m_pm.display(tout, m_level_ps.get(m_rel.m_rfunc[m_l_rf].ps_idx));
tout << "\n";
}
if (is_set(m_u_rf)) {
tout << " upper poly idx=" << m_rel.m_rfunc[m_u_rf].ps_idx << ": ";
m_pm.display(tout, m_level_ps.get(m_rel.m_rfunc[m_u_rf].ps_idx));
tout << "\n";
}
);
tout << "\n add_level_projections_sector at level " << m_level << "\n";
tout << " Lower bound rf=" << m_l_rf << ", Upper bound rf=" << m_u_rf << "\n";
if (is_set(m_l_rf)) {
tout << " lower poly idx=" << m_rel.m_rfunc[m_l_rf].ps_idx << ": ";
m_pm.display(tout, m_level_ps.get(m_rel.m_rfunc[m_l_rf].ps_idx));
tout << "\n";
}
if (is_set(m_u_rf)) {
tout << " upper poly idx=" << m_rel.m_rfunc[m_u_rf].ps_idx << ": ";
m_pm.display(tout, m_level_ps.get(m_rel.m_rfunc[m_u_rf].ps_idx));
tout << "\n";
}
);
// Lines 11-12 (Algorithm 1): add projections for each p
// Note: Algorithm 1 adds disc + ldcf for ALL polynomials (classical delineability)
// We additionally omit leading coefficients for rootless polynomials when possible
@ -1135,12 +1073,12 @@ namespace nlsat {
bool is_upper_bound = is_set(m_u_rf) && i == m_rel.m_rfunc[m_u_rf].ps_idx;
TRACE(lws,
tout << " poly[" << i << "] is_lower=" << is_lower_bound << " is_upper=" << is_upper_bound;
tout << " omit_lc[i]=" << (i < m_omit_lc.size() ? (m_omit_lc[i] ? "true" : "false") : "N/A");
tout << " add_lc=" << add_lc << "\n";
);
tout << " poly[" << i << "] is_lower=" << is_lower_bound << " is_upper=" << is_upper_bound;
tout << " omit_lc[i]=" << (i < m_omit_lc.size() ? (m_omit_lc[i] ? "true" : "false") : "N/A");
tout << " add_lc=" << add_lc << "\n";
);
if (add_lc && i < usize(m_poly_has_roots) && !m_poly_has_roots[i])
if (add_lc && !poly_has_roots(i))
if (lc && !is_zero(lc) && m_am.eval_sign_at(lc, sample()) != 0)
add_lc = false;
@ -1159,22 +1097,21 @@ namespace nlsat {
if (add_disc && get_req(i) != inv_req::ord && i < m_omit_disc.size() && m_omit_disc[i])
add_disc = false;
add_projections_for(p, m_level, witness, add_lc, add_disc);
add_projection_for_poly(p, m_level, witness, add_lc, add_disc);
}
}
bool is_section() { return m_I[m_level].is_section();}
void add_level_projections_section() {
SASSERT(m_I[m_level].section);
poly* section_p = m_I[m_level].l.get();
SASSERT(is_section());
compute_omit_lc_section();
// m_omit_disc is computed by compute_omit_disc_for_spanning_tree() in process_level_section()
poly* section_p = m_I[m_level].l.get();
SASSERT(section_p);
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
polynomial_ref p(m_level_ps.get(i), m_pm);
bool is_section_poly = section_p && p.get() == section_p;
bool has_roots = i < usize(m_poly_has_roots) && m_poly_has_roots[i];
bool is_section_poly = p.get() == section_p;
polynomial_ref lc(m_pm);
unsigned deg = m_pm.degree(p, m_level);
lc = m_pm.coeff(p, m_level, deg);
@ -1184,7 +1121,7 @@ namespace nlsat {
if (!is_section_poly && m_rel_mode == spanning_tree) {
if (i < m_omit_lc.size() && m_omit_lc[i])
add_lc = false;
if (add_lc && !has_roots)
if (add_lc && !poly_has_roots(i))
if (lc && !is_zero(lc) && m_am.eval_sign_at(lc, sample()) != 0)
add_lc = false;
}
@ -1205,7 +1142,7 @@ namespace nlsat {
if (lc && !is_zero(lc) && m_am.eval_sign_at(lc, sample()) != 0)
witness = polynomial_ref(m_pm);
add_projections_for(p, m_level, witness, add_lc, add_disc);
add_projection_for_poly(p, m_level, witness, add_lc, add_disc);
}
}
@ -1214,15 +1151,22 @@ namespace nlsat {
// Need at least 2 polynomials for a spanning tree to have edges
if (m_spanning_tree_threshold < 2)
return false;
if (m_rel.m_rfunc.size() < 2)
if (m_rel.m_rfunc.size() < 4) // 2 different bounds + at least 2 same out of bounds
return false;
compute_side_mask();
if (m_side_mask.size() == 0)
return false;
if (!is_set(m_l_rf) || !is_set(m_u_rf))
return false;
if (same_boundary_poly())
return false; //spanning tree does not make sense: all same pairs will be obtained with the biggest cell
unsigned both_count = 0;
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
for (unsigned i = 0; i < m_level_ps.size(); ++i)
if (m_side_mask[i] == 3)
if (++both_count >= m_spanning_tree_threshold)
return true;
}
return false;
}
@ -1239,18 +1183,9 @@ namespace nlsat {
void process_level_section(bool have_interval) {
SASSERT(m_I[m_level].section);
clear_level_state(); // Clear stale state from previous level
m_rel_mode = biggest_cell; // default
if (have_interval) {
// Check spanning tree threshold first
if (should_use_spanning_tree()) {
fill_relation_pairs_for_section_spanning_tree();
compute_omit_lc_both_sides(true);
compute_omit_disc_for_spanning_tree();
m_rel_mode = spanning_tree;
} else {
fill_relation_pairs_for_section();
}
fill_relation_pairs_for_section_biggest_cell();
SASSERT(relation_invariant());
}
upgrade_bounds_to_ord();
@ -1260,18 +1195,17 @@ namespace nlsat {
void process_level_sector(bool have_interval) {
SASSERT(!m_I[m_level].section);
clear_level_state(); // Clear stale state from previous level
m_rel_mode = biggest_cell; // default
if (have_interval) {
// Check spanning tree threshold first
if (should_use_spanning_tree()) {
fill_relation_with_spanning_tree_heuristic();
fill_relation_sector_spanning_tree();
compute_omit_lc_both_sides(true);
compute_omit_disc_for_spanning_tree();
m_rel_mode = spanning_tree;
} else {
fill_relation_pairs_for_sector();
}
} else
fill_relation_sector_biggest_cell();
SASSERT(relation_invariant());
}
upgrade_bounds_to_ord();
@ -1279,17 +1213,7 @@ namespace nlsat {
add_relation_resultants();
}
void process_level() {
TRACE(lws,
tout << "\n--- process_level: level=" << m_level << " ---\n";
tout << "Polynomials at this level (" << m_level_ps.size() << "):\n";
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
tout << " ps[" << i << "]: ";
m_pm.display(tout, m_level_ps.get(i));
tout << "\n";
}
);
void collect_non_null_witnesses() {
// Line 10/11: detect nullification + pick a non-zero coefficient witness per p.
m_witnesses.clear();
m_witnesses.reserve(m_level_ps.size());
@ -1300,22 +1224,33 @@ namespace nlsat {
fail();
m_witnesses.push_back(w);
}
}
void display_polys_at_level(std::ostream& out) {
TRACE(lws, out << "Polynomials at level " << m_level << "\n";
for (unsigned i = 0; i < m_level_ps.size(); ++i)
m_pm.display(out, m_level_ps.get(i)) << "\n";);
}
void process_level() {
TRACE(lws, display_polys_at_level(tout););
clear_level_state(); // Clear stale state from previous level
collect_non_null_witnesses();
// Lines 3-8: Θ + I_level + relation ≼
bool have_interval = build_interval();
TRACE(lws,
tout << "Interval: ";
display(tout, m_solver, m_I[m_level]);
tout << "\n";
tout << "Section? " << (m_I[m_level].section ? "yes" : "no") << "\n";
tout << "have_interval=" << have_interval << ", rfunc.size=" << m_rel.m_rfunc.size() << "\n";
for (unsigned i = 0; i < m_rel.m_rfunc.size(); ++i) {
tout << " rfunc[" << i << "]: ps_idx=" << m_rel.m_rfunc[i].ps_idx << ", val=";
m_am.display(tout, m_rel.m_rfunc[i].val);
tout << "\n";
}
);
tout << "Interval: ";
display(tout, m_solver, m_I[m_level]);
tout << "\n";
tout << "Section? " << (m_I[m_level].section ? "yes" : "no") << "\n";
tout << "have_interval=" << have_interval << ", rfunc.size=" << m_rel.m_rfunc.size() << "\n";
for (unsigned i = 0; i < m_rel.m_rfunc.size(); ++i) {
tout << " rfunc[" << i << "]: ps_idx=" << m_rel.m_rfunc[i].ps_idx << ", val=";
m_am.display(tout, m_rel.m_rfunc[i].val);
tout << "\n";
}
);
if (m_I[m_level].section)
process_level_section(have_interval);
@ -1323,74 +1258,50 @@ namespace nlsat {
process_level_sector(have_interval);
}
bool poly_has_roots(unsigned i) { return i < usize(m_poly_has_roots) && m_poly_has_roots[i]; }
void process_top_level() {
TRACE(lws,
tout << "\n--- process_top_level: level=" << m_n << " ---\n";
tout << "Polynomials at top level (" << m_level_ps.size() << "):\n";
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
tout << " ps[" << i << "]: ";
m_pm.display(tout, m_level_ps.get(i));
tout << "\n";
}
);
m_witnesses.clear();
m_witnesses.reserve(m_level_ps.size());
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
polynomial_ref p(m_level_ps.get(i), m_pm);
polynomial_ref w = choose_nonzero_coeff(p, m_n);
if (!w)
fail();
m_witnesses.push_back(w);
}
// Resultants between adjacent root functions (a lightweight ordering for the top level).
TRACE(lws, display_polys_at_level(tout););
collect_non_null_witnesses();
add_adjacent_root_resultants();
// Projections (coeff witness, disc, leading coeff).
// Note: SMT-RAT's levelwise implementation additionally has dedicated heuristics for
// selecting resultants and selectively adding leading coefficients (see OneCellCAD.h,
// sectionHeuristic/sectorHeuristic). Z3's levelwise currently uses the relation_mode
// ordering heuristics instead of these specialized cases.
for (unsigned i = 0; i < m_level_ps.size(); ++i) {
polynomial_ref p(m_level_ps.get(i), m_pm);
polynomial_ref lc(m_pm);
unsigned deg = m_pm.degree(p, m_n);
lc = m_pm.coeff(p, m_n, deg);
bool add_lc = true; // Let todo_set handle duplicates if witness == lc
if (i < usize(m_poly_has_roots) && !m_poly_has_roots[i]) {
bool add_lc = true;
if (!poly_has_roots(i))
if (lc && !is_zero(lc) && m_am.eval_sign_at(lc, sample()) != 0)
add_lc = false;
}
// SMT-RAT-style coeffNonNull: if the leading coefficient is already non-zero at the sample
// if the leading coefficient is already non-zero at the sample
// AND we're adding lc, we do not need to project an additional non-null coefficient witness.
polynomial_ref witness = m_witnesses[i];
if (add_lc && witness && !is_const(witness) && add_lc)
if (add_lc && witness && !is_const(witness))
if (lc && !is_zero(lc) && m_am.eval_sign_at(lc, sample()) != 0)
witness = polynomial_ref(m_pm);
add_projections_for(p, m_n, witness, add_lc, true);
witness = polynomial_ref(m_pm); // zero the witnsee as lc will be the witness
add_projection_for_poly(p, m_n, witness, add_lc, true); //true to add the discriminant
}
}
std_vector<root_function_interval> single_cell_work() {
TRACE(lws,
tout << "Input polynomials (" << m_P.size() << "):\n";
for (unsigned i = 0; i < m_P.size(); ++i) {
tout << " p[" << i << "]: ";
m_pm.display(tout, m_P.get(i));
tout << "\n";
}
tout << "Sample values:\n";
for (unsigned j = 0; j < m_n; ++j) {
tout << " x" << j << " = ";
m_am.display(tout, sample().value(j));
tout << "\n";
}
);
tout << "Input polynomials (" << m_P.size() << "):\n";
for (unsigned i = 0; i < m_P.size(); ++i) {
tout << " p[" << i << "]: ";
m_pm.display(tout, m_P.get(i)) << "\n";
}
tout << "Sample values:\n";
for (unsigned j = 0; j < m_n; ++j) {
tout << " x" << j << " = ";
m_am.display(tout, sample().value(j)) << "\n";
}
);
if (m_n == 0)
return m_I;
if (m_n == 0) return m_I;
m_todo.reset();
@ -1402,8 +1313,7 @@ namespace nlsat {
});
}
if (m_todo.empty())
return m_I;
if (m_todo.empty()) return m_I;
// Process top level m_n (we project from x_{m_n} and do not return an interval for it).
if (m_todo.max_var() == m_n) {
@ -1417,19 +1327,19 @@ namespace nlsat {
SASSERT(m_level < m_n);
process_level();
TRACE(lws,
tout << "After level " << m_level << ": m_todo.empty()=" << m_todo.empty();
if (!m_todo.empty()) tout << ", m_todo.max_var()=" << m_todo.max_var();
tout << "\n";
);
tout << "After level " << m_level << ": m_todo.empty()=" << m_todo.empty();
if (!m_todo.empty()) tout << ", m_todo.max_var()=" << m_todo.max_var();
tout << "\n";
);
}
TRACE(lws,
for (unsigned i = 0; i < m_I.size(); ++i) {
tout << "I[" << i << "]: ";
display(tout, m_solver, m_I[i]);
tout << "\n";
}
);
for (unsigned i = 0; i < m_I.size(); ++i) {
tout << "I[" << i << "]: ";
display(tout, m_solver, m_I[i]);
tout << "\n";
}
);
return m_I;
}

20
src/nlsat/nlsat_explain.cpp
View file

@ -63,6 +63,9 @@ namespace nlsat {
// Lower-stage polynomials encountered during normalization that need to be projected
polynomial_ref_vector m_lower_stage_polys;
// Store last levelwise input for debugging unsound lemmas
polynomial_ref_vector m_last_lws_input_polys;
// temporary fields for storing the result
scoped_literal_vector * m_result = nullptr;
@ -88,6 +91,7 @@ namespace nlsat {
m_core1(s),
m_core2(s),
m_lower_stage_polys(m_pm),
m_last_lws_input_polys(m_pm),
m_evaluator(ev) {
m_simplify_cores = false;
m_full_dimensional = false;
@ -1012,6 +1016,11 @@ namespace nlsat {
*/
bool levelwise_single_cell(polynomial_ref_vector & ps, var max_x, polynomial::cache & cache) {
// Store polynomials for debugging unsound lemmas
m_last_lws_input_polys.reset();
for (unsigned i = 0; i < ps.size(); i++)
m_last_lws_input_polys.push_back(ps.get(i));
levelwise lws(m_solver, ps, max_x, sample(), m_pm, m_am, cache);
auto cell = lws.single_cell();
if (lws.failed())
@ -1062,6 +1071,7 @@ namespace nlsat {
ps.push_back(p);
bool use_all_coeffs = false;
if (m_solver.apply_levelwise()) {
bool levelwise_ok = levelwise_single_cell(ps, max_x, m_cache);
m_solver.record_levelwise_result(levelwise_ok);
@ -1849,6 +1859,16 @@ namespace nlsat {
m_imp->maximize(x, n, ls, val, unbounded);
}
void explain::display_last_lws_input(std::ostream& out) {
out << "=== POLYNOMIALS PASSED TO LEVELWISE ===\n";
for (unsigned i = 0; i < m_imp->m_last_lws_input_polys.size(); i++) {
out << " p[" << i << "]: ";
m_imp->m_pm.display(out, m_imp->m_last_lws_input_polys.get(i));
out << "\n";
}
out << "=== END LEVELWISE INPUT (" << m_imp->m_last_lws_input_polys.size() << " polynomials) ===\n";
}
void explain::test_root_literal(atom::kind k, var y, unsigned i, poly* p, scoped_literal_vector & result) {
m_imp->test_root_literal(k, y, i, p, result);
}

5
src/nlsat/nlsat_explain.h
View file

@ -103,6 +103,11 @@ namespace nlsat {
*/
void maximize(var x, unsigned n, literal const * ls, scoped_anum& val, bool& unbounded);
/**
Print the polynomials that were passed to levelwise in the last call (for debugging).
*/
void display_last_lws_input(std::ostream& out);
/**
Unit test routine.
*/

2
src/nlsat/nlsat_params.pyg
View file

@ -23,6 +23,6 @@ def_module_params('nlsat',
('zero_disc', BOOL, False, "add_zero_assumption to the vanishing discriminant."),
('known_sat_assignment_file_name', STRING, "", "the file name of a known solution: used for debugging only"),
('lws', BOOL, True, "apply levelwise."),
('lws_spt_threshold', UINT, 3, "minimum both-side polynomial count to apply spanning tree optimization; < 2 disables spanning tree"),
('lws_spt_threshold', UINT, 2, "minimum both-side polynomial count to apply spanning tree optimization; < 2 disables spanning tree"),
('canonicalize', BOOL, True, "canonicalize polynomials.")
))

40
src/nlsat/nlsat_solver.cpp
View file

@ -54,7 +54,6 @@ Revision History:
namespace nlsat {
typedef chashtable<ineq_atom*, ineq_atom::hash_proc, ineq_atom::eq_proc> ineq_atom_table;
typedef chashtable<root_atom*, root_atom::hash_proc, root_atom::eq_proc> root_atom_table;
@ -250,6 +249,7 @@ namespace nlsat {
stats m_stats;
std::string m_debug_known_solution_file_name;
bool m_apply_lws;
bool m_last_conflict_used_lws = false; // Track if last conflict explanation used levelwise
unsigned m_lws_spt_threshold = 3;
imp(solver& s, ctx& c):
m_ctx(c),
@ -1022,7 +1022,8 @@ namespace nlsat {
}
// Helper: Setup checker solver and translate atoms/clauses
void setup_checker(imp& checker, scoped_bool_vars& tr, unsigned n, literal const* cls, assumption_set a) {
// Returns false if the lemma cannot be properly translated for checking
bool setup_checker(imp& checker, scoped_bool_vars& tr, unsigned n, literal const* cls, assumption_set a) {
auto pconvert = [&](poly* p) {
return convert(m_pm, p, checker.m_pm);
};
@ -1058,7 +1059,9 @@ namespace nlsat {
}
else {
// root atom cannot be translated due to variable ordering
bv = checker.mk_bool_var();
// Skip lemma check in this case
TRACE(nlsat, tout << "check_lemma: skipping due to untranslatable root atom\n";);
return false;
}
}
else {
@ -1085,12 +1088,20 @@ namespace nlsat {
literal nlit(tr[lit.var()], !lit.sign());
checker.mk_external_clause(1, &nlit, nullptr);
}
return true; // Successfully set up the checker
}
// Helper: Display unsound lemma failure information
void display_unsound_lemma(imp& checker, scoped_bool_vars& tr, unsigned n, literal const* cls) {
verbose_stream() << "\n";
verbose_stream() << "========== UNSOUND LEMMA DETECTED ==========\n";
verbose_stream() << "Levelwise used for this conflict: " << (m_last_conflict_used_lws ? "YES" : "NO") << "\n";
// Print polynomials passed to levelwise
if (m_last_conflict_used_lws) {
m_explain.display_last_lws_input(verbose_stream());
}
verbose_stream() << "The following lemma is NOT implied by the original clauses:\n";
display(verbose_stream() << " Lemma: ", n, cls) << "\n\n";
verbose_stream() << "Reason: Found a satisfying assignment where:\n";
@ -1130,7 +1141,7 @@ namespace nlsat {
// Create a separate reslimit for the checker with 10 second timeout
reslimit checker_rlimit;
cancel_eh<reslimit> eh(checker_rlimit);
scoped_timer timer(1000, &eh);
scoped_timer timer(1000, &eh); // one second
ctx c(checker_rlimit, m_ctx.m_params, m_ctx.m_incremental);
solver solver2(c);
@ -1139,14 +1150,20 @@ namespace nlsat {
checker.m_log_lemmas = false;
checker.m_dump_mathematica = false;
checker.m_inline_vars = false;
checker.m_apply_lws = false; // Disable levelwise for checker to avoid recursive issues
scoped_bool_vars tr(checker);
setup_checker(checker, tr, n, cls, a);
if (!setup_checker(checker, tr, n, cls, a))
return; // Lemma contains untranslatable atoms, skip check
lbool r = checker.check();
if (r == l_undef) // Checker timed out - skip this lemma check
return;
if (r == l_true) {
// Before reporting unsound, dump the lemma to see what we're checking
verbose_stream() << "Dumping lemma that internal checker thinks is not a tautology:\n";
verbose_stream() << "Checker levelwise calls: " << checker.m_stats.m_levelwise_calls << "\n";
log_lemma(verbose_stream(), n, cls, true, "internal-check-fail");
display_unsound_lemma(checker, tr, n, cls);
exit(1);
}
@ -2125,9 +2142,8 @@ namespace nlsat {
collect(assumptions, m_learned);
del_clauses(m_valids);
if (m_check_lemmas)
for (clause* c : m_learned)
check_lemma(c->size(), c->data(), nullptr);
// Note: Don't check learned clauses here - they are the result of resolution
// and may not be tautologies. Conflict lemmas are checked in resolve_lazy_justification.
assumptions.reset();
assumptions.append(result);
@ -2364,6 +2380,8 @@ namespace nlsat {
}
if (m_check_lemmas) {
TRACE(nlsat, tout << "Checking lazy clause with " << m_lazy_clause.size() << " literals:\n";
display(tout, m_lazy_clause.size(), m_lazy_clause.data()) << "\n";);
check_lemma(m_lazy_clause.size(), m_lazy_clause.data(), nullptr);
m_valids.push_back(mk_clause_core(m_lazy_clause.size(), m_lazy_clause.data(), false, nullptr));
}
@ -2609,9 +2627,8 @@ namespace nlsat {
TRACE(nlsat, tout << "new lemma:\n"; display(tout, m_lemma.size(), m_lemma.data()); tout << "\n";
tout << "found_decision: " << found_decision << "\n";);
if (m_check_lemmas) {
check_lemma(m_lemma.size(), m_lemma.data(), m_lemma_assumptions.get());
}
// Note: Don't check m_lemma here - it's the result of resolution
// and may not be a tautology. Conflict lemmas are checked in resolve_lazy_justification.
// if (m_log_lemmas)
// log_lemma(std::cout, m_lemma.size(), m_lemma.data(), false);
@ -4645,6 +4662,7 @@ namespace nlsat {
void solver::record_levelwise_result(bool success) {
m_imp->m_stats.m_levelwise_calls++;
m_imp->m_last_conflict_used_lws = success; // Track for unsound lemma reporting
if (!success) {
m_imp->m_stats.m_levelwise_failures++;
// m_imp->m_apply_lws = false; // is it useful to throttle

281
src/test/nlsat.cpp
View file

@ -1164,7 +1164,288 @@ static void tst_unsound_lws2380() {
run_test(true); // Levelwise projection
}
// Test case for unsound lemma - levelwise produces cell that's too large
// Input: 5 polynomials with max_var=x3, sample x0=-7, x1=-1, x2=1
// Counterexample: x0=-4, x1=-8, x2=5, x3=6
static void tst_unsound_lws_x3() {
std::cout << "=== tst_unsound_lws_x3 ===\n";
params_ref ps;
ps.set_bool("lws", true);
reslimit rlim;
nlsat::solver s(rlim, ps, false);
anum_manager & am = s.am();
nlsat::pmanager & pm = s.pm();
nlsat::assignment sample_as(am);
nlsat::assignment counter_as(am);
polynomial::cache cache(pm);
// Create 4 variables x0, x1, x2, x3
nlsat::var x0 = s.mk_var(false);
nlsat::var x1 = s.mk_var(false);
nlsat::var x2 = s.mk_var(false);
nlsat::var x3 = s.mk_var(false);
polynomial_ref _x0(pm), _x1(pm), _x2(pm), _x3(pm);
_x0 = pm.mk_polynomial(x0);
_x1 = pm.mk_polynomial(x1);
_x2 = pm.mk_polynomial(x2);
_x3 = pm.mk_polynomial(x3);
// p[0]: x3 + x0
polynomial_ref p0(pm);
p0 = _x3 + _x0;
// p[1]: x1
polynomial_ref p1(pm);
p1 = _x1;
// p[2]: 9 x1^2 x3^2 - 6 x0 x1 x2 x3 + 3 x0 x1^2 x3 - 12 x1^2 x3 + x0^2 x2^2 - x0^2 x1 x2 - 4 x0 x1 x2 + 2 x0 x1^2 + 4 x1^2
polynomial_ref p2(pm);
p2 = 9 * (_x1^2) * (_x3^2)
- 6 * _x0 * _x1 * _x2 * _x3
+ 3 * _x0 * (_x1^2) * _x3
- 12 * (_x1^2) * _x3
+ (_x0^2) * (_x2^2)
- (_x0^2) * _x1 * _x2
- 4 * _x0 * _x1 * _x2
+ 2 * _x0 * (_x1^2)
+ 4 * (_x1^2);
// p[3]: 9 x1^2 x3^2 - 6 x0 x1 x2 x3 + 6 x0 x1^2 x3 - 12 x1^2 x3 + x0^2 x2^2 - 2 x0^2 x1 x2 - 4 x0 x1 x2 + x0^2 x1^2 + 4 x0 x1^2 + 4 x1^2
polynomial_ref p3(pm);
p3 = 9 * (_x1^2) * (_x3^2)
- 6 * _x0 * _x1 * _x2 * _x3
+ 6 * _x0 * (_x1^2) * _x3
- 12 * (_x1^2) * _x3
+ (_x0^2) * (_x2^2)
- 2 * (_x0^2) * _x1 * _x2
- 4 * _x0 * _x1 * _x2
+ (_x0^2) * (_x1^2)
+ 4 * _x0 * (_x1^2)
+ 4 * (_x1^2);
// p[4]: x3 + x1 + x0
polynomial_ref p4(pm);
p4 = _x3 + _x1 + _x0;
std::cout << "p0: " << p0 << "\n";
std::cout << "p1: " << p1 << "\n";
std::cout << "p2: " << p2 << "\n";
std::cout << "p3: " << p3 << "\n";
std::cout << "p4: " << p4 << "\n\n";
// Set sample: x0=-7, x1=-1, x2=1, x3=? (need to pick a value in the cell)
// For the sample, we need an x3 value. Let's use x3=8 (which is > -x0 = 7, so p0 > 0)
scoped_anum val(am);
am.set(val, -7); sample_as.set(x0, val);
am.set(val, -1); sample_as.set(x1, val);
am.set(val, 1); sample_as.set(x2, val);
am.set(val, 8); sample_as.set(x3, val);
// Counterexample: x0=-4, x1=-8, x2=5, x3=6
am.set(val, -4); counter_as.set(x0, val);
am.set(val, -8); counter_as.set(x1, val);
am.set(val, 5); counter_as.set(x2, val);
am.set(val, 6); counter_as.set(x3, val);
std::cout << "Sample point: x0=-7, x1=-1, x2=1, x3=8\n";
std::cout << "Counterexample: x0=-4, x1=-8, x2=5, x3=6\n\n";
// Evaluate polynomials at sample
std::cout << "Polynomial signs at SAMPLE:\n";
std::cout << " p0 sign: " << am.eval_sign_at(p0, sample_as) << "\n";
std::cout << " p1 sign: " << am.eval_sign_at(p1, sample_as) << "\n";
std::cout << " p2 sign: " << am.eval_sign_at(p2, sample_as) << "\n";
std::cout << " p3 sign: " << am.eval_sign_at(p3, sample_as) << "\n";
std::cout << " p4 sign: " << am.eval_sign_at(p4, sample_as) << "\n\n";
// Evaluate polynomials at counterexample
std::cout << "Polynomial signs at COUNTEREXAMPLE:\n";
std::cout << " p0 sign: " << am.eval_sign_at(p0, counter_as) << "\n";
std::cout << " p1 sign: " << am.eval_sign_at(p1, counter_as) << "\n";
std::cout << " p2 sign: " << am.eval_sign_at(p2, counter_as) << "\n";
std::cout << " p3 sign: " << am.eval_sign_at(p3, counter_as) << "\n";
std::cout << " p4 sign: " << am.eval_sign_at(p4, counter_as) << "\n\n";
// Set solver assignment for levelwise (without x3)
am.set(val, -7); sample_as.set(x0, val);
am.set(val, -1); sample_as.set(x1, val);
am.set(val, 1); sample_as.set(x2, val);
s.set_rvalues(sample_as);
// Build polynomial vector
polynomial_ref_vector polys(pm);
polys.push_back(p0);
polys.push_back(p1);
polys.push_back(p2);
polys.push_back(p3);
polys.push_back(p4);
unsigned max_x = x3;
// Print roots of each polynomial at sample
std::cout << "Roots of polynomials at sample (in x3):\n";
for (unsigned i = 0; i < polys.size(); ++i) {
polynomial_ref p(polys.get(i), pm);
if (pm.max_var(p) != x3) {
std::cout << " p" << i << ": max_var is not x3, skipping\n";
continue;
}
scoped_anum_vector roots(am);
am.isolate_roots(p, nlsat::undef_var_assignment(sample_as, x3), roots);
std::cout << " p" << i << " roots: ";
if (roots.empty()) {
std::cout << "(none)";
} else {
for (unsigned j = 0; j < roots.size(); ++j) {
if (j > 0) std::cout << ", ";
am.display_decimal(std::cout, roots[j], 5);
}
}
std::cout << "\n";
}
std::cout << "\n";
// Compute and evaluate resultant of p3 and p4
std::cout << "Resultant of p3 and p4 (in x3):\n";
polynomial_ref res_p3_p4(pm);
{
pm.resultant(p3, p4, x3, res_p3_p4);
std::cout << " Res(p3, p4) = ";
pm.display(std::cout, res_p3_p4);
std::cout << "\n";
std::cout << " Sign at sample (x0=-7, x1=-1, x2=1): " << am.eval_sign_at(res_p3_p4, sample_as) << "\n";
std::cout << " Sign at counter (x0=-4, x1=-8, x2=5): " << am.eval_sign_at(res_p3_p4, counter_as) << "\n";
// Check roots of the resultant at x2 level (parametric in x0, x1)
std::cout << " Roots at sample x0,x1 (in x2): ";
scoped_anum_vector res_roots(am);
nlsat::assignment partial_sample(am);
scoped_anum val(am);
am.set(val, -7); partial_sample.set(x0, val);
am.set(val, -1); partial_sample.set(x1, val);
am.isolate_roots(res_p3_p4, nlsat::undef_var_assignment(partial_sample, x2), res_roots);
for (unsigned j = 0; j < res_roots.size(); ++j) {
if (j > 0) std::cout << ", ";
am.display_decimal(std::cout, res_roots[j], 5);
}
std::cout << "\n";
// Check roots at counterexample x0,x1
std::cout << " Roots at counter x0,x1 (in x2): ";
nlsat::assignment partial_counter(am);
am.set(val, -4); partial_counter.set(x0, val);
am.set(val, -8); partial_counter.set(x1, val);
scoped_anum_vector res_roots_counter(am);
am.isolate_roots(res_p3_p4, nlsat::undef_var_assignment(partial_counter, x2), res_roots_counter);
for (unsigned j = 0; j < res_roots_counter.size(); ++j) {
if (j > 0) std::cout << ", ";
am.display_decimal(std::cout, res_roots_counter[j], 5);
}
std::cout << "\n";
// Compute and check discriminant of Res(p3,p4) in x2
std::cout << "\n Discriminant of Res(p3,p4) in x2:\n";
polynomial_ref disc_res(pm);
pm.discriminant(res_p3_p4, x2, disc_res);
std::cout << " Disc = ";
pm.display(std::cout, disc_res);
std::cout << "\n";
std::cout << " Sign at sample (x0=-7, x1=-1): " << am.eval_sign_at(disc_res, sample_as) << "\n";
std::cout << " Sign at counter (x0=-4, x1=-8): " << am.eval_sign_at(disc_res, counter_as) << "\n";
}
std::cout << "\n";
std::cout << "Running levelwise with max_x = x3\n";
// Run levelwise
nlsat::levelwise lws(s, polys, max_x, s.sample(), pm, am, cache);
auto cell = lws.single_cell();
std::cout << "Levelwise " << (lws.failed() ? "FAILED" : "succeeded") << "\n";
std::cout << "Cell intervals (count=" << cell.size() << "):\n";
for (auto const& interval : cell) {
nlsat::display(std::cout << " ", s, interval) << "\n";
}
// Evaluate cell bounds at counterexample to check if counterexample is in cell
std::cout << "\n--- Checking if counterexample is in cell ---\n";
std::cout << "For a SECTOR (lower_root, upper_root), variable x satisfies:\n";
std::cout << " x > lower_root AND x < upper_root\n\n";
// For univariate evaluation, we need to substitute lower vars
// Level 0: x0 interval, evaluate at x0=-4
// Level 1: x1 interval (parametric in x0), evaluate at (x0=-4, x1=-8)
// Level 2: x2 interval (parametric in x0,x1), evaluate at (x0=-4,x1=-8,x2=5)
bool counterexample_outside_cell = false;
for (unsigned i = 0; i < cell.size(); ++i) {
auto const& interval = cell[i];
nlsat::var level = i;
std::cout << "Level " << level << ":\n";
// Build assignment up to this level (exclusive) for root isolation
nlsat::assignment partial_as(am);
scoped_anum val(am);
if (level > 0) { am.set(val, -4); partial_as.set(x0, val); }
if (level > 1) { am.set(val, -8); partial_as.set(x1, val); }
if (level > 2) { am.set(val, 5); partial_as.set(x2, val); }
scoped_anum counter_val(am);
if (level == 0) am.set(counter_val, -4);
else if (level == 1) am.set(counter_val, -8);
else if (level == 2) am.set(counter_val, 5);
if (interval.is_section()) {
std::cout << " Section case\n";
} else {
// Isolate roots and check bounds
if (!interval.l_inf()) {
polynomial_ref lower_p(interval.l, pm);
scoped_anum_vector lower_roots(am);
am.isolate_roots(lower_p, nlsat::undef_var_assignment(partial_as, level), lower_roots);
if (lower_roots.size() >= interval.l_index) {
std::cout << " Lower root (root[" << interval.l_index << "]): ";
am.display_decimal(std::cout, lower_roots[interval.l_index - 1], 10);
std::cout << "\n";
std::cout << " Counter x" << level << " = ";
am.display_decimal(std::cout, counter_val, 10);
int cmp = am.compare(counter_val, lower_roots[interval.l_index - 1]);
std::cout << " -> " << (cmp > 0 ? "ABOVE" : (cmp < 0 ? "BELOW" : "EQUAL")) << " lower bound\n";
if (cmp <= 0) counterexample_outside_cell = true;
}
}
if (!interval.u_inf()) {
polynomial_ref upper_p(interval.u, pm);
scoped_anum_vector upper_roots(am);
am.isolate_roots(upper_p, nlsat::undef_var_assignment(partial_as, level), upper_roots);
if (upper_roots.size() >= interval.u_index) {
std::cout << " Upper root (root[" << interval.u_index << "]): ";
am.display_decimal(std::cout, upper_roots[interval.u_index - 1], 10);
std::cout << "\n";
std::cout << " Counter x" << level << " = ";
am.display_decimal(std::cout, counter_val, 10);
int cmp = am.compare(counter_val, upper_roots[interval.u_index - 1]);
std::cout << " -> " << (cmp > 0 ? "ABOVE" : (cmp < 0 ? "BELOW" : "EQUAL")) << " upper bound\n";
if (cmp >= 0) counterexample_outside_cell = true;
}
}
}
std::cout << "\n";
}
// The counterexample has different polynomial signs than the sample.
// For a sound cell, the counterexample must be OUTSIDE the cell.
ENSURE(counterexample_outside_cell);
std::cout << "SUCCESS: Counterexample is OUTSIDE the cell (cell is sound)\n";
std::cout << "=== END tst_unsound_lws_x3 ===\n\n";
}
void tst_nlsat() {
tst_unsound_lws_x3();
std::cout << "------------------\n";
tst_unsound_lws2380();
std::cout << "------------------\n";
tst_polynomial_cache_mk_unique();