fix #4116

delta has to be computed based on Simplex tableau not on difference graph.

Author: NikolajBjorner

src/math/simplex/simplex.cpp

@@ -20,7 +20,50 @@ Module Name:
 #include "math/simplex/simplex.h"
 #include "math/simplex/sparse_matrix_def.h"
 #include "math/simplex/simplex_def.h"
+#include "util/rational.h"
+#include "util/inf_rational.h"
+
 namespace simplex {
     template class simplex<mpz_ext>;
     template class simplex<mpq_ext>;
+
+    static void refine_delta(rational& delta, inf_rational const& l, inf_rational const& u) {
+        if (l.get_rational() < u.get_rational() && l.get_infinitesimal() > u.get_infinitesimal()) {
+            rational new_delta = (u.get_rational() - l.get_rational()) / (l.get_infinitesimal() - u.get_infinitesimal());
+            if (new_delta < delta) { 
+                delta = new_delta;
+            }
+        }
+    }
+
+
+    void ensure_rational_solution(simplex<mpq_ext>& S) {
+        rational delta(1);
+        for (unsigned i = 0; i < S.get_num_vars(); ++i) {
+            auto const& _value = S.get_value(i);
+            inf_rational value(rational(_value.first), rational(_value.second));
+            if (S.lower_valid(i)) {
+                auto const& _bound = S.get_lower(i);
+                inf_rational bound(rational(_bound.first), rational(_bound.second));
+                refine_delta(delta, bound, value);
+            }
+            if (S.upper_valid(i)) {
+                auto const& _bound = S.get_upper(i);
+                inf_rational bound(rational(_bound.first), rational(_bound.second));
+                refine_delta(delta, value, bound);
+            }
+        }
+        unsynch_mpq_inf_manager inf_mgr;
+        scoped_mpq_inf q(inf_mgr);
+        for (unsigned i = 0; i < S.get_num_vars(); ++i) {
+            auto const& _value = S.get_value(i);
+            rational inf(_value.second);
+            if (!inf.is_zero()) {
+                rational fin = rational(_value.first) + inf * delta;
+                inf = 0;
+                inf_mgr.set(q, fin.to_mpq(), inf.to_mpq());
+                S.set_value(i, q);
+            }
+        }
+    }
 };

src/math/simplex/simplex.h

@@ -133,6 +133,8 @@ namespace simplex {
         void  set_upper(var_t var, eps_numeral const& b);
         void  get_lower(var_t var, scoped_eps_numeral& b) const { b = m_vars[var].m_lower; }
         void  get_upper(var_t var, scoped_eps_numeral& b) const { b = m_vars[var].m_upper; }
+        eps_numeral const& get_lower(var_t var) const { return m_vars[var].m_lower; }
+        eps_numeral const& get_upper(var_t var) const { return m_vars[var].m_upper; }
         bool  above_lower(var_t var, eps_numeral const& b) const;
         bool  below_upper(var_t var, eps_numeral const& b) const;
         bool  below_lower(var_t v) const;
@@ -198,6 +200,7 @@ namespace simplex {
         bool is_feasible() const;
     };
 
+    void ensure_rational_solution(simplex<mpq_ext>& s);
 };
 
 #endif

src/math/simplex/simplex_def.h

@@ -5,8 +5,6 @@ Module Name:
 
     simplex_def.h
 
-Abstract:
-
 Author:
 
     Nikolaj Bjorner (nbjorner) 2014-01-15

src/smt/diff_logic.h

@@ -414,6 +414,7 @@ class dl_graph {
                         }
                         break;
                     case DL_PROCESSED:
+                        TRACE("arith", display_edge(tout << "processed twice: ", e););
                         // if two edges with the same source/target occur in the graph.
                         break;
                     default:
@@ -514,10 +515,8 @@ class dl_graph {
     // The method assumes the graph is feasible before the invocation.
     bool enable_edge(edge_id id) {
         edge& e = m_edges[id];
+        SASSERT(is_feasible());
         bool r = true;
-        if (!is_feasible()) {
-            return false;
-        }
         if (!e.is_enabled()) {
             e.enable(m_timestamp);
             m_last_enabled_edge = id;

src/smt/theory_dense_diff_logic_def.h

@@ -985,6 +985,8 @@ namespace smt {
         TRACE("opt", S.display(tout); );    
         SASSERT(is_sat != l_false);
         lbool is_fin = S.minimize(w);
+
+        ensure_rational_solution(S);
 
         switch (is_fin) {
         case l_true: {
@@ -1010,10 +1012,10 @@ namespace smt {
             }
             for (unsigned i = 0; i < num_nodes; ++i) {
                 mpq_inf const& val = S.get_value(i);
-                rational q(val.first), eps(val.second);
+                rational q(val.first);
+                SASSERT(rational(val.second).is_zero());
                 numeral  a(q);
                 m_assignment[i] = a;
-                // TBD: if epsilon is != 0, then adjust a by some small fraction.
             }
             inf_eps result(rational(0), r);
             blocker = mk_gt(v, result);

src/smt/theory_diff_logic_def.h

@@ -341,7 +341,9 @@ void theory_diff_logic<Ext>::pop_scope_eh(unsigned num_scopes) {
     m_scopes.shrink(new_lvl);
     unsigned num_edges = m_graph.get_num_edges();
     m_graph.pop(num_scopes);
-    if (num_edges != m_graph.get_num_edges() && m_num_simplex_edges > 0) {
+    TRACE("arith", m_graph.display(tout););
+    SASSERT(m_graph.is_feasible());
+    if (true || (num_edges != m_graph.get_num_edges() && m_num_simplex_edges > 0)) {
         m_S.reset();
         m_num_simplex_edges = 0;
         m_objective_rows.reset();
@@ -565,6 +567,7 @@ bool theory_diff_logic<Ext>::propagate_atom(atom* a) {
     }
     int edge_id = a->get_asserted_edge();
     if (!m_graph.enable_edge(edge_id)) {
+        TRACE("arith", display(tout););
         set_neg_cycle_conflict();
         return false;
     }
@@ -685,11 +688,8 @@ void theory_diff_logic<Ext>::set_neg_cycle_conflict() {
     context & ctx = get_context();
     TRACE("arith_conflict", 
           tout << "conflict: ";
-          for (unsigned i = 0; i < lits.size(); ++i) {
-              ctx.display_literal_info(tout, lits[i]);
-          }
-          tout << "\n";
-          );
+          for (literal lit : lits) ctx.display_literal_info(tout, lit);          
+          tout << "\n";);
 
     if (dump_lemmas()) {
         symbol logic(m_lia_or_lra == is_lia ? "QF_LIA" : "QF_LRA");
@@ -911,6 +911,7 @@ void theory_diff_logic<Ext>::compute_delta() {
 }
 
 
+
 template<typename Ext>
 void theory_diff_logic<Ext>::init_model(smt::model_generator & m) {
     m_factory = alloc(arith_factory, get_manager());
@@ -935,9 +936,11 @@ model_value_proc * theory_diff_logic<Ext>::mk_value(enode * n, model_generator &
 
 template<typename Ext>
 void theory_diff_logic<Ext>::display(std::ostream & out) const {
-    for (unsigned i = 0; i < m_atoms.size(); ++i) {
-        m_atoms[i]->display(*this, out);
+    out << "atoms\n";
+    for (atom* a : m_atoms) {
+        a->display(*this, out) << "\n";
     }
+    out << "graph\n";
     m_graph.display(out);
 }
 
@@ -1150,6 +1153,8 @@ void theory_diff_logic<Ext>::update_simplex(Simplex& S) {
     }
     S.set_lower(node2simplex(get_zero(true)), mpq_inf(mpq(0), mpq(0)));
     S.set_upper(node2simplex(get_zero(true)), mpq_inf(mpq(0), mpq(0)));
+    S.set_lower(node2simplex(get_zero(false)), mpq_inf(mpq(0), mpq(0)));
+    S.set_upper(node2simplex(get_zero(false)), mpq_inf(mpq(0), mpq(0)));
     svector<unsigned> vars;
     scoped_mpq_vector coeffs(mgr);
     coeffs.push_back(mpq(1));
@@ -1215,6 +1220,8 @@ typename theory_diff_logic<Ext>::inf_eps theory_diff_logic<Ext>::value(theory_va
      return r;
 }
 
+
+
 template<typename Ext>
 typename theory_diff_logic<Ext>::inf_eps 
 theory_diff_logic<Ext>::maximize(theory_var v, expr_ref& blocker, bool& has_shared) {
@@ -1224,6 +1231,9 @@ theory_diff_logic<Ext>::maximize(theory_var v, expr_ref& blocker, bool& has_shar
     Simplex& S = m_S;
     ast_manager& m = get_manager();
 
+    CTRACE("arith",!m_graph.is_feasible(), m_graph.display(tout););
+    SASSERT(m_graph.is_feasible());
+
     update_simplex(S);
 
     TRACE("arith",
@@ -1235,7 +1245,12 @@ theory_diff_logic<Ext>::maximize(theory_var v, expr_ref& blocker, bool& has_shar
           tout << "Free coefficient " << m_objective_consts[v] << "\n";
           );
 
-    TRACE("opt", S.display(tout); display(tout););
+    TRACE("opt", 
+          S.display(tout); 
+          for (unsigned i = 0; i < m_graph.get_num_nodes(); ++i)
+              tout << "$" << i << ": " << node2simplex(i) << "\n";
+          display(tout);
+          );
 
     // optimize    
     lbool is_sat = S.make_feasible();
@@ -1252,8 +1267,6 @@ theory_diff_logic<Ext>::maximize(theory_var v, expr_ref& blocker, bool& has_shar
         simplex::mpq_ext::eps_numeral const& val = S.get_value(w);
         inf_rational r(-rational(val.first), -rational(val.second));
         Simplex::row row = m_objective_rows[v];
-        TRACE("opt", tout << r << " " << "\n"; 
-              S.display_row(tout, row, true););
         Simplex::row_iterator it = S.row_begin(row), end = S.row_end(row);
         expr_ref_vector& core = m_objective_assignments[v];
         expr_ref tmp(m);
@@ -1269,13 +1282,21 @@ theory_diff_logic<Ext>::maximize(theory_var v, expr_ref& blocker, bool& has_shar
                 }
             }
         }
-        compute_delta();
+        ensure_rational_solution(S);
+        TRACE("opt", tout << r << " " << "\n"; 
+              S.display_row(tout, row, true);
+              S.display(tout);
+              );
+
         for (unsigned i = 0; i < m_graph.get_num_nodes(); ++i) {
             unsigned w = node2simplex(i);
-            simplex::mpq_ext::eps_numeral const& val = S.get_value(w);
-            rational r = rational(val.first) + m_delta*rational(val.second);
+            auto const& val = S.get_value(w);
+            SASSERT(rational(val.second).is_zero());
+            rational r = rational(val.first);
             m_graph.set_assignment(i, numeral(r));
         }
+        CTRACE("arith",!m_graph.is_feasible(), m_graph.display(tout););
+        SASSERT(m_graph.is_feasible());
         inf_eps r1(rational(0), r);
         blocker = mk_gt(v, r1);
         return inf_eps(rational(0), r + m_objective_consts[v]);

src/util/mpq_inf.h

@@ -286,4 +286,8 @@ typedef mpq_inf_manager<false> synch_mpq_inf_manager;
 #endif
 typedef mpq_inf_manager<false> unsynch_mpq_inf_manager;
 
+typedef _scoped_numeral<unsynch_mpq_inf_manager> scoped_mpq_inf;
+typedef _scoped_numeral<synch_mpq_inf_manager> scoped_synch_mpq_inf;
+typedef _scoped_numeral_vector<unsynch_mpq_inf_manager> scoped_mpq_inf_vector;
+
 #endif