Merge branch 'opt' of https://git01.codeplex.com/z3 into opt

src/smt/diff_logic.h

@@ -25,6 +25,7 @@ Revision History:
 #include"trace.h"
 #include"warning.h"
 #include"uint_set.h"
+#include<deque>
 
 typedef int dl_var;
 
@@ -281,6 +282,8 @@ public:
 
     unsigned get_num_edges() const { return m_edges.size(); }
 
+    unsigned get_num_nodes() const { return m_out_edges.size(); }
+
     dl_var get_source(edge_id id) const {  return m_edges[id].get_source(); }
 
     dl_var get_target(edge_id id) const {  return m_edges[id].get_target(); }
@@ -931,6 +934,120 @@ public:
         return found;
     }
 
+    // Return true if there is an edge source --> target (also counting disabled edges).
+    // If there is such edge, return its edge_id in parameter id.
+    bool get_edge_id(dl_var source, dl_var target, edge_id & id) const {
+        edge_id_vector const & edges = m_out_edges[source];
+        typename edge_id_vector::const_iterator it  = edges.begin();
+        typename edge_id_vector::const_iterator end = edges.end();
+        for (; it != end; ++it) {
+            id = *it;
+            edge const & e = m_edges[id];
+            if (e.get_target() == target) {
+                return true;
+            }
+        }
+        return false;
+    }
+
+    edges const & get_all_edges() const {
+        return m_edges;
+    }
+
+    void get_neighbours_undirected(dl_var current, svector<dl_var> & neighbours) {
+	    neighbours.reset();
+	    edge_id_vector & out_edges = m_out_edges[current];
+        typename edge_id_vector::iterator it = out_edges.begin(), end = out_edges.end();
+        for (; it != end; ++it) {
+            edge_id e_id = *it;
+            edge & e     = m_edges[e_id];
+            SASSERT(e.get_source() == current);
+            dl_var neighbour = e.get_target();
+            neighbours.push_back(neighbour);
+        }
+        edge_id_vector & in_edges = m_in_edges[current];
+        typename edge_id_vector::iterator it2 = in_edges.begin(), end2 = in_edges.end();
+	    for (; it2 != end2; ++it2) {
+            edge_id e_id = *it2;
+            edge & e     = m_edges[e_id];
+            SASSERT(e.get_target() == current);
+            dl_var neighbour = e.get_source();
+            neighbours.push_back(neighbour);
+        }
+    }
+
+    void dfs_undirected(dl_var start, svector<dl_var> & threads) {
+        threads.reset();
+        threads.resize(get_num_nodes());
+	    uint_set discovered, explored;
+	    svector<dl_var> nodes;
+        discovered.insert(start);
+	    nodes.push_back(start);
+	    dl_var prev = start;
+	    while(!nodes.empty()) {
+		    dl_var current = nodes.back();
+            SASSERT(discovered.contains(current) && !explored.contains(current));
+		    svector<dl_var> neighbours;
+		    get_neighbours_undirected(current, neighbours);
+            SASSERT(!neighbours.empty());
+            bool found = false;
+		    for (unsigned i = 0; i < neighbours.size(); ++i) {
+                dl_var next = neighbours[i];
+                DEBUG_CODE(
+                edge_id id;
+                SASSERT(get_edge_id(current, next, id) || get_edge_id(next, current, id)););                
+                if (!discovered.contains(next) && !explored.contains(next)) {
+                    TRACE("diff_logic", tout << "thread[" << prev << "] --> " << next << std::endl;);
+                    threads[prev] = next;
+                    prev = next;
+                    discovered.insert(next);
+	                nodes.push_back(next);
+                    found = true;
+                    break;
+                }
+		    }            
+            SASSERT(!nodes.empty());
+            if (!found) {
+                explored.insert(current);
+                nodes.pop_back();
+            }
+	    }
+	    threads[prev] = start;
+    }
+
+    void bfs_undirected(dl_var start, svector<dl_var> & parents, svector<dl_var> & depths) {
+        parents.reset();
+        parents.resize(get_num_nodes());
+        parents[start] = -1;
+        depths.reset();
+        depths.resize(get_num_nodes());
+	    uint_set visited;
+	    std::deque<dl_var> nodes;
+	    visited.insert(start);
+	    nodes.push_front(start);
+	    while(!nodes.empty()) {
+            dl_var current = nodes.back();
+            nodes.pop_back();
+		    SASSERT(visited.contains(current));
+            svector<dl_var> neighbours;
+		    get_neighbours_undirected(current, neighbours);
+            SASSERT(!neighbours.empty());
+		    for (unsigned i = 0; i < neighbours.size(); ++i) {
+			    dl_var next = neighbours[i];
+                DEBUG_CODE(
+                edge_id id;
+                SASSERT(get_edge_id(current, next, id) || get_edge_id(next, current, id)););
+                if (!visited.contains(next)) {
+                    TRACE("diff_logic", tout << "parents[" << next << "] --> " << current << std::endl;);
+	                parents[next] = current;
+	                depths[next] = depths[current] + 1;
+	                visited.insert(next);
+	                nodes.push_front(next);
+                }
+		    }
+	    }
+    }
+
     template<typename Functor>
     void enumerate_edges(dl_var source, dl_var target, Functor& f) {
         edge_id_vector & edges = m_out_edges[source];
@@ -1076,6 +1193,10 @@ public:
         return m_assignment[v]; 
     }
 
+    void set_assignment(dl_var v, numeral const & n) {
+        m_assignment[v] = n; 
+    }
+
     unsigned get_timestamp() const {
         return m_timestamp;
     }

src/smt/elim_term_ite.cpp

@@ -33,16 +33,16 @@ void elim_term_ite::operator()(expr * n,
     proof * pr2;
     get_cached(n, r2, pr2);
     r = r2;
-    switch (m_manager.proof_mode()) {
+    switch (m.proof_mode()) {
     case PGM_DISABLED:
-        pr = m_manager.mk_undef_proof();
+        pr = m.mk_undef_proof();
         break;
     case PGM_COARSE:
         remove_duplicates(m_coarse_proofs);
-        pr = n == r2 ? m_manager.mk_oeq_reflexivity(n) : m_manager.mk_apply_defs(n, r, m_coarse_proofs.size(), m_coarse_proofs.c_ptr());
+        pr = n == r2 ? m.mk_oeq_reflexivity(n) : m.mk_apply_defs(n, r, m_coarse_proofs.size(), m_coarse_proofs.c_ptr());
         break;
     case PGM_FINE:
-        pr = pr2 == 0 ? m_manager.mk_oeq_reflexivity(n) : pr2;
+        pr = pr2 == 0 ? m.mk_oeq_reflexivity(n) : pr2;
         break;
     }
     m_coarse_proofs.reset();
@@ -107,36 +107,36 @@ void elim_term_ite::reduce1_app(app * n) {
     m_args.reset();
     
     func_decl * decl = n->get_decl();
-    proof_ref p1(m_manager);
+    proof_ref p1(m);
     get_args(n, m_args, p1);
-    if (!m_manager.fine_grain_proofs())
+    if (!m.fine_grain_proofs())
         p1 = 0;
 
-    expr_ref r(m_manager);
-    r = m_manager.mk_app(decl, m_args.size(), m_args.c_ptr());
-    if (m_manager.is_term_ite(r)) {
-        expr_ref   new_def(m_manager);
-        proof_ref  new_def_pr(m_manager);
-        app_ref   new_r(m_manager);
-        proof_ref  new_pr(m_manager);
+    expr_ref r(m);
+    r = m.mk_app(decl, m_args.size(), m_args.c_ptr());
+    if (m.is_term_ite(r)) {
+        expr_ref   new_def(m);
+        proof_ref  new_def_pr(m);
+        app_ref   new_r(m);
+        proof_ref  new_pr(m);
         if (m_defined_names.mk_name(r, new_def, new_def_pr, new_r, new_pr)) {
-            CTRACE("elim_term_ite_bug", new_def.get() == 0, tout << mk_ismt2_pp(r, m_manager) << "\n";);
+            CTRACE("elim_term_ite_bug", new_def.get() == 0, tout << mk_ismt2_pp(r, m) << "\n";);
             SASSERT(new_def.get() != 0);
             m_new_defs->push_back(new_def);
-            if (m_manager.fine_grain_proofs()) {
+            if (m.fine_grain_proofs()) {
                 m_new_def_proofs->push_back(new_def_pr);
-                new_pr = m_manager.mk_transitivity(p1, new_pr);
+                new_pr = m.mk_transitivity(p1, new_pr);
             }
             else {
                 // [Leo] This looks fishy... why do we add 0 into m_coarse_proofs when fine_grain_proofs are disabled? 
                 new_pr = 0;
-                if (m_manager.proofs_enabled())
+                if (m.proofs_enabled())
                     m_coarse_proofs.push_back(new_pr);
             }
         }
         else {
             SASSERT(new_def.get() == 0);
-            if (!m_manager.fine_grain_proofs())
+            if (!m.fine_grain_proofs())
                 new_pr = 0;
         }
         cache_result(n, new_r, new_pr);
@@ -151,8 +151,8 @@ void elim_term_ite::reduce1_quantifier(quantifier * q) {
     proof * new_body_pr;
     get_cached(q->get_expr(), new_body, new_body_pr);
     
-    quantifier * new_q = m_manager.update_quantifier(q, new_body);
-    proof *      p     = q == new_q ? 0 : m_manager.mk_oeq_quant_intro(q, new_q, new_body_pr);   
+    quantifier * new_q = m.update_quantifier(q, new_body);
+    proof *      p     = q == new_q ? 0 : m.mk_oeq_quant_intro(q, new_q, new_body_pr);   
     cache_result(q, new_q, p);
 }
 

src/smt/old_interval.cpp

@@ -72,7 +72,11 @@ ext_numeral & ext_numeral::operator-=(ext_numeral const & other) {
 }
 
 ext_numeral & ext_numeral::operator*=(ext_numeral const & other) {
-    if (is_zero() || other.is_zero()) {
+    if (is_zero()) {
+        m_kind = FINITE;
+        return *this;
+    }
+    if (other.is_zero()) {
         m_kind = FINITE;
         m_value.reset();
         return *this;
@@ -295,6 +299,8 @@ interval & interval::operator*=(interval const & other) {
     }
     if (other.is_zero()) {
         *this = other;
+        m_lower_dep = m_manager.mk_join(m_lower_dep, m_upper_dep);
+        m_upper_dep = m_lower_dep;
         return *this;
     }
 

src/smt/params/smt_params.cpp

@@ -49,6 +49,7 @@ void smt_params::updt_params(params_ref const & p) {
     qi_params::updt_params(p);
     theory_arith_params::updt_params(p);
     theory_bv_params::updt_params(p);
+    theory_pb_params::updt_params(p);
     // theory_array_params::updt_params(p);
     updt_local_params(p);
 }

src/smt/params/smt_params.h

@@ -25,6 +25,7 @@ Revision History:
 #include"theory_arith_params.h"
 #include"theory_array_params.h"
 #include"theory_bv_params.h"
+#include"theory_pb_params.h"
 #include"theory_datatype_params.h"
 #include"preprocessor_params.h"
 #include"context_params.h"
@@ -73,7 +74,8 @@ struct smt_params : public preprocessor_params,
                     public qi_params, 
                     public theory_arith_params, 
                     public theory_array_params, 
-                    public theory_bv_params, 
+                    public theory_bv_params,
+                    public theory_pb_params,
                     public theory_datatype_params {
     bool             m_display_proof;
     bool             m_display_dot_proof;

src/smt/params/smt_params_helper.pyg

@@ -44,6 +44,11 @@ def_module_params(module_name='smt',
                           ('arith.branch_cut_ratio', UINT, 2, 'branch/cut ratio for linear integer arithmetic'),
                           ('arith.int_eq_branch', BOOL, False, 'branching using derived integer equations'),
                           ('arith.ignore_int', BOOL, False, 'treat integer variables as real'),
+                          ('arith.dump_lemmas', BOOL, False, 'dump arithmetic theory lemmas to files'),   
+                          ('pb.conflict_frequency', UINT, 1000, 'conflict frequency for Pseudo-Boolean theory'),
+                          ('pb.learn_complements', BOOL, True, 'learn complement literals for Pseudo-Boolean theory'),
+                          ('pb.enable_compilation', BOOL, True, 'enable compilation into sorting circuits for Pseudo-Boolean'),
+                          ('pb.enable_simplex', BOOL, False, 'enable simplex to check rational feasibility'),
                           ('array.weak', BOOL, False, 'weak array theory'),
                           ('array.extensional', BOOL, True, 'extensional array theory')
                           ))

src/smt/params/theory_arith_params.cpp

@@ -34,6 +34,7 @@ void theory_arith_params::updt_params(params_ref const & _p) {
     m_arith_int_eq_branching = p.arith_int_eq_branch();
     m_arith_ignore_int = p.arith_ignore_int();
     m_arith_bound_prop = static_cast<bound_prop_mode>(p.arith_propagation_mode());
+    m_arith_dump_lemmas = p.arith_dump_lemmas();
 }
 
 

src/smt/params/theory_arith_params.h

@@ -27,7 +27,8 @@ enum arith_solver_id {
     AS_DIFF_LOGIC,
     AS_ARITH,
     AS_DENSE_DIFF_LOGIC,
-    AS_UTVPI
+    AS_UTVPI,
+    AS_OPTINF
 };
 
 enum bound_prop_mode {

src/smt/params/theory_pb_params.cpp

@@ -0,0 +1,28 @@
+/*++
+Copyright (c) 2014 Microsoft Corporation
+
+Module Name:
+
+    theory_pb_params.cpp
+
+Abstract:
+
+    <abstract>
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2014-01-01
+
+Revision History:
+
+--*/
+#include"theory_pb_params.h"
+#include"smt_params_helper.hpp"
+
+void theory_pb_params::updt_params(params_ref const & _p) {
+    smt_params_helper p(_p);
+    m_pb_conflict_frequency = p.pb_conflict_frequency();
+    m_pb_learn_complements = p.pb_learn_complements();
+    m_pb_enable_compilation = p.pb_enable_compilation();
+    m_pb_enable_simplex = p.pb_enable_simplex();
+}

src/smt/params/theory_pb_params.h

@@ -0,0 +1,41 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    theory_pb_params.h
+
+Abstract:
+
+    <abstract>
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2014-01-01
+
+Revision History:
+
+--*/
+#ifndef _THEORY_PB_PARAMS_H_
+#define _THEORY_PB_PARAMS_H_
+
+#include"params.h"
+
+
+struct theory_pb_params {
+    unsigned m_pb_conflict_frequency;
+    bool     m_pb_learn_complements;
+    bool     m_pb_enable_compilation;
+    bool     m_pb_enable_simplex;
+    theory_pb_params(params_ref const & p = params_ref()):
+        m_pb_conflict_frequency(1000),
+        m_pb_learn_complements(true),
+        m_pb_enable_compilation(true),
+        m_pb_enable_simplex(false)
+    {}
+    
+    void updt_params(params_ref const & p);
+};
+
+#endif /* _THEORY_PB_PARAMS_H_ */
+

src/smt/smt_b_justification.h

@@ -54,7 +54,7 @@ namespace smt {
             m_data(BOXTAGINT(void*, l.index(), BIN_CLAUSE)) {
         }
         
-        explicit b_justification(justification * js):
+        explicit b_justification(justification * js):            
             m_data(TAG(void*, js, JUSTIFICATION)) {
             SASSERT(js);
         }

src/smt/smt_context.cpp

@@ -199,6 +199,9 @@ namespace smt {
     bool context::bcp() {
         SASSERT(!inconsistent());
         while (m_qhead < m_assigned_literals.size()) {
+            if (m_cancel_flag) {
+                return true;
+            }
             literal l      = m_assigned_literals[m_qhead];
             SASSERT(get_assignment(l) == l_true);
             m_qhead++;
@@ -225,9 +228,6 @@ namespace smt {
                     case l_true:
                         break;
                     }
-                    if (m_cancel_flag) {
-                        return true;
-                    }
                 }
             }
 
@@ -1334,6 +1334,7 @@ namespace smt {
         TRACE("propagate_bool_var_enode_bug", tout << "var: " << v << " #" << bool_var2expr(v)->get_id() << "\n";);
         SASSERT(v < static_cast<int>(m_b_internalized_stack.size()));
         enode * n  = bool_var2enode(v);
+        CTRACE("mk_bool_var", !n, tout << "No enode for " << v << "\n";);
         bool sign  = val == l_false;
         if (n->merge_tf())
             add_eq(n, sign ? m_false_enode : m_true_enode, eq_justification(literal(v, sign)));
@@ -2853,7 +2854,7 @@ namespace smt {
 
     void context::init_assumptions(unsigned num_assumptions, expr * const * assumptions) {
         reset_assumptions();
-        m_bool_var2assumption.reset();
+        m_literal2assumption.reset();
         m_unsat_core.reset();
         if (num_assumptions > 0) {
             // We must give a chance to the theories to propagate before we create a new scope...
@@ -2869,16 +2870,16 @@ namespace smt {
                 proof * pr = m_manager.mk_asserted(curr_assumption);
                 internalize_assertion(curr_assumption, pr, 0);
                 literal l = get_literal(curr_assumption);
-                m_bool_var2assumption.insert(l.var(), curr_assumption);
+                m_literal2assumption.insert(l.index(), curr_assumption);
                 // mark_as_relevant(l); <<< not needed
                 // internalize_assertion marked l as relevant.
                 SASSERT(is_relevant(l));
-                TRACE("assumptions", tout << mk_pp(curr_assumption, m_manager) << "\n";);
+                TRACE("assumptions", tout << l << ":" << mk_pp(curr_assumption, m_manager) << "\n";);
                 if (m_manager.proofs_enabled())
                     assign(l, mk_justification(justification_proof_wrapper(*this, pr)));
                 else
                     assign(l, b_justification::mk_axiom());
-                m_assumptions.push_back(l.var());
+                m_assumptions.push_back(l);
                 get_bdata(l.var()).m_assumption = true;
             }
         }
@@ -2888,10 +2889,10 @@ namespace smt {
     }
 
     void context::reset_assumptions() {
-        bool_var_vector::iterator it  = m_assumptions.begin();
-        bool_var_vector::iterator end = m_assumptions.end();
+        literal_vector::iterator it  = m_assumptions.begin();
+        literal_vector::iterator end = m_assumptions.end();
         for (; it != end; ++it)
-            get_bdata(*it).m_assumption = false;
+            get_bdata(it->var()).m_assumption = false;
         m_assumptions.reset();
     }
 
@@ -2908,10 +2909,9 @@ namespace smt {
             literal l = *it;
             TRACE("unsat_core_bug", tout << "answer literal: " << l << "\n";);
             SASSERT(get_bdata(l.var()).m_assumption);
-            SASSERT(m_bool_var2assumption.contains(l.var()));
-            expr * a = 0;
-            m_bool_var2assumption.find(l.var(), a);
-            SASSERT(a);
+            if (!m_literal2assumption.contains(l.index())) l.neg();
+            SASSERT(m_literal2assumption.contains(l.index()));
+            expr * a = m_literal2assumption[l.index()];
             if (!already_found_assumptions.contains(a)) {
                 already_found_assumptions.insert(a);
                 m_unsat_core.push_back(a);
@@ -2952,7 +2952,7 @@ namespace smt {
        \brief Execute some finalization code after performing the search.
     */
     void context::check_finalize(lbool r) {
-        TRACE("after_search", display(tout););
+        TRACE("after_search", display(tout << "result: " << r << "\n"););
         display_profile(verbose_stream());
     }
 
@@ -3239,7 +3239,7 @@ namespace smt {
                 for (; it != end; ++it)
                     (*it)->restart_eh();
                 TRACE("mbqi_bug_detail", tout << "before instantiating quantifiers...\n";); 
-                m_qmanager->restart_eh();
+                m_qmanager->restart_eh();                
             }
             if (m_fparams.m_simplify_clauses)
                 simplify_clauses();
@@ -3939,6 +3939,20 @@ namespace smt {
         return th->get_value(n, value);
     }
 
+    bool context::update_model(bool refinalize) {
+        final_check_status fcs = FC_DONE;
+        if (refinalize) {
+            fcs = final_check();
+        }
+        TRACE("opt", tout << (refinalize?"refinalize":"no-op") << " " << fcs << "\n";);
+        if (fcs == FC_DONE) {
+            mk_proto_model(l_true);
+            m_model = m_proto_model->mk_model();
+        }
+            
+        return fcs == FC_DONE;
+    }
+
     void context::mk_proto_model(lbool r) {
         TRACE("get_model", 
               display(tout); 

src/smt/smt_context.h

@@ -206,9 +206,9 @@ namespace smt {
         // Unsat core extraction
         //
         // -----------------------------------
-        typedef u_map<expr *>  bool_var2assumption;
-        bool_var_vector             m_assumptions;
-        bool_var2assumption         m_bool_var2assumption; // maps an expression associated with a literal to the original assumption
+        typedef u_map<expr *>  literal2assumption;
+        literal_vector             m_assumptions;
+        literal2assumption         m_literal2assumption; // maps an expression associated with a literal to the original assumption
         expr_ref_vector             m_unsat_core;
 
         // -----------------------------------
@@ -337,6 +337,10 @@ namespace smt {
             return get_assignment(literal(v));
         }
 
+        literal_vector const & assigned_literals() const { 
+            return m_assigned_literals; 
+        }
+
         lbool get_assignment(expr * n) const;
 
         // Similar to get_assignment, but returns l_undef if n is not internalized.
@@ -1391,6 +1395,8 @@ namespace smt {
         
         void get_model(model_ref & m) const;
 
+        bool update_model(bool refinalize);
+
         void get_proto_model(proto_model_ref & m) const;
         
         unsigned get_num_asserted_formulas() const { return m_asserted_formulas.get_num_formulas(); }

src/smt/smt_context_pp.cpp

@@ -427,7 +427,7 @@ namespace smt {
         }
         expr_ref n(m_manager);
         literal2expr(~consequent, n);
-        pp.display(out, n);
+        pp.display_smt2(out, n);
     }
 
     static unsigned g_lemma_id = 0;
@@ -437,9 +437,9 @@ namespace smt {
     void context::display_lemma_as_smt_problem(unsigned num_antecedents, literal const * antecedents, literal consequent, const char * logic) const {
         char buffer[BUFFER_SZ];
 #ifdef _WINDOWS
-        sprintf_s(buffer, BUFFER_SZ, "lemma_%d.smt", g_lemma_id);
+        sprintf_s(buffer, BUFFER_SZ, "lemma_%d.smt2", g_lemma_id);
 #else
-        sprintf(buffer, "lemma_%d.smt", g_lemma_id);
+        sprintf(buffer, "lemma_%d.smt2", g_lemma_id);
 #endif
         std::ofstream out(buffer);
         display_lemma_as_smt_problem(out, num_antecedents, antecedents, consequent, logic);
@@ -468,7 +468,7 @@ namespace smt {
         }
         expr_ref n(m_manager);
         literal2expr(~consequent, n);
-        pp.display(out, n);
+        pp.display_smt2(out, n);
     }
 
     void context::display_lemma_as_smt_problem(unsigned num_antecedents, literal const * antecedents, 
@@ -476,9 +476,9 @@ namespace smt {
                                                literal consequent, const char * logic) const {
         char buffer[BUFFER_SZ];
 #ifdef _WINDOWS
-        sprintf_s(buffer, BUFFER_SZ, "lemma_%d.smt", g_lemma_id);
+        sprintf_s(buffer, BUFFER_SZ, "lemma_%d.smt2", g_lemma_id);
 #else
-        sprintf(buffer, "lemma_%d.smt", g_lemma_id);
+        sprintf(buffer, "lemma_%d.smt2", g_lemma_id);
 #endif
         std::ofstream out(buffer);
         display_lemma_as_smt_problem(out, num_antecedents, antecedents, num_eq_antecedents, eq_antecedents, consequent, logic);

src/smt/smt_farkas_util.cpp

@@ -0,0 +1,356 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    smt_farkas_util.cpp
+
+Abstract:
+
+    Utility for combining inequalities using coefficients obtained from Farkas lemmas.
+    
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2013-11-2.
+
+Revision History:
+
+    NB. This utility is specialized to proofs generated by the arithmetic solvers.
+
+--*/
+
+#include "smt_farkas_util.h"
+#include "ast_pp.h"
+#include "th_rewriter.h"
+#include "bool_rewriter.h"
+
+
+namespace smt {
+
+    farkas_util::farkas_util(ast_manager& m):
+        m(m),
+        a(m),
+        m_ineqs(m),
+        m_split_literals(false),
+        m_time(0) {
+    }
+
+    void farkas_util::mk_coerce(expr*& e1, expr*& e2) {
+        if (a.is_int(e1) && a.is_real(e2)) {
+            e1 = a.mk_to_real(e1);
+        }
+        else if (a.is_int(e2) && a.is_real(e1)) {
+            e2 = a.mk_to_real(e2);
+        }
+    }
+
+    // TBD: arith_decl_util now supports coercion, so this should be deprecated.
+    app* farkas_util::mk_add(expr* e1, expr* e2) {
+        mk_coerce(e1, e2);
+        return a.mk_add(e1, e2);
+    }
+
+    app* farkas_util::mk_mul(expr* e1, expr* e2) {
+        mk_coerce(e1, e2);
+        return a.mk_mul(e1, e2);
+    }
+
+    app* farkas_util::mk_le(expr* e1, expr* e2) {
+        mk_coerce(e1, e2);
+        return a.mk_le(e1, e2);
+    }
+
+    app* farkas_util::mk_ge(expr* e1, expr* e2) {
+        mk_coerce(e1, e2);
+        return a.mk_gt(e1, e2);
+    }
+
+    app* farkas_util::mk_gt(expr* e1, expr* e2) {
+        mk_coerce(e1, e2);
+        return a.mk_gt(e1, e2);
+    }
+
+    app* farkas_util::mk_lt(expr* e1, expr* e2) {
+        mk_coerce(e1, e2);
+        return a.mk_lt(e1, e2);
+    }
+
+    void farkas_util::mul(rational const& c, expr* e, expr_ref& res) {
+        expr_ref tmp(m);
+        if (c.is_one()) {
+            tmp = e;
+        }
+        else {
+            tmp = mk_mul(a.mk_numeral(c, c.is_int() && a.is_int(e)), e);
+        }
+        res = mk_add(res, tmp);
+    }
+
+    bool farkas_util::is_int_sort(app* c) {
+        SASSERT(m.is_eq(c) || a.is_le(c) || a.is_lt(c) || a.is_gt(c) || a.is_ge(c));
+        SASSERT(a.is_int(c->get_arg(0)) || a.is_real(c->get_arg(0)));
+        return a.is_int(c->get_arg(0));
+    }
+
+    bool farkas_util::is_int_sort() {
+        SASSERT(!m_ineqs.empty());
+        return is_int_sort(m_ineqs[0].get());
+    }
+
+    void farkas_util::normalize_coeffs() {
+        rational l(1);
+        for (unsigned i = 0; i < m_coeffs.size(); ++i) {
+            l = lcm(l, denominator(m_coeffs[i]));
+        }
+        if (!l.is_one()) {
+            for (unsigned i = 0; i < m_coeffs.size(); ++i) {
+                m_coeffs[i] *= l;
+            }
+        }
+        m_normalize_factor = l;
+    }
+
+    app* farkas_util::mk_one() {
+        return a.mk_numeral(rational(1), true);
+    }
+
+    app* farkas_util::fix_sign(bool is_pos, app* c) {
+        expr* x, *y;
+        SASSERT(m.is_eq(c) || a.is_le(c) || a.is_lt(c) || a.is_gt(c) || a.is_ge(c));
+        bool is_int = is_int_sort(c);
+        if (is_int && is_pos && (a.is_lt(c, x, y) || a.is_gt(c, y, x))) {
+            return mk_le(mk_add(x, mk_one()), y);
+        }
+        if (is_int && !is_pos && (a.is_le(c, x, y) || a.is_ge(c, y, x))) {
+            // !(x <= y) <=> x > y <=> x >= y + 1
+            return mk_ge(x, mk_add(y, mk_one()));
+        }
+        if (is_pos) {
+            return c;
+        }
+        if (a.is_le(c, x, y)) return mk_gt(x, y);
+        if (a.is_lt(c, x, y)) return mk_ge(x, y);
+        if (a.is_ge(c, x, y)) return mk_lt(x, y);
+        if (a.is_gt(c, x, y)) return mk_le(x, y);
+        UNREACHABLE();
+        return c;
+    }
+
+    void farkas_util::partition_ineqs() {
+        m_reps.reset();
+        m_his.reset();
+        ++m_time;
+        for (unsigned i = 0; i < m_ineqs.size(); ++i) {
+            m_reps.push_back(process_term(m_ineqs[i].get()));
+        }
+        unsigned head = 0;
+        while (head < m_ineqs.size()) {
+            unsigned r = find(m_reps[head]);
+            unsigned tail = head;
+            for (unsigned i = head+1; i < m_ineqs.size(); ++i) {
+                if (find(m_reps[i]) == r) {
+                    ++tail;
+                    if (tail != i) {
+                        SASSERT(tail < i);
+                        std::swap(m_reps[tail], m_reps[i]);
+                        app_ref tmp(m);
+                        tmp = m_ineqs[i].get();
+                        m_ineqs[i] = m_ineqs[tail].get();
+                        m_ineqs[tail] = tmp;
+                        std::swap(m_coeffs[tail], m_coeffs[i]);
+                    }
+                }
+            }
+            head = tail + 1;
+            m_his.push_back(head);
+        }
+    }
+
+    unsigned farkas_util::find(unsigned idx) {
+        if (m_ts.size() <= idx) {
+            m_roots.resize(idx+1);
+            m_size.resize(idx+1);
+            m_ts.resize(idx+1); 
+            m_roots[idx] = idx;
+            m_ts[idx] = m_time;
+            m_size[idx] = 1;                
+            return idx;
+        }
+        if (m_ts[idx] != m_time) {
+            m_size[idx] = 1;
+            m_ts[idx]    = m_time;
+            m_roots[idx] = idx;
+            return idx;
+        }
+        while (true) {
+            if (m_roots[idx] == idx) {
+                return idx;
+            }
+            idx = m_roots[idx];
+        }
+    }
+
+    void farkas_util::merge(unsigned i, unsigned j) {
+        i = find(i);
+        j = find(j);
+        if (i == j) {
+            return;
+        }
+        if (m_size[i] > m_size[j]) {
+            std::swap(i, j);
+        }
+        m_roots[i] = j;
+        m_size[j] += m_size[i];
+    }
+    unsigned farkas_util::process_term(expr* e) {
+        unsigned r = e->get_id();
+        ptr_vector<expr> todo;
+        ast_mark mark;
+        todo.push_back(e);
+        while (!todo.empty()) {
+            e = todo.back();
+            todo.pop_back();
+            if (mark.is_marked(e)) {
+                continue;
+            }
+            mark.mark(e, true);
+            if (is_uninterp(e)) {
+                merge(r, e->get_id());
+            }
+            if (is_app(e)) {
+                app* a = to_app(e);
+                for (unsigned i = 0; i < a->get_num_args(); ++i) {
+                    todo.push_back(a->get_arg(i));
+                }
+            }
+        }
+        return r;
+    }
+    expr_ref farkas_util::extract_consequence(unsigned lo, unsigned hi) {
+        bool is_int = is_int_sort();
+        app_ref zero(a.mk_numeral(rational::zero(), is_int), m);
+        expr_ref res(m);
+        res = zero;
+        bool is_strict = false;
+        bool is_eq     = true;
+        expr* x, *y;
+        for (unsigned i = lo; i < hi; ++i) {
+            app* c = m_ineqs[i].get();
+            if (m.is_eq(c, x, y)) {
+                mul(m_coeffs[i],  x, res);
+                mul(-m_coeffs[i], y, res);
+            }
+            if (a.is_lt(c, x, y) || a.is_gt(c, y, x)) {
+                mul(m_coeffs[i],  x, res);
+                mul(-m_coeffs[i], y, res);
+                is_strict = true;
+                is_eq = false;
+            }
+            if (a.is_le(c, x, y) || a.is_ge(c, y, x)) {
+                mul(m_coeffs[i],  x, res);
+                mul(-m_coeffs[i], y, res);
+                is_eq = false;
+            }
+        }
+        
+        zero = a.mk_numeral(rational::zero(), a.is_int(res));
+        if (is_eq) {
+            res = m.mk_eq(res, zero);
+        }
+        else if (is_strict) {
+            res = mk_lt(res, zero);
+        }
+        else {
+            res = mk_le(res, zero);
+        }            
+        res = m.mk_not(res);
+        th_rewriter rw(m);
+        params_ref params;
+        params.set_bool("gcd_rounding", true);
+        rw.updt_params(params);
+        proof_ref pr(m);
+        expr_ref result(m);
+        rw(res, result, pr);
+        fix_dl(result);
+        return result;            
+    }
+
+    void farkas_util::fix_dl(expr_ref& r) {
+        expr* e;
+        if (m.is_not(r, e)) {
+            r = e;
+            fix_dl(r);
+            r = m.mk_not(r);
+            return;
+        }
+        expr* e1, *e2, *e3, *e4;
+        if ((m.is_eq(r, e1, e2) || a.is_lt(r, e1, e2) || a.is_gt(r, e1, e2) || 
+             a.is_le(r, e1, e2) || a.is_ge(r, e1, e2))) {
+            if (a.is_add(e1, e3, e4) && a.is_mul(e3)) {
+                r = m.mk_app(to_app(r)->get_decl(), a.mk_add(e4,e3), e2);
+            }
+        }
+    }
+
+    void farkas_util::reset() {
+        m_ineqs.reset();
+        m_coeffs.reset();        
+    }
+    
+    void farkas_util::add(rational const & coef, app * c) {
+        bool is_pos = true;
+        expr* e;
+        while (m.is_not(c, e)) {
+            is_pos = !is_pos;
+            c = to_app(e);
+        }
+        
+        if (!coef.is_zero() && !m.is_true(c)) {
+            m_coeffs.push_back(coef);                
+            m_ineqs.push_back(fix_sign(is_pos, c));                
+        }
+    }
+    
+    expr_ref farkas_util::get() {
+        m_normalize_factor = rational::one();
+        expr_ref res(m);
+        if (m_coeffs.empty()) {
+            res = m.mk_false();
+            return res;
+        }
+        bool is_int = is_int_sort();
+        if (is_int) {                
+            normalize_coeffs();
+        }
+        
+        if (m_split_literals) {
+            // partition equalities into variable disjoint sets.
+            // take the conjunction of these instead of the
+            // linear combination.
+            partition_ineqs();
+            expr_ref_vector lits(m);
+            unsigned lo = 0;
+            for (unsigned i = 0; i < m_his.size(); ++i) {
+                unsigned hi = m_his[i];
+                lits.push_back(extract_consequence(lo, hi));
+                lo = hi;
+            }
+            bool_rewriter(m).mk_or(lits.size(), lits.c_ptr(), res);
+            IF_VERBOSE(2, { if (lits.size() > 1) { verbose_stream() << "combined lemma: " << mk_pp(res, m) << "\n"; } });
+        }
+        else {
+            res = extract_consequence(0, m_coeffs.size());
+        }
+
+        TRACE("arith", 
+              for (unsigned i = 0; i < m_coeffs.size(); ++i) {
+                  tout << m_coeffs[i] << " * (" << mk_pp(m_ineqs[i].get(), m) << ") ";
+              }
+              tout << "\n";
+              tout << res << "\n";
+              );
+
+        return res;
+    }
+}
+

src/smt/smt_farkas_util.h

@@ -0,0 +1,96 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    smt_farkas_util.h
+
+Abstract:
+
+    Utility for combining inequalities using coefficients obtained from Farkas lemmas.
+    
+Author:
+
+    Nikolaj Bjorner (nbjorne) 2013-11-2.
+
+Revision History:
+
+    NB. This utility is specialized to proofs generated by the arithmetic solvers.
+
+--*/
+
+#ifndef __FARKAS_UTIL_H_
+#define __FARKAS_UTIL_H_
+
+#include "arith_decl_plugin.h"
+
+namespace smt {
+
+    class farkas_util {
+        ast_manager&     m;
+        arith_util       a;
+        app_ref_vector   m_ineqs;
+        vector<rational> m_coeffs;
+        rational         m_normalize_factor;
+        // utilities for separating coefficients
+        bool             m_split_literals;
+        unsigned         m_time;
+        unsigned_vector  m_roots, m_size, m_his, m_reps, m_ts;
+
+        void mk_coerce(expr*& e1, expr*& e2);
+        app* mk_add(expr* e1, expr* e2);
+        app* mk_mul(expr* e1, expr* e2);
+        app* mk_le(expr* e1, expr* e2);
+        app* mk_ge(expr* e1, expr* e2);
+        app* mk_gt(expr* e1, expr* e2);
+        app* mk_lt(expr* e1, expr* e2);
+        void mul(rational const& c, expr* e, expr_ref& res);
+        bool is_int_sort(app* c);
+        bool is_int_sort();
+        void normalize_coeffs();
+        app* mk_one();
+        app* fix_sign(bool is_pos, app* c);
+        void partition_ineqs();
+        unsigned find(unsigned idx);
+        void merge(unsigned i, unsigned j);
+        unsigned process_term(expr* e);
+        expr_ref extract_consequence(unsigned lo, unsigned hi);
+        void fix_dl(expr_ref& r);
+
+    public:
+        farkas_util(ast_manager& m);
+
+        /**
+           \brief Reset state
+         */
+        void reset();
+
+        /** 
+            \brief add a multiple of constraint c to the current state 
+         */
+        void add(rational const & coef, app * c);
+
+        /**
+           \brief Extract the complement of premises multiplied by Farkas coefficients.
+        */
+        expr_ref get();
+
+        /**
+           \brief Coefficients are normalized for integer problems.
+           Retrieve multiplicant for normalization.
+         */
+        rational const& get_normalize_factor() const { return m_normalize_factor; }
+
+        /**
+           \brief extract one or multiple consequences based on literal partition.
+           Multiple consequences are strongst modulo a partition of variables.
+           Consequence generation under literal partitioning maintains difference logic constraints.
+           That is, if the original constraints are difference logic, then the consequent 
+           produced by literal partitioning is also difference logic. 
+         */
+        void set_split_literals(bool f) { m_split_literals = f; }
+
+    };
+}
+
+#endif

src/smt/smt_internalizer.cpp

@@ -814,6 +814,7 @@ namespace smt {
     */
     bool_var context::mk_bool_var(expr * n) {
         SASSERT(!b_internalized(n));
+        //SASSERT(!m_manager.is_not(n));
         unsigned id = n->get_id();
         bool_var v  = m_b_internalized_stack.size();
 #ifndef _EXTERNAL_RELEASE 
@@ -828,7 +829,7 @@ namespace smt {
         }
 #endif
         TRACE("mk_bool_var", tout << "creating boolean variable: " << v << " for:\n" << mk_pp(n, m_manager) << "\n";);
-        TRACE("mk_var_bug", tout << "mk_bool: " << v << "\n";);
+        TRACE("mk_var_bug", tout << "mk_bool: " << v << "\n";);        
         set_bool_var(id, v);
         m_bdata.reserve(v+1);
         m_activity.reserve(v+1);

src/smt/smt_justification.h

@@ -225,6 +225,7 @@ namespace smt {
         virtual proof * mk_proof(conflict_resolution & cr) = 0;
 
         virtual char const * get_name() const { return "simple"; }
+
     };
 
     class simple_theory_justification : public simple_justification {
@@ -274,6 +275,7 @@ namespace smt {
 
 
         virtual char const * get_name() const { return "theory-propagation"; }
+        
     };
      
     class theory_conflict_justification : public simple_theory_justification {

src/smt/smt_setup.cpp

@@ -30,6 +30,7 @@ Revision History:
 #include"theory_dummy.h"
 #include"theory_dl.h"
 #include"theory_seq_empty.h"
+#include"theory_pb.h"
 
 namespace smt {
 
@@ -300,7 +301,7 @@ namespace smt {
     }
 
     void setup::setup_QF_IDL() {
-        TRACE("setup", tout << "setup_QF_IDL(st)\n";);
+        TRACE("setup", tout << "setup_QF_IDL()\n";);
         m_params.m_relevancy_lvl       = 0;
         m_params.m_arith_expand_eqs    = true;
         m_params.m_arith_reflect       = false;
@@ -354,6 +355,10 @@ namespace smt {
                 m_context.register_plugin(alloc(smt::theory_dense_si, m_manager, m_params));
             else
                 m_context.register_plugin(alloc(smt::theory_dense_i, m_manager, m_params));
+
+        }
+        else if (!m_params.m_arith_auto_config_simplex && !is_dense(st)) {
+            m_context.register_plugin(alloc(smt::theory_idl, m_manager, m_params));            
         }
         else {
             // if (st.m_arith_k_sum < rational(INT_MAX / 8)) {
@@ -373,6 +378,7 @@ namespace smt {
         m_params.m_arith_reflect    = false;
         m_params.m_nnf_cnf          = false;
         m_params.m_arith_eq_bounds  = true;
+        m_params.m_arith_expand_eqs = true;
         m_params.m_phase_selection  = PS_ALWAYS_FALSE;
         m_params.m_restart_strategy = RS_GEOMETRIC;
         m_params.m_restart_factor   = 1.5;
@@ -405,6 +411,15 @@ namespace smt {
                 return;
             }
         }
+#if 0
+        switch (m_params.m_arith_mode) {
+        case AS_DIFF_LOGIC:
+        case AS_DENSE_DIFF_LOGIC:
+        case AS_UTVPI:
+            setup_arith();
+            return;
+        }
+#endif
         m_params.m_arith_eq_bounds  = true;
         m_params.m_phase_selection  = PS_ALWAYS_FALSE;
         m_params.m_restart_strategy = RS_GEOMETRIC;
@@ -414,7 +429,7 @@ namespace smt {
             m_context.register_plugin(alloc(smt::theory_mi_arith, m_manager, m_params));
         }
         // else if (st.m_arith_k_sum < rational(INT_MAX / 8))
-        //    m_context.register_plugin(alloc(smt::theory_si_arith, m_manager, m_params));
+        //   m_context.register_plugin(alloc(smt::theory_dense_si, m_manager, m_params));
         else
             m_context.register_plugin(alloc(smt::theory_i_arith, m_manager, m_params));
     }
@@ -688,7 +703,12 @@ namespace smt {
     }
 
     void setup::setup_mi_arith() {
-        m_context.register_plugin(alloc(smt::theory_mi_arith, m_manager, m_params));
+        if (m_params.m_arith_mode == AS_OPTINF) {
+            m_context.register_plugin(alloc(smt::theory_inf_arith, m_manager, m_params));            
+        }
+        else {
+            m_context.register_plugin(alloc(smt::theory_mi_arith, m_manager, m_params));
+        }
     }
 
     void setup::setup_arith() {
@@ -697,6 +717,7 @@ namespace smt {
             m_context.register_plugin(alloc(smt::theory_dummy, m_manager.mk_family_id("arith"), "no arithmetic"));
             break;
         case AS_DIFF_LOGIC:
+            m_params.m_arith_expand_eqs  = true;
             if (m_params.m_arith_fixnum) {
                 if (m_params.m_arith_int_only)
                     m_context.register_plugin(alloc(smt::theory_fidl, m_manager, m_params));
@@ -711,6 +732,7 @@ namespace smt {
             }
             break;
         case AS_DENSE_DIFF_LOGIC:
+            m_params.m_arith_expand_eqs  = true;
             if (m_params.m_arith_fixnum) {
                 if (m_params.m_arith_int_only)
                     m_context.register_plugin(alloc(smt::theory_dense_si, m_manager, m_params));
@@ -725,11 +747,15 @@ namespace smt {
             }
             break;
         case AS_UTVPI:
+            m_params.m_arith_expand_eqs  = true;
             if (m_params.m_arith_int_only)
                 m_context.register_plugin(alloc(smt::theory_iutvpi, m_manager));
             else
                 m_context.register_plugin(alloc(smt::theory_rutvpi, m_manager));          
             break;
+        case AS_OPTINF:
+            m_context.register_plugin(alloc(smt::theory_inf_arith, m_manager, m_params));            
+            break;
         default:
             if (m_params.m_arith_int_only)
                 m_context.register_plugin(alloc(smt::theory_i_arith, m_manager, m_params));
@@ -780,6 +806,10 @@ namespace smt {
         m_context.register_plugin(alloc(theory_seq_empty, m_manager));
     }
 
+    void setup::setup_card() {
+        m_context.register_plugin(alloc(theory_pb, m_manager, m_params));
+    }
+
     void setup::setup_unknown() {
         setup_arith();
         setup_arrays();
@@ -787,6 +817,7 @@ namespace smt {
         setup_datatypes();
         setup_dl();
         setup_seq();
+        setup_card();
     }
 
     void setup::setup_unknown(static_features & st) {

src/smt/smt_setup.h

@@ -42,7 +42,7 @@ namespace smt {
     class setup {
         context &          m_context;
         ast_manager &      m_manager;
-        smt_params & m_params;
+        smt_params &       m_params;
         symbol             m_logic;
         bool               m_already_configured;
         void setup_auto_config();
@@ -92,7 +92,7 @@ namespace smt {
         void setup_arith();
         void setup_dl();
         void setup_seq();
-        void setup_instgen();
+        void setup_card();
         void setup_i_arith();
         void setup_mi_arith();
     public:

src/smt/spanning_tree.h

@@ -0,0 +1,82 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    spanning_tree.h
+
+Abstract:
+
+    Represent spanning trees with needed operations for Network Simplex
+
+Author:
+
+    Anh-Dung Phan (t-anphan) 2013-11-06
+
+Notes:
+   
+--*/
+#ifndef _SPANNING_TREE_H_
+#define _SPANNING_TREE_H_
+
+#include "diff_logic.h"
+#include "spanning_tree_base.h"
+
+namespace smt {
+
+    template<typename Ext>
+    class thread_spanning_tree : public spanning_tree_base, protected Ext {
+    protected:        
+        typedef dl_edge<Ext> edge;
+        typedef dl_graph<Ext> graph;     
+        typedef typename Ext::numeral numeral;
+        typedef typename Ext::fin_numeral fin_numeral;
+
+        // Store the parent of a node i in the spanning tree
+        svector<node_id> m_pred;
+        // Store the number of edge on the path from node i to the root
+        svector<int> m_depth;
+        svector<node_id> m_thread;           // Store the pointer from node i to the next node in depth-first search order
+
+        svector<edge_id> m_tree;          // i |-> edge between (i, m_pred[i])
+
+        node_id m_root_t2;
+
+        graph & m_graph;
+
+        void swap_order(node_id q, node_id v);
+        node_id find_rev_thread(node_id n) const;
+        void fix_depth(node_id start, node_id after_end);
+        node_id get_final(int start);
+        bool is_preorder_traversal(node_id start, node_id end);   
+        node_id get_common_ancestor(node_id u, node_id v);
+        bool is_forward_edge(edge_id e_id) const;
+        bool is_ancestor_of(node_id ancestor, node_id child);
+
+    public:      
+        thread_spanning_tree(graph & g);
+
+        virtual void initialize(svector<edge_id> const & tree);
+        void get_descendants(node_id start, svector<node_id> & descendants);
+        
+        virtual void update(edge_id enter_id, edge_id leave_id);                
+        void get_path(node_id start, node_id end, svector<edge_id> & path, svector<bool> & against);              
+        bool in_subtree_t2(node_id child);
+
+        bool check_well_formed();        
+    };
+
+    template<typename Ext>
+    class basic_spanning_tree : public thread_spanning_tree<Ext> {
+    private:
+        graph * m_tree_graph;
+
+    public:
+        basic_spanning_tree(graph & g);
+        void initialize(svector<edge_id> const & tree);
+        void update(edge_id enter_id, edge_id leave_id);
+    };
+
+}
+
+#endif

src/smt/spanning_tree_base.h

@@ -0,0 +1,53 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    spanning_tree_base.h
+
+Abstract:
+
+    Represent spanning trees with needed operations for Network Simplex
+
+Author:
+
+    Anh-Dung Phan (t-anphan) 2013-11-06
+
+Notes:
+   
+--*/
+
+#ifndef _SPANNING_TREE_BASE_H_
+#define _SPANNING_TREE_BASE_H_
+
+#include "util.h"
+#include "vector.h"
+
+namespace smt {
+    template<typename TV>
+    inline std::string pp_vector(std::string const & label, TV v) {
+        std::ostringstream oss;
+        oss << label << " ";
+        for (unsigned i = 0; i < v.size(); ++i) {
+            oss << v[i] << " ";
+        }
+        oss << std::endl;
+        return oss.str();
+    }
+
+    class spanning_tree_base { 
+    public:
+        typedef int node_id;
+        typedef int edge_id;
+        virtual void initialize(svector<edge_id> const & tree) = 0;
+        virtual void get_descendants(node_id start, svector<node_id> & descendants) = 0;
+        
+        virtual void update(edge_id enter_id, edge_id leave_id) = 0;
+        virtual void get_path(node_id start, node_id end, svector<edge_id> & path, svector<bool> & against) = 0;
+        virtual bool in_subtree_t2(node_id child) = 0;
+
+        virtual bool check_well_formed() = 0;
+    };
+}
+
+#endif

src/smt/spanning_tree_def.h

@@ -0,0 +1,489 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    spanning_tree_def.h
+
+Abstract:
+
+
+Author:
+
+    Anh-Dung Phan (t-anphan) 2013-11-06
+
+Notes:
+   
+--*/
+
+#ifndef _SPANNING_TREE_DEF_H_
+#define _SPANNING_TREE_DEF_H_
+
+#include "spanning_tree.h"
+
+namespace smt {
+
+    template<typename Ext>
+    thread_spanning_tree<Ext>::thread_spanning_tree(graph & g) :
+        m_graph(g) {
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::initialize(svector<edge_id> const & tree) {
+        m_tree = tree;
+
+        unsigned num_nodes = m_graph.get_num_nodes();
+        m_pred.resize(num_nodes);
+        m_depth.resize(num_nodes);
+        m_thread.resize(num_nodes);
+        
+        node_id root = num_nodes - 1;
+        m_pred[root] = -1;
+        m_depth[root] = 0;
+        m_thread[root] = 0;
+            
+        // Create artificial edges from/to root node to/from other nodes and initialize the spanning tree
+        for (int i = 0; i < root; ++i) {
+            m_pred[i] = root;
+            m_depth[i] = 1;
+            m_thread[i] = i + 1;
+        }
+
+        TRACE("network_flow", {
+            tout << pp_vector("Predecessors", m_pred) << pp_vector("Threads", m_thread); 
+            tout << pp_vector("Depths", m_depth) << pp_vector("Tree", m_tree);
+        });
+    }
+
+    template<typename Ext>
+    typename thread_spanning_tree<Ext>::node_id thread_spanning_tree<Ext>::get_common_ancestor(node_id u, node_id v) {
+        while (u != v) {
+            if (m_depth[u] > m_depth[v])
+                u = m_pred[u];
+            else 
+                v = m_pred[v];                    
+        }
+        return u;
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::get_path(node_id start, node_id end, svector<edge_id> & path, svector<bool> & against) {
+        node_id join = get_common_ancestor(start, end);
+        path.reset();
+        while (start != join) {
+            edge_id e_id = m_tree[start];
+            path.push_back(e_id);
+            against.push_back(is_forward_edge(e_id));
+            start = m_pred[start];
+        }
+        while (end != join) {
+            edge_id e_id = m_tree[end];
+            path.push_back(e_id);
+            against.push_back(!is_forward_edge(e_id));
+            end = m_pred[end];
+        }
+    }
+
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::is_forward_edge(edge_id e_id) const {
+        node_id start = m_graph.get_source(e_id);
+        node_id end = m_graph.get_target(e_id);
+        SASSERT(m_pred[start] == end || m_pred[end] == start);
+        return m_pred[start] == end;
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::get_descendants(node_id start, svector<node_id> & descendants) {
+        descendants.reset();        
+        descendants.push_back(start);
+        node_id u = m_thread[start];
+        while (m_depth[u] > m_depth[start]) {         
+            descendants.push_back(u);
+            u = m_thread[u];
+        }
+        
+    }
+
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::in_subtree_t2(node_id child) {
+        if (m_depth[child] < m_depth[m_root_t2]) {
+            return false;
+        }
+        return is_ancestor_of(m_root_t2, child);
+    }
+
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::is_ancestor_of(node_id ancestor, node_id child) {
+        for (node_id n = child; n != -1; n = m_pred[n]) {
+            if (n == ancestor) {
+                return true;
+            }
+        }
+        return false;
+    }
+        
+    /**
+        \brief add entering_edge, remove leaving_edge from spanning tree.
+
+        Old tree:                    New tree:
+                       root                      root
+                     /      \                  /      \
+                    x       y                  x       y
+                  /   \    / \               /   \    / \
+                       u     s                   u      s
+                       |    /                          /
+                       v    w                    v    w
+                      / \    \                  / \    \
+                         z   p                     z   p
+                          \                         \ /
+                           q                         q
+        */
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::update(edge_id enter_id, edge_id leave_id) {    
+        node_id p = m_graph.get_source(enter_id);
+        node_id q = m_graph.get_target(enter_id);
+        node_id u = m_graph.get_source(leave_id);
+        node_id v = m_graph.get_target(leave_id);
+        
+        if (m_pred[u] == v) {
+            std::swap(u, v);
+        }
+
+        SASSERT(m_pred[v] == u);         
+
+        if (is_ancestor_of(v, p)) {
+            std::swap(p, q);    
+        }
+
+        SASSERT(is_ancestor_of(v, q));
+        
+        TRACE("network_flow", { 
+            tout << "update_spanning_tree: (" << p << ", " << q << ") enters, (";
+            tout << u << ", " << v << ") leaves\n";
+        });
+
+        // Old threads: alpha -> v -*-> f(v) -> beta | p -*-> f(p) -> gamma
+        // New threads: alpha -> beta                | p -*-> f(p) -> v -*-> f(v) -> gamma
+
+        node_id f_p = get_final(p);
+        node_id f_v = get_final(v);
+        node_id alpha = find_rev_thread(v);
+        node_id beta = m_thread[f_v];
+        node_id gamma = m_thread[f_p];
+
+        if (v != gamma) {
+            m_thread[alpha] = beta;
+            m_thread[f_p] = v;
+            m_thread[f_v] = gamma;
+        }
+
+        node_id old_pred = m_pred[q]; 
+        // Update stem nodes from q to v
+        if (q != v) {
+            for (node_id n = q; n != v; ) {
+                SASSERT(old_pred != u); // the last processed node_id is v
+                SASSERT(-1 != m_pred[old_pred]);
+                int next_old_pred = m_pred[old_pred];  
+                swap_order(n, old_pred);
+                m_tree[old_pred] = m_tree[n];
+                n = old_pred;
+                old_pred = next_old_pred;
+            }     
+        }
+        
+        m_pred[q] = p; 
+        m_tree[q] = enter_id;
+        m_root_t2 = q;
+
+        node_id after_final_q = (v == gamma) ? beta : gamma;
+        fix_depth(q, after_final_q);
+
+        SASSERT(!in_subtree_t2(p));
+        SASSERT(in_subtree_t2(q));
+        SASSERT(!in_subtree_t2(u));
+        SASSERT(in_subtree_t2(v));
+
+        TRACE("network_flow", {
+            tout << pp_vector("Predecessors", m_pred) << pp_vector("Threads", m_thread); 
+            tout << pp_vector("Depths", m_depth) << pp_vector("Tree", m_tree);
+            });
+    }
+
+    /**
+        swap v and q in tree.
+        - fixup m_thread
+        - fixup m_pred
+
+        Case 1: final(q) == final(v)
+        -------
+        Old thread: prev -> v -*-> alpha -> q -*-> final(q) -> next
+        New thread: prev -> q -*-> final(q) -> v -*-> alpha -> next
+
+        Case 2: final(q) != final(v)
+        -------
+        Old thread: prev -> v -*-> alpha -> q -*-> final(q) -> beta -*-> final(v) -> next
+        New thread: prev -> q -*-> final(q) -> v -*-> alpha -> beta -*-> final(v) -> next
+                 
+    */
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::swap_order(node_id q, node_id v) {
+        SASSERT(q != v);
+        SASSERT(m_pred[q] == v);        
+        SASSERT(is_preorder_traversal(v, get_final(v)));
+        node_id prev = find_rev_thread(v);
+        node_id f_q = get_final(q);
+        node_id f_v = get_final(v);
+        node_id next = m_thread[f_v];
+        node_id alpha = find_rev_thread(q);
+
+        if (f_q == f_v) {
+            SASSERT(f_q != v && alpha != next);
+            m_thread[f_q] = v;
+            m_thread[alpha] = next;
+            f_q = alpha;
+        }
+        else {            
+            node_id beta = m_thread[f_q];
+            SASSERT(f_q != v && alpha != beta);
+            m_thread[f_q] = v;
+            m_thread[alpha] = beta;
+            f_q = f_v;
+        }
+        SASSERT(prev != q);
+        m_thread[prev] = q;
+        m_pred[v] = q;
+        // Notes: f_q has to be used since m_depth hasn't been updated yet.
+        SASSERT(is_preorder_traversal(q, f_q));
+    }
+
+    /**
+        \brief Check invariants of main data-structures.
+
+        Spanning tree of m_graph + root is represented using:
+        
+        svector<edge_state> m_states;      edge_id |-> edge_state
+        svector<node_id> m_pred;              node_id |-> node
+        svector<int>     m_depth;             node_id |-> int
+        svector<node_id> m_thread;            node_id |-> node
+
+        Tree is determined by m_pred:
+        - m_pred[root] == -1
+        - m_pred[n] = m != n     for each node_id n, acyclic until reaching root.
+        - m_depth[m_pred[n]] + 1 == m_depth[n] for each n != root        
+
+        m_thread is a linked list traversing all nodes.
+        Furthermore, the nodes linked in m_thread follows a 
+        depth-first traversal order.
+        
+    */
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::check_well_formed() {
+        node_id root = m_pred.size()-1;
+
+        // Check that m_thread traverses each node.
+        // This gets checked using union-find as well.
+        svector<bool> found(m_thread.size(), false);
+        found[root] = true;
+        for (node_id x = m_thread[root]; x != root; x = m_thread[x]) {
+            SASSERT(x != m_thread[x]);
+            found[x] = true;
+        }
+        for (unsigned i = 0; i < found.size(); ++i) {
+            SASSERT(found[i]);
+        }
+
+        // m_pred is acyclic, and points to root.
+        SASSERT(m_pred[root] == -1);
+        SASSERT(m_depth[root] == 0);
+        for (node_id i = 0; i < root; ++i) {
+            SASSERT(m_depth[m_pred[i]] < m_depth[i]);
+        }
+
+        // m_depth[x] denotes distance from x to the root node
+        for (node_id x = m_thread[root]; x != root; x = m_thread[x]) {
+            SASSERT(m_depth[x] > 0);
+            SASSERT(m_depth[x] == m_depth[m_pred[x]] + 1);
+        }
+
+        // m_thread forms a spanning tree over [0..root]
+        // Union-find structure
+        svector<int> roots(m_pred.size(), -1);
+                      
+        for (node_id x = m_thread[root]; x != root; x = m_thread[x]) {
+            node_id y = m_pred[x];
+            // We are now going to check the edge between x and y
+            SASSERT(find(roots, x) != find(roots, y));
+            merge(roots, x, y);
+        }
+
+        // All nodes belong to the same spanning tree
+        for (unsigned i = 0; i < roots.size(); ++i) {            
+            SASSERT(roots[i] + roots.size() == 0 || roots[i] >= 0);            
+        }   
+
+        for (unsigned i = 0; i < m_tree.size(); ++i) {           
+            node_id src = m_graph.get_source(m_tree[i]);
+            node_id tgt = m_graph.get_target(m_tree[i]);
+            SASSERT(m_pred[src] == tgt || m_pred[tgt] == src);            
+        }  
+
+        return true;
+    }
+
+    static unsigned find(svector<int>& roots, unsigned x) {
+        unsigned old_x = x;
+        while (roots[x] >= 0) {
+            x = roots[x];
+        }
+        SASSERT(roots[x] < 0);
+        if (old_x != x) {
+            roots[old_x] = x;
+        }
+        return x;
+    }
+
+    static void merge(svector<int>& roots, unsigned x, unsigned y) {
+        x = find(roots, x);
+        y = find(roots, y);
+        SASSERT(roots[x] < 0 && roots[y] < 0);
+        if (x == y) {
+            return;
+        }
+        if (roots[x] > roots[y]) {
+            std::swap(x, y);
+        }
+        SASSERT(roots[x] <= roots[y]);
+        roots[x] += roots[y];
+        roots[y] = x;
+    }
+
+    /**
+        \brief find node_id that points to 'n' in m_thread
+    */
+    template<typename Ext>
+    typename thread_spanning_tree<Ext>::node_id thread_spanning_tree<Ext>::find_rev_thread(node_id n) const {     
+        node_id ancestor = m_pred[n];
+        SASSERT(ancestor != -1);
+        while (m_thread[ancestor] != n) {
+            ancestor = m_thread[ancestor];
+        }
+        return ancestor;
+    }
+
+    template<typename Ext>
+    void thread_spanning_tree<Ext>::fix_depth(node_id start, node_id after_end) {
+        while (start != after_end) {
+            SASSERT(m_pred[start] != -1);
+            m_depth[start] = m_depth[m_pred[start]]+1;
+            start = m_thread[start];
+        }
+    }
+        
+    template<typename Ext>
+    typename thread_spanning_tree<Ext>::node_id thread_spanning_tree<Ext>::get_final(int start) {
+        int n = start;
+        while (m_depth[m_thread[n]] > m_depth[start]) {
+            n = m_thread[n];
+        }
+        return n;
+    }
+
+    template<typename Ext>
+    bool thread_spanning_tree<Ext>::is_preorder_traversal(node_id start, node_id end) {
+        // get children of start
+        uint_set children;
+        children.insert(start);
+        node_id root = m_pred.size()-1;
+        for (int i = 0; i < root; ++i) {
+            for (int j = 0; j < root; ++j) {
+                if (children.contains(m_pred[j])) {
+                    children.insert(j);
+                }
+            }
+        }
+        // visit children using m_thread
+        children.remove(start);
+        do {
+            start = m_thread[start];
+            SASSERT(children.contains(start));
+            children.remove(start);
+        }
+        while (start != end);
+        SASSERT(children.empty());
+        return true;
+    }
+
+    // Basic spanning tree
+    template<typename Ext>
+    basic_spanning_tree<Ext>::basic_spanning_tree(graph & g) : thread_spanning_tree<Ext>(g) {        
+    }
+
+    template<typename Ext>
+    void basic_spanning_tree<Ext>::initialize(svector<edge_id> const & tree) {        
+        m_tree_graph = alloc(graph);
+        m_tree = tree;
+        unsigned num_nodes = m_graph.get_num_nodes();
+        for (unsigned i = 0; i < num_nodes; ++i) {
+            m_tree_graph->init_var(i);
+        }
+
+        vector<edge> const & es = m_graph.get_all_edges();
+        svector<edge_id>::const_iterator it = m_tree.begin(), end = m_tree.end();
+        for(; it != end; ++it) {
+            edge const & e = es[*it];
+            m_tree_graph->add_edge(e.get_source(), e.get_target(), e.get_weight(), explanation());
+        }
+
+        node_id root = num_nodes - 1;
+        m_tree_graph->bfs_undirected(root, m_pred, m_depth);
+        m_tree_graph->dfs_undirected(root, m_thread);
+    }
+
+    template<typename Ext>
+    void basic_spanning_tree<Ext>::update(edge_id enter_id, edge_id leave_id) {
+        if (m_tree_graph) dealloc(m_tree_graph);
+        m_tree_graph = alloc(graph);
+        unsigned num_nodes = m_graph.get_num_nodes();
+        for (unsigned i = 0; i < num_nodes; ++i) {
+            m_tree_graph->init_var(i);
+        }
+
+        vector<edge> const & es = m_graph.get_all_edges();
+        svector<edge_id>::const_iterator it = m_tree.begin(), end = m_tree.end();
+        for(; it != end; ++it) {
+            edge const & e = es[*it];
+            if (leave_id != *it) {
+                m_tree_graph->add_edge(e.get_source(), e.get_target(), e.get_weight(), explanation());
+            }
+        }
+        edge const & e = es[enter_id];
+        m_tree_graph->add_edge(e.get_source(), e.get_target(), e.get_weight(), explanation());
+
+        node_id root = num_nodes - 1;
+        m_tree_graph->bfs_undirected(root, m_pred, m_depth);
+        m_tree_graph->dfs_undirected(root, m_thread);
+
+        vector<edge> const & tree_edges = m_tree_graph->get_all_edges();
+        for (unsigned i = 0; i < tree_edges.size(); ++i) {
+            edge const & e = tree_edges[i];
+            dl_var src = e.get_source();
+            dl_var tgt = e.get_target();
+            edge_id id;
+            VERIFY(m_graph.get_edge_id(src, tgt, id));
+            SASSERT(tgt == m_pred[src] || src == m_pred[tgt]);
+            if (tgt == m_pred[src]) {
+                m_tree[src] = id;
+            }
+            else {
+                m_tree[tgt] = id;
+            }
+        }
+
+        node_id p = m_graph.get_source(enter_id);
+        node_id q = m_graph.get_target(enter_id);
+        m_root_t2 = p == m_pred[q] ? q : p;
+    }
+
+}
+
+#endif

src/smt/theory_arith.cpp

@@ -25,5 +25,6 @@ namespace smt {
     template class theory_arith<i_ext>;
     // template class theory_arith<si_ext>;
     // template class theory_arith<smi_ext>;
-    
+
+    template class smt::theory_arith<inf_ext>;    
 };

src/smt/theory_arith.h

@@ -36,6 +36,8 @@ Revision History:
 #include"grobner.h"
 #include"arith_simplifier_plugin.h"
 #include"arith_eq_solver.h"
+#include"theory_opt.h"
+#include"uint_set.h"
 
 namespace smt {
     
@@ -80,7 +82,7 @@ namespace smt {
     */
 
     template<typename Ext>
-    class theory_arith : public theory, private Ext {
+    class theory_arith : public theory, public theory_opt, private Ext {
     public:
         typedef typename Ext::numeral     numeral;
         typedef typename Ext::inf_numeral inf_numeral;
@@ -90,6 +92,7 @@ namespace smt {
         static const int    dead_row_id = -1;
     protected:
         bool proofs_enabled() const { return get_manager().proofs_enabled(); }
+        bool coeffs_enabled() const { return proofs_enabled() || m_bound_watch != null_bool_var; }
         
         struct linear_monomial {
             numeral     m_coeff;
@@ -262,6 +265,7 @@ namespace smt {
             inf_numeral const & get_value() const { return m_value; }
             virtual bool has_justification() const { return false; }
             virtual void push_justification(antecedents& antecedents, numeral const& coeff, bool proofs_enabled) {}
+            virtual void display(theory_arith const& th, std::ostream& out) const;
         };
 
 
@@ -273,14 +277,14 @@ namespace smt {
         class atom : public bound {
         protected:
             bool_var    m_bvar;
-            numeral     m_k;
+            inf_numeral m_k;
             unsigned    m_atom_kind:2;   // atom kind
             unsigned    m_is_true:1;     // cache: true if the atom was assigned to true.
         public:
-            atom(bool_var bv, theory_var v, numeral const & k, atom_kind kind);
+            atom(bool_var bv, theory_var v, inf_numeral const & k, atom_kind kind);
             atom_kind get_atom_kind() const { return static_cast<atom_kind>(m_atom_kind); }
             virtual ~atom() {}
-            numeral const & get_k() const { return m_k; }
+            inline inf_numeral const & get_k() const { return m_k; }
             bool_var get_bool_var() const { return m_bvar; }
             bool is_true() const { return m_is_true; }
             void assign_eh(bool is_true, inf_numeral const & epsilon);
@@ -288,6 +292,7 @@ namespace smt {
             virtual void push_justification(antecedents& a, numeral const& coeff, bool proofs_enabled) { 
                 a.push_lit(literal(get_bool_var(), !m_is_true), coeff, proofs_enabled); 
             }
+            virtual void display(theory_arith const& th, std::ostream& out) const;
         };
 
         class eq_bound : public bound { 
@@ -306,6 +311,7 @@ namespace smt {
                 SASSERT(m_lhs->get_root() == m_rhs->get_root());
                 a.push_eq(enode_pair(m_lhs, m_rhs), coeff, proofs_enabled); 
             }
+            virtual void display(theory_arith const& th, std::ostream& out) const;
         };
 
         class derived_bound : public bound {
@@ -320,6 +326,7 @@ namespace smt {
             virtual void push_justification(antecedents& a, numeral const& coeff, bool proofs_enabled); 
             virtual void push_lit(literal l, numeral const&) { m_lits.push_back(l); }
             virtual void push_eq(enode_pair const& p, numeral const&) { m_eqs.push_back(p); }
+            virtual void display(theory_arith const& th, std::ostream& out) const;
         };
     
         class justified_derived_bound : public derived_bound {
@@ -433,6 +440,7 @@ namespace smt {
         bool                    m_eager_gcd; // true if gcd should be applied at every add_row
         unsigned                m_final_check_idx;
 
+
         // backtracking
         svector<bound_trail>    m_bound_trail;
         svector<unsigned>       m_unassigned_atoms_trail;
@@ -869,13 +877,29 @@ namespace smt {
         row               m_tmp_row;
 
         void add_tmp_row(row & r1, numeral const & coeff, row const & r2);
-        theory_var pick_var_to_leave(theory_var x_j, bool inc, numeral & a_ij, inf_numeral & gain);
+        theory_var pick_var_to_leave(bool has_int, theory_var x_j, bool inc, numeral & a_ij, inf_numeral & gain, bool& skiped_row);
+        bool is_safe_to_leave(theory_var x, bool inc, bool& has_int, bool& is_shared);
+        bool move_to_bound(theory_var x_i, bool inc);
         template<bool invert>
         void add_tmp_row_entry(row & r, numeral const & coeff, theory_var v);
-        bool max_min(row & r, bool max);
-        bool max_min(theory_var v, bool max);
+        enum max_min_t { UNBOUNDED, AT_BOUND, OPTIMIZED, BEST_EFFORT};
+        max_min_t max_min(theory_var v, bool max, bool& has_shared);
+        max_min_t max_min_orig(row & r, bool max, bool& has_shared);
         bool max_min(svector<theory_var> const & vars);
 
+        max_min_t max_min_new(row& r, bool max, bool& has_shared);
+        bool unbounded_gain(inf_numeral const & max_gain) const;
+        bool safe_gain(inf_numeral const& min_gain, inf_numeral const & max_gain) const;
+        void normalize_gain(numeral const& divisor, inf_numeral & max_gain) const;
+        void init_gains(theory_var x, bool inc, inf_numeral& min_gain, inf_numeral& max_gain);
+        bool update_gains(bool inc, theory_var x_i, numeral const& a_ij, 
+                          inf_numeral& min_gain, inf_numeral& max_gain);
+        bool move_to_bound_new(theory_var x_i, bool inc, unsigned& best_efforts, bool& has_shared);
+        bool pick_var_to_leave(
+            theory_var x_j, bool inc, numeral & a_ij, 
+            inf_numeral& min_gain, inf_numeral& max_gain, 
+            bool& shared, theory_var& x_i);
+
         // -----------------------------------
         //
         // Non linear
@@ -1004,6 +1028,27 @@ namespace smt {
         virtual bool get_value(enode * n, expr_ref & r);
 
 
+        // -----------------------------------
+        //
+        // Optimization
+        //
+        // -----------------------------------
+        virtual inf_eps_rational<inf_rational> maximize(theory_var v, expr_ref& blocker, bool& has_shared);
+        virtual inf_eps_rational<inf_rational> value(theory_var v);
+        virtual theory_var add_objective(app* term);
+        virtual expr_ref mk_ge(filter_model_converter& fm, theory_var v, inf_numeral const& val);
+        void enable_record_conflict(expr* bound);
+        void record_conflict(unsigned num_lits, literal const * lits, 
+                          unsigned num_eqs, enode_pair const * eqs,
+                          unsigned num_params, parameter* params);
+        inf_eps_rational<inf_rational> conflict_minimize();
+    private:
+        virtual expr_ref mk_gt(theory_var v);
+
+        bool_var m_bound_watch;
+        inf_eps_rational<inf_rational> m_upper_bound;
+        bool get_theory_vars(expr * n, uint_set & vars);
+    public:
         // -----------------------------------
         //
         // Pretty Printing
@@ -1034,6 +1079,8 @@ namespace smt {
         void display_bounds_in_smtlib() const;
         void display_nl_monomials(std::ostream & out) const;
         void display_coeff_exprs(std::ostream & out, sbuffer<coeff_expr> const & p) const;
+        void display_interval(std::ostream& out, interval const& i);
+        void display_deps(std::ostream& out, v_dependency* dep);
 
     protected:
         // -----------------------------------
@@ -1049,9 +1096,11 @@ namespace smt {
         bool wf_rows() const;
         bool wf_column(theory_var v) const;
         bool wf_columns() const;
+        bool valid_assignment() const;
         bool valid_row_assignment() const;
         bool valid_row_assignment(row const & r) const;
         bool satisfy_bounds() const;
+        bool satisfy_integrality() const;
 #endif
     };
     
@@ -1071,6 +1120,7 @@ namespace smt {
         static inf_numeral mk_inf_numeral(numeral const & n, numeral const & r) {
             return inf_numeral(n, r);
         }
+        static bool is_infinite(inf_numeral const& ) { return false; }
         mi_ext() : m_int_epsilon(rational(1)), m_real_epsilon(rational(0), true) {}
     };
 
@@ -1087,6 +1137,8 @@ namespace smt {
             UNREACHABLE();
             return inf_numeral(n);
         }
+        static bool is_infinite(inf_numeral const& ) { return false; }
+
         i_ext() : m_int_epsilon(1), m_real_epsilon(1) {}
     };
 
@@ -1103,6 +1155,8 @@ namespace smt {
             UNREACHABLE();
             return inf_numeral(n);
         }
+        static bool is_infinite(inf_numeral const& ) { return false; }
+
         si_ext(): m_int_epsilon(s_integer(1)), m_real_epsilon(s_integer(1)) {}
     };
     
@@ -1123,13 +1177,41 @@ namespace smt {
         static inf_numeral mk_inf_numeral(numeral const& n, numeral const& i) {
             return inf_numeral(n, i);
         }
+        static bool is_infinite(inf_numeral const& ) { return false; }
+
         smi_ext() : m_int_epsilon(s_integer(1)), m_real_epsilon(s_integer(0), true) {}
     };
+
+    class inf_ext {
+    public:
+        typedef rational     numeral;
+        typedef inf_eps_rational<inf_rational>      inf_numeral;
+        inf_numeral   m_int_epsilon;
+        inf_numeral   m_real_epsilon;
+        numeral fractional_part(inf_numeral const& n) {
+            SASSERT(n.is_rational());
+            return n.get_rational() - floor(n);
+        }
+        static numeral fractional_part(numeral const & n) {
+            return n - floor(n);
+        }
+        static inf_numeral mk_inf_numeral(numeral const & n, numeral const & r) {
+            return inf_numeral(inf_rational(n, r));
+        }
+        static bool is_infinite(inf_numeral const& n) { 
+            return !n.get_infinity().is_zero(); 
+        }
+
+        inf_ext() : m_int_epsilon(inf_rational(rational(1))), m_real_epsilon(inf_rational(rational(0), true)) {}
+    };
+
     
     typedef theory_arith<mi_ext> theory_mi_arith;
     typedef theory_arith<i_ext> theory_i_arith;
+    typedef smt::theory_arith<inf_ext> theory_inf_arith;
     // typedef theory_arith<si_ext> theory_si_arith;
     // typedef theory_arith<smi_ext> theory_smi_arith;
+
     
 };
 

src/smt/theory_arith_core.h

@@ -821,7 +821,7 @@ namespace smt {
             m_new_atoms.push_back(a1);
             return;
         }
-        numeral const & k1(a1->get_k());
+        inf_numeral const & k1(a1->get_k());
         atom_kind kind1 = a1->get_atom_kind();
         TRACE("mk_bound_axioms", tout << "making bound axioms for v" << v << " " << kind1 << " " << k1 << "\n";);
         typename atoms::iterator it  = occs.begin();
@@ -831,7 +831,7 @@ namespace smt {
         typename atoms::iterator hi_inf = end, hi_sup = end;
         for (; it != end; ++it) {
             atom * a2 = *it;            
-            numeral const & k2(a2->get_k());
+            inf_numeral const & k2(a2->get_k());
             atom_kind kind2 = a2->get_atom_kind();
             SASSERT(k1 != k2 || kind1 != kind2);
             if (kind2 == A_LOWER) {
@@ -864,12 +864,13 @@ namespace smt {
         theory_var v = a1->get_var();
         literal   l1(a1->get_bool_var()); 
         literal   l2(a2->get_bool_var()); 
-        numeral const & k1(a1->get_k());
-        numeral const & k2(a2->get_k());
+        inf_numeral const & k1(a1->get_k());
+        inf_numeral const & k2(a2->get_k());
         atom_kind kind1 = a1->get_atom_kind();
         atom_kind kind2 = a2->get_atom_kind();
         bool v_is_int = is_int(v);
         SASSERT(v == a2->get_var());
+        if (k1 == k2 && kind1 == kind2) return;
         SASSERT(k1 != k2 || kind1 != kind2);
         parameter coeffs[3] = { parameter(symbol("farkas")), 
                                 parameter(rational(1)), parameter(rational(1)) };
@@ -890,7 +891,7 @@ namespace smt {
             else {
                 // k1 > hi_inf, k1 <= x => ~(x <= hi_inf)
                 mk_clause(~l1, ~l2, 3, coeffs);
-                if (v_is_int && k1 == k2 + numeral(1)) {
+                if (v_is_int && k1 == k2 + inf_numeral(1)) {
                     // k1 <= x or x <= k1-1
                     mk_clause(l1, l2, 3, coeffs);
                 }
@@ -904,7 +905,7 @@ namespace smt {
             else {
                 // k1 < k2, k2 <= x => ~(x <= k1)
                 mk_clause(~l1, ~l2, 3, coeffs); 
-                if (v_is_int && k1 == k2 - numeral(1)) {
+                if (v_is_int && k1 == k2 - inf_numeral(1)) {
                     // x <= k1 or k1+l <= x
                     mk_clause(l1, l2, 3, coeffs);
                 }
@@ -1000,14 +1001,14 @@ namespace smt {
         typename atoms::iterator it, 
         typename atoms::iterator end,
         bool& found_compatible) {
-        numeral const & k1(a1->get_k());
+        inf_numeral const & k1(a1->get_k());
         typename atoms::iterator result = end;
         found_compatible = false;
         for (; it != end; ++it) {
             atom * a2 = *it;            
             if (a1 == a2) continue;
             if (a2->get_atom_kind() != kind) continue;
-            numeral const & k2(a2->get_k());
+            inf_numeral const & k2(a2->get_k());
             found_compatible = true;
             if (k2 <= k1) {
                 result = it;
@@ -1027,13 +1028,13 @@ namespace smt {
         typename atoms::iterator it, 
         typename atoms::iterator end,
         bool& found_compatible) {
-        numeral const & k1(a1->get_k());
+        inf_numeral const & k1(a1->get_k());
         found_compatible = false;
         for (; it != end; ++it) {
             atom * a2 = *it;            
             if (a1 == a2) continue;
             if (a2->get_atom_kind() != kind) continue;
-            numeral const & k2(a2->get_k());
+            inf_numeral const & k2(a2->get_k());
             found_compatible = true;
             if (k1 < k2) {
                 return it;
@@ -1081,7 +1082,7 @@ namespace smt {
         ctx.set_var_theory(bv, get_id());
         rational _k;
         VERIFY(m_util.is_numeral(rhs, _k));
-        numeral   k(_k);
+        inf_numeral   k(_k);
         atom * a = alloc(atom, bv, v, k, kind);
         mk_bound_axioms(a);
         m_unassigned_atoms[v]++;
@@ -1269,7 +1270,9 @@ namespace smt {
                 result = FC_GIVEUP;
                 break;
             case FC_CONTINUE:
-                TRACE("final_check_arith", tout << "continue arith...\n";);
+                TRACE("final_check_arith", 
+                      tout << "continue arith..." 
+                      << (get_context().inconsistent()?"inconsistent\n":"\n"););
                 return FC_CONTINUE;
             }
         }
@@ -1519,7 +1522,8 @@ namespace smt {
         m_assume_eq_head(0),
         m_nl_rounds(0),
         m_nl_gb_exhausted(false),
-        m_nl_new_exprs(m) {
+        m_nl_new_exprs(m),
+        m_bound_watch(null_bool_var) {
     }
 
     template<typename Ext>
@@ -2184,7 +2188,7 @@ namespace smt {
               tout << "is_below_lower: " << below_lower(x_i) << ", is_above_upper: " << above_upper(x_i) << "\n";);
         antecedents& ante = get_antecedents();
         explain_bound(r, idx, !is_below, delta, ante);
-        b->push_justification(ante, numeral(1), proofs_enabled());
+        b->push_justification(ante, numeral(1), coeffs_enabled());
        
 
         set_conflict(ante.lits().size(), ante.lits().c_ptr(), 
@@ -2327,11 +2331,12 @@ namespace smt {
     void theory_arith<Ext>::sign_bound_conflict(bound * b1, bound * b2) {
         SASSERT(b1->get_var() == b2->get_var());
         antecedents& ante = get_antecedents();
-        b1->push_justification(ante, numeral(1), proofs_enabled());
-        b2->push_justification(ante, numeral(1), proofs_enabled());
+        b1->push_justification(ante, numeral(1), coeffs_enabled());
+        b2->push_justification(ante, numeral(1), coeffs_enabled());
 
         set_conflict(ante.lits().size(), ante.lits().c_ptr(), ante.eqs().size(), ante.eqs().c_ptr(), ante, is_int(b1->get_var()), "farkas");
-        TRACE("arith_conflict", tout << "bound conflict\n";);
+        TRACE("arith_conflict", tout << "bound conflict v" << b1->get_var() << "\n";
+              tout << "bounds: " << b1 << " " << b2 << "\n";);
     }
 
     // -----------------------------------
@@ -2587,7 +2592,7 @@ namespace smt {
                 if (!b->has_justification())
                     continue;
                 if (!relax_bounds() || delta.is_zero()) {
-                    b->push_justification(ante, it->m_coeff, proofs_enabled());
+                    b->push_justification(ante, it->m_coeff, coeffs_enabled());
                     continue;
                 }
                 numeral coeff = it->m_coeff;
@@ -2649,7 +2654,7 @@ namespace smt {
                 SASSERT(!is_b_lower || k_2 <= k_1);
                 SASSERT(is_b_lower  || k_2 >= k_1);
                 if (new_atom == 0) {
-                    b->push_justification(ante, coeff, proofs_enabled());
+                    b->push_justification(ante, coeff, coeffs_enabled());
                     continue;
                 }
                 SASSERT(!is_b_lower || k_2 < k_1);
@@ -2663,7 +2668,7 @@ namespace smt {
                     delta -= coeff*(k_2 - k_1);
                 }
                 TRACE("propagate_bounds", tout << "delta (after replace): " << delta << "\n";);
-                new_atom->push_justification(ante, coeff, proofs_enabled());
+                new_atom->push_justification(ante, coeff, coeffs_enabled());
                 SASSERT(delta >= inf_numeral::zero());
             }
         }
@@ -2682,7 +2687,7 @@ namespace smt {
             bool_var bv = a->get_bool_var();
             literal  l(bv);
             if (get_context().get_assignment(bv) == l_undef) {
-                numeral const & k2 = a->get_k();
+                inf_numeral const & k2 = a->get_k();
                 delta.reset();
                 if (a->get_atom_kind() == A_LOWER) {
                     // v >= k  k >= k2  |-  v >= k2
@@ -2863,13 +2868,13 @@ namespace smt {
               for (unsigned i = 0; i < num_literals; i++) {
                   ctx.display_detailed_literal(tout, lits[i]);
                   tout << " ";
-                  if (proofs_enabled()) {
+                  if (coeffs_enabled()) {
                       tout << "bound: " << bounds.lit_coeffs()[i] << "\n";
                   }
               }
               for (unsigned i = 0; i < num_eqs; i++) {
                   tout << "#" << eqs[i].first->get_owner_id() << "=#" << eqs[i].second->get_owner_id() << " ";
-                  if (proofs_enabled()) {
+                  if (coeffs_enabled()) {
                       tout << "bound: " << bounds.eq_coeffs()[i] << "\n";
                   }
               }
@@ -2877,6 +2882,7 @@ namespace smt {
                   tout << bounds.params(proof_rule)[i] << "\n";
               }
               tout << "\n";);
+        record_conflict(num_literals, lits, num_eqs, eqs, bounds.num_params(), bounds.params(proof_rule));
         ctx.set_conflict(
             ctx.mk_justification(
                 ext_theory_conflict_justification(get_id(), r, num_literals, lits, num_eqs, eqs, 
@@ -2893,8 +2899,8 @@ namespace smt {
         typename vector<row_entry>::const_iterator end = r.end_entries();
         for (; it != end; ++it) {
             if (!it->is_dead() && is_fixed(it->m_var)) {
-                lower(it->m_var)->push_justification(antecedents, it->m_coeff, proofs_enabled());
-                upper(it->m_var)->push_justification(antecedents, it->m_coeff, proofs_enabled());                
+                lower(it->m_var)->push_justification(antecedents, it->m_coeff, coeffs_enabled());
+                upper(it->m_var)->push_justification(antecedents, it->m_coeff, coeffs_enabled());                
             }
         }
     }
@@ -2997,6 +3003,9 @@ namespace smt {
                 if (!get_context().is_shared(get_enode(v)))
                     continue;
                 inf_numeral const & val = get_value(v);
+                if (Ext::is_infinite(val)) {
+                    continue;
+                }
                 rational value = val.get_rational().to_rational() + m_epsilon.to_rational() * val.get_infinitesimal().to_rational();
                 theory_var v2;
                 if (mapping.find(value, v2)) {
@@ -3220,11 +3229,13 @@ namespace smt {
                 case QUASI_BASE:
                     SASSERT(m_columns[v].size() == 1);
                     del_row(get_var_row(v));
+                    TRACE("arith_make_feasible", tout << "del row v" << v << "\n";);
                     break; 
                 case BASE:
                     SASSERT(lazy_pivoting_lvl() != 0 || m_columns[v].size() == 1);
                     if (lazy_pivoting_lvl() > 0)
                         eliminate<false>(v, false);
+                    TRACE("arith_make_feasible", tout << "del row v" << v << "\n";);
                     del_row(get_var_row(v));
                     break;
                 case NON_BASE: {
@@ -3236,6 +3247,10 @@ namespace smt {
                         pivot<false>(r.get_base_var(), v, r[entry->m_row_idx].m_coeff, false);
                         SASSERT(is_base(v));
                         del_row(get_var_row(v));
+                        TRACE("arith_make_feasible", tout << "del row v" << v << "\n";);
+                    }
+                    else {
+                        TRACE("arith_make_feasible", tout << "no row v" << v << "\n";);
                     }
                     break;
                 } }

src/smt/theory_arith_int.h

@@ -33,11 +33,15 @@ namespace smt {
     // Integrality
     //
     // -----------------------------------
+       
     
     /**
        \brief Move non base variables to one of its bounds.
        If the variable does not have bounds, it is integer, but it is not assigned to an integer value, then the 
        variable is set to an integer value.
+       In mixed integer/real problems moving a real variable to a bound could cause an integer value to 
+       have an infinitesimal. Such an assignment would disable mk_gomory_cut, and Z3 would loop.
+       
     */
     template<typename Ext>
     void theory_arith<Ext>::move_non_base_vars_to_bounds() {
@@ -413,10 +417,10 @@ namespace smt {
         for (; it != end; ++it) {
             // All non base variables must be at their bounds and assigned to rationals (that is, infinitesimals are not allowed).
             if (!it->is_dead() && it->m_var != b && (!at_bound(it->m_var) || !get_value(it->m_var).is_rational())) {
-                TRACE("gomory_cut", tout << "row is gomory cut target:\n";
+                TRACE("gomory_cut", tout << "row is not gomory cut target:\n";
                       display_var(tout, it->m_var);
                       tout << "at_bound:      " << at_bound(it->m_var) << "\n";
-                      tout << "infinitesimal: " << get_value(it->m_var).is_rational() << "\n";);
+                      tout << "infinitesimal: " << !get_value(it->m_var).is_rational() << "\n";);
                 return false;
             }
         }
@@ -473,7 +477,8 @@ namespace smt {
     */
     template<typename Ext>
     bool theory_arith<Ext>::mk_gomory_cut(row const & r) {
-        SASSERT(!all_coeff_int(r));
+        // The following assertion is wrong. It may be violated in mixed-integer problems.
+        // SASSERT(!all_coeff_int(r));
         theory_var x_i = r.get_base_var();
         
         SASSERT(is_int(x_i));
@@ -525,7 +530,7 @@ namespace smt {
                         }
                         // k += new_a_ij * lower_bound(x_j).get_rational();
                         k.addmul(new_a_ij, lower_bound(x_j).get_rational());
-                        lower(x_j)->push_justification(ante, numeral::zero(), proofs_enabled());
+                        lower(x_j)->push_justification(ante, numeral::zero(), coeffs_enabled());
                     }
                     else {
                         SASSERT(at_upper(x_j));
@@ -541,7 +546,7 @@ namespace smt {
                         }
                         // k += new_a_ij * upper_bound(x_j).get_rational();
                         k.addmul(new_a_ij, upper_bound(x_j).get_rational());
-                        upper(x_j)->push_justification(ante, numeral::zero(), proofs_enabled());
+                        upper(x_j)->push_justification(ante, numeral::zero(), coeffs_enabled());
                     }
                     pol.push_back(row_entry(new_a_ij, x_j));
                 }
@@ -566,7 +571,7 @@ namespace smt {
                             }
                             // k += new_a_ij * lower_bound(x_j).get_rational();
                             k.addmul(new_a_ij, lower_bound(x_j).get_rational());
-                            lower(x_j)->push_justification(ante, numeral::zero(), proofs_enabled());
+                            lower(x_j)->push_justification(ante, numeral::zero(), coeffs_enabled());
                         }
                         else {
                             SASSERT(at_upper(x_j));
@@ -579,7 +584,7 @@ namespace smt {
                             new_a_ij.neg(); // the upper terms are inverted
                             // k += new_a_ij * upper_bound(x_j).get_rational();
                             k.addmul(new_a_ij, upper_bound(x_j).get_rational());
-                            upper(x_j)->push_justification(ante, numeral::zero(), proofs_enabled());
+                            upper(x_j)->push_justification(ante, numeral::zero(), coeffs_enabled());
                         }
                         TRACE("gomory_cut_detail", tout << "new_a_ij: " << new_a_ij << "\n";);
                         pol.push_back(row_entry(new_a_ij, x_j));
@@ -772,8 +777,8 @@ namespace smt {
                         // u += ncoeff * lower_bound(v).get_rational();
                         u.addmul(ncoeff, lower_bound(v).get_rational());
                     }
-                    lower(v)->push_justification(ante, numeral::zero(), proofs_enabled());
-                    upper(v)->push_justification(ante, numeral::zero(), proofs_enabled());
+                    lower(v)->push_justification(ante, numeral::zero(), coeffs_enabled());
+                    upper(v)->push_justification(ante, numeral::zero(), coeffs_enabled());
                 }
                 else if (gcds.is_zero()) {
                     gcds = abs_ncoeff; 
@@ -1379,6 +1384,7 @@ namespace smt {
         m_branch_cut_counter++;
         // TODO: add giveup code
         if (m_branch_cut_counter % m_params.m_arith_branch_cut_ratio == 0) {
+            TRACE("opt", display(tout););
             move_non_base_vars_to_bounds();
             if (!make_feasible()) {
                 TRACE("arith_int", tout << "failed to move variables to bounds.\n";);
@@ -1390,7 +1396,9 @@ namespace smt {
                 TRACE("arith_int", tout << "v" << int_var << " does not have an integer assignment: " << get_value(int_var) << "\n";);
                 SASSERT(is_base(int_var));
                 row const & r = m_rows[get_var_row(int_var)];
-                mk_gomory_cut(r);
+                if (!mk_gomory_cut(r)) {
+                    // silent failure
+                }
                 return FC_CONTINUE;
             }
         }
@@ -1400,7 +1408,7 @@ namespace smt {
             }
             theory_var int_var = find_infeasible_int_base_var();
             if (int_var != null_theory_var) {
-                TRACE("arith_int", tout << "v" << int_var << " does not have and integer assignment: " << get_value(int_var) << "\n";);
+                TRACE("arith_int", tout << "v" << int_var << " does not have an integer assignment: " << get_value(int_var) << "\n";);
                 // apply branching 
                 branch_infeasible_int_var(int_var);
                 return FC_CONTINUE;

src/smt/theory_arith_inv.h

@@ -199,6 +199,27 @@ namespace smt {
         }
         return true;
     }
+
+    template<typename Ext>
+    bool theory_arith<Ext>::satisfy_integrality() const {
+        int num = get_num_vars();
+        for (theory_var v = 0; v < num; v++) {
+            if (is_int(v) && !get_value(v).is_int()) {
+                TRACE("bound_bug", display_var(tout, v); display(tout););
+                return false;
+            }
+        }
+        return true;
+    }
+
+    template<typename Ext>
+    bool theory_arith<Ext>::valid_assignment() const {
+        return 
+            valid_row_assignment() &&
+            satisfy_bounds() &&
+            satisfy_integrality();
+    }
+
 #endif
 
 };

src/smt/theory_arith_nl.h

@@ -333,12 +333,15 @@ namespace smt {
     void theory_arith<Ext>::mul_bound_of(expr * var, unsigned power, interval & target) {
         theory_var v  = expr2var(var);
         interval i = mk_interval_for(v);
-        TRACE("non_linear", tout << "bound: " << i << "\n" << mk_pp(var, get_manager()) << "\n";
+        
+        TRACE("non_linear", 
+              display_interval(tout << "bound: ",i); tout << i << "\n";
+              tout << mk_pp(var, get_manager()) << "\n";
               tout << "power " << power << ": " << expt(i, power) << "\n";
-              tout << "target before: " << target << "\n";); 
+              display_interval(tout << "target before: ", target); tout << "\n";); 
         i.expt(power);
         target *= i;
-        TRACE("non_linear", tout << "target after: " << target << "\n";);
+        TRACE("non_linear", display_interval(tout << "target after: ", target); tout << "\n";);
     }
 
     /**
@@ -427,12 +430,12 @@ namespace smt {
     template<typename Ext>
     void theory_arith<Ext>::mk_derived_nl_bound(theory_var v, inf_numeral const & coeff, bound_kind k, v_dependency * dep) {
         inf_numeral coeff_norm = normalize_bound(v, coeff, k);
-        TRACE("buggy_bound", tout << "v" << v << " " << coeff << " " << coeff_norm << " " << k << "\n";);
         derived_bound * new_bound = alloc(derived_bound, v, coeff_norm, k);
         m_bounds_to_delete.push_back(new_bound);
         m_asserted_bounds.push_back(new_bound);
         // copy justification to new bound
         dependency2new_bound(dep, *new_bound);
+        TRACE("buggy_bound", new_bound->display(*this, tout); tout << "\n";);
     }
 
     /**
@@ -449,7 +452,8 @@ namespace smt {
                 new_lower += get_epsilon(v);
             bound * old_lower = lower(v);
             if (old_lower == 0 || new_lower > old_lower->get_value()) {
-                TRACE("non_linear", tout << "NEW lower bound for v" << v << " " << new_lower << "\n";);
+                TRACE("non_linear", tout << "NEW lower bound for v" << v << " " << new_lower << "\n";
+                      display_interval(tout, i); tout << "\n";);
                 mk_derived_nl_bound(v, new_lower, B_LOWER, i.get_lower_dependencies());
                 r = true;
             }
@@ -460,7 +464,8 @@ namespace smt {
                 new_upper -= get_epsilon(v);
             bound * old_upper = upper(v);
             if (old_upper == 0 || new_upper < old_upper->get_value()) {
-                TRACE("non_linear", tout << "NEW upper bound for v" << v << " " << new_upper << "\n";);
+                TRACE("non_linear", tout << "NEW upper bound for v" << v << " " << new_upper << "\n";
+                      display_interval(tout, i); tout << "\n";);
                 mk_derived_nl_bound(v, new_upper, B_UPPER, i.get_upper_dependencies());
                 r = true;
             }
@@ -897,31 +902,54 @@ namespace smt {
         m_tmp_lit_set.reset();
         m_tmp_eq_set.reset();
 
-        bool found_zero = false;
         SASSERT(is_pure_monomial(m));
-        for (unsigned i = 0; i < to_app(m)->get_num_args(); i++) {
-            expr * arg = to_app(m)->get_arg(i);
-            if (!found_zero) {                                                                  
-                theory_var _var = expr2var(arg);                                                
-                if (is_fixed(_var)) {                                                           
-                    bound * l  = lower(_var);                                                   
-                    bound * u  = upper(_var);                                                   
-                    if (l->get_value().is_zero()) {                                             
-                        /* if zero was found, then it is the explanation */                     
-                        SASSERT(k.is_zero());                                                   
-                        found_zero = true;                                                      
-                        m_tmp_lit_set.reset();                                                  
-                        m_tmp_eq_set.reset();                                                   
-                        new_lower->m_lits.reset();                                              
-                        new_lower->m_eqs.reset();                                               
-                    }                                                                           
-                    accumulate_justification(*l, *new_lower, numeral::zero(), m_tmp_lit_set, m_tmp_eq_set);      
-                    accumulate_justification(*u, *new_lower, numeral::zero(), m_tmp_lit_set, m_tmp_eq_set);      
-                }                                                                               
-            }                                                                                   
-        }
+        bool found_zero = false;
+        for (unsigned i = 0; !found_zero && i < to_app(m)->get_num_args(); i++) {            
+            expr * arg = to_app(m)->get_arg(i);            
+            theory_var _var = expr2var(arg);                                                
+            if (is_fixed(_var)) {                                                           
+                bound * l  = lower(_var);                                                   
+                bound * u  = upper(_var);                                                   
+                if (l->get_value().is_zero()) {                                             
+                    /* if zero was found, then it is the explanation */                     
+                    SASSERT(k.is_zero());                                                   
+                    found_zero = true;                                                      
+                    m_tmp_lit_set.reset();                                                  
+                    m_tmp_eq_set.reset();                                                   
+                    new_lower->m_lits.reset();                                              
+                    new_lower->m_eqs.reset();                                               
+                }                                                                           
+                accumulate_justification(*l, *new_lower, numeral::zero(), m_tmp_lit_set, m_tmp_eq_set);
+                
+                TRACE("non_linear", 
+                      for (unsigned j = 0; j < new_lower->m_lits.size(); ++j) {
+                          ctx.display_detailed_literal(tout, new_lower->m_lits[j]);
+                          tout << " ";
+                      }
+                      tout << "\n";);
+      
+                accumulate_justification(*u, *new_lower, numeral::zero(), m_tmp_lit_set, m_tmp_eq_set);
+      
+                TRACE("non_linear", 
+                      for (unsigned j = 0; j < new_lower->m_lits.size(); ++j) {
+                          ctx.display_detailed_literal(tout, new_lower->m_lits[j]);
+                          tout << " ";
+                      }
+                      tout << "\n";);
+                    
+            }                                                                               
+        } 
         new_upper->m_lits.append(new_lower->m_lits);
         new_upper->m_eqs.append(new_lower->m_eqs);
+
+        TRACE("non_linear", 
+              tout << "lower: " << new_lower << " upper: " << new_upper << "\n";
+              for (unsigned j = 0; j < new_upper->m_lits.size(); ++j) {
+                  ctx.display_detailed_literal(tout, new_upper->m_lits[j]);
+                  tout << " ";
+              }
+              tout << "\n";);
+
         return true;
     }
 
@@ -1887,6 +1915,10 @@ namespace smt {
         derived_bound  b(null_theory_var, inf_numeral(0), B_LOWER);
         dependency2new_bound(d, b);
         set_conflict(b.m_lits.size(), b.m_lits.c_ptr(), b.m_eqs.size(), b.m_eqs.c_ptr(), ante, is_lia, "arith_nl");
+        TRACE("non_linear", 
+              for (unsigned i = 0; i < b.m_lits.size(); ++i) {
+                  tout << b.m_lits[i] << " ";
+              });
     }
 
     /**

src/smt/theory_arith_pp.h

@@ -395,9 +395,30 @@ namespace smt {
     template<typename Ext>
     void theory_arith<Ext>::display_bound(std::ostream & out, bound * b, unsigned indent) const {
         for (unsigned i = 0; i < indent; i++) out << "  ";
-        theory_var v = b->get_var();
-        enode * e    = get_enode(v);
-        out << "v" << v << " #" << e->get_owner_id() << " " << (b->get_bound_kind() == B_LOWER ? ">=" : "<=") << " " << b->get_value() << "\n";
+        b->display(*this, out);
+        out << "\n";
+    }
+
+    template<typename Ext>
+    void theory_arith<Ext>::display_deps(std::ostream & out, v_dependency* dep) {
+        ptr_vector<void> bounds;
+        m_dep_manager.linearize(dep, bounds);
+        m_tmp_lit_set.reset();
+        m_tmp_eq_set.reset();
+        ptr_vector<void>::const_iterator it  = bounds.begin();
+        ptr_vector<void>::const_iterator end = bounds.end();
+        for (; it != end; ++it) {
+            bound * b = static_cast<bound*>(*it);
+            out << " ";
+            b->display(*this, out);
+        }
+    }
+
+    template<typename Ext>
+    void theory_arith<Ext>::display_interval(std::ostream & out, interval const& i) {
+        i.display(out);
+        display_deps(out << " lo:", i.get_lower_dependencies());
+        display_deps(out << " hi:", i.get_upper_dependencies());
     }
 
     template<typename Ext>
@@ -428,7 +449,7 @@ namespace smt {
     template<typename Ext>
     void theory_arith<Ext>::display_atom(std::ostream & out, atom * a, bool show_sign) const {
         theory_var      v = a->get_var();
-        numeral const & k = a->get_k();
+        inf_numeral const & k = a->get_k();
         enode *         e = get_enode(v);
         if (show_sign) {
             if (!a->is_true()) 

src/smt/theory_array_base.cpp

@@ -362,7 +362,6 @@ namespace smt {
          (store A i v) <--- v is used as an value
     */
     bool theory_array_base::is_shared(theory_var v) const {
-        context & ctx  = get_context();
         enode * n      = get_enode(v);
         enode * r      = n->get_root();
         bool is_array  = false;

src/smt/theory_array_full.h

@@ -19,8 +19,9 @@ Revision History:
 #ifndef _THEORY_ARRAY_FULL_H_
 #define _THEORY_ARRAY_FULL_H_
 
-#include"theory_array.h"
+#include "theory_array.h"
 #include "simplifier.h"
+#include "ast_trail.h"
 
 namespace smt {
 

src/smt/theory_datatype.cpp

@@ -61,6 +61,13 @@ namespace smt {
             if (antecedent == null_literal) {
                 ctx.assign_eq(lhs, ctx.get_enode(rhs), eq_justification::mk_axiom());
             }
+            else if (ctx.get_assignment(antecedent) != l_true) {
+                literal l(mk_eq(lhs->get_owner(), rhs, true));
+                ctx.mark_as_relevant(l);
+                ctx.mark_as_relevant(antecedent);
+                literal lits[2] = {l, ~antecedent};
+                ctx.mk_th_axiom(get_id(), 2, lits);
+            }
             else {
                 SASSERT(ctx.get_assignment(antecedent) == l_true);
                 region & r   = ctx.get_region();
@@ -143,6 +150,48 @@ namespace smt {
         ctx.set_conflict(ctx.mk_justification(ext_theory_conflict_justification(get_id(), reg, 1, &l, 1, &p)));
     }
 
+    /**
+       \brief Given a field update n := { r with field := v } for constructor C, assert the axioms:
+       (=> (is-C r) (= (acc_j n) (acc_j r))) for acc_j != field
+       (=> (is-C r) (= (field n) v))         for acc_j != field
+       (=> (not (is-C r)) (= n r))
+    */
+    void theory_datatype::assert_update_field_axioms(enode * n) {
+        m_stats.m_assert_update_field++;
+        SASSERT(is_update_field(n));
+        context & ctx = get_context();
+        ast_manager & m  = get_manager();
+        app*        own  = n->get_owner();
+        expr*       arg1 = own->get_arg(0);
+        expr*       arg2 = own->get_arg(1);
+        func_decl * upd  = n->get_decl();
+        func_decl * acc  = to_func_decl(upd->get_parameter(0).get_ast());
+        func_decl * con  = m_util.get_accessor_constructor(acc);
+        func_decl * rec  = m_util.get_constructor_recognizer(con);
+        ptr_vector<func_decl> const * accessors   = m_util.get_constructor_accessors(con);
+        ptr_vector<func_decl>::const_iterator it  = accessors->begin();
+        ptr_vector<func_decl>::const_iterator end = accessors->end();
+        app_ref rec_app(m.mk_app(rec, arg1), m);
+        ctx.internalize(rec_app, false);
+        literal is_con(ctx.get_bool_var(rec_app));
+        for (; it != end; ++it) {
+            enode* arg;
+            func_decl * acc1   = *it;
+            if (acc1 == acc) {
+                arg = n->get_arg(1);
+            }
+            else {
+                app* acc_app = m.mk_app(acc1, arg1);
+                ctx.internalize(acc_app, false);
+                arg = ctx.get_enode(acc_app);
+            }
+            app * acc_own = m.mk_app(acc1, own);
+            assert_eq_axiom(arg, acc_own, is_con); 
+        }
+        // update_field is identity if 'n' is not created by a matching constructor.        
+        assert_eq_axiom(n, arg1, ~is_con);
+    }
+
     theory_var theory_datatype::mk_var(enode * n) {
         theory_var r  = theory::mk_var(n);
         theory_var r2 = m_find.mk_var();
@@ -150,15 +199,17 @@ namespace smt {
         SASSERT(r == static_cast<int>(m_var_data.size()));
         m_var_data.push_back(alloc(var_data));
         var_data * d  = m_var_data[r];
+        context & ctx   = get_context();
+        ctx.attach_th_var(n, this, r);
         if (is_constructor(n)) {
             d->m_constructor = n;
-            get_context().attach_th_var(n, this, r);
             assert_accessor_axioms(n);
         }
+        else if (is_update_field(n)) {
+            assert_update_field_axioms(n);
+        }
         else {
             ast_manager & m = get_manager();
-            context & ctx   = get_context();
-            ctx.attach_th_var(n, this, r);
             sort * s      = m.get_sort(n->get_owner());
             if (m_util.get_datatype_num_constructors(s) == 1) {
                 func_decl * c = m_util.get_datatype_constructors(s)->get(0);
@@ -192,7 +243,7 @@ namespace smt {
             ctx.set_var_theory(bv, get_id());
             ctx.set_enode_flag(bv, true);
         }
-        if (is_constructor(term)) {
+        if (is_constructor(term) || is_update_field(term)) {
             SASSERT(!is_attached_to_var(e));
             // *** We must create a theory variable for each argument that has sort datatype ***
             //
@@ -478,6 +529,7 @@ namespace smt {
         st.update("datatype splits", m_stats.m_splits);
         st.update("datatype constructor ax", m_stats.m_assert_cnstr);
         st.update("datatype accessor ax", m_stats.m_assert_accessor);
+        st.update("datatype update ax", m_stats.m_assert_update_field);
     }
     
     void theory_datatype::display_var(std::ostream & out, theory_var v) const {

src/smt/theory_datatype.h

@@ -41,7 +41,7 @@ namespace smt {
 
         struct stats {
             unsigned   m_occurs_check, m_splits;
-            unsigned   m_assert_cnstr, m_assert_accessor;
+            unsigned   m_assert_cnstr, m_assert_accessor, m_assert_update_field;
             void reset() { memset(this, 0, sizeof(stats)); }
             stats() { reset(); }
         };
@@ -58,14 +58,17 @@ namespace smt {
         bool is_constructor(app * f) const { return m_util.is_constructor(f); }
         bool is_recognizer(app * f) const { return m_util.is_recognizer(f); }
         bool is_accessor(app * f) const { return m_util.is_accessor(f); }
-        
+        bool is_update_field(app * f) const { return m_util.is_update_field(f); }
+
         bool is_constructor(enode * n) const { return is_constructor(n->get_owner()); }
         bool is_recognizer(enode * n) const { return is_recognizer(n->get_owner()); }
         bool is_accessor(enode * n) const { return is_accessor(n->get_owner()); }
+        bool is_update_field(enode * n) const { return m_util.is_update_field(n->get_owner()); }
 
         void assert_eq_axiom(enode * lhs, expr * rhs, literal antecedent);
         void assert_is_constructor_axiom(enode * n, func_decl * c, literal antecedent);
         void assert_accessor_axioms(enode * n);
+        void assert_update_field_axioms(enode * n);
         void add_recognizer(theory_var v, enode * recognizer);
         void propagate_recognizer(theory_var v, enode * r);
         void sign_recognizer_conflict(enode * c, enode * r);

src/smt/theory_dense_diff_logic.h

@@ -25,6 +25,8 @@ TODO: eager equality propagation
 #include"theory_arith_params.h"
 #include"arith_decl_plugin.h"
 #include"arith_eq_adapter.h"
+#include"theory_opt.h"
+
 
 namespace smt {
 
@@ -41,7 +43,7 @@ namespace smt {
     };
     
     template<typename Ext>
-    class theory_dense_diff_logic : public theory, private Ext {
+    class theory_dense_diff_logic : public theory, public theory_opt, private Ext {
     public:
         theory_dense_diff_logic_statistics      m_stats;
 
@@ -127,6 +129,12 @@ namespace smt {
         svector<scope>        m_scopes;
         bool                  m_non_diff_logic_exprs;
 
+        // For optimization purpose
+        typedef vector <std::pair<theory_var, rational> > objective_term;
+        vector<objective_term>         m_objectives;
+        vector<rational>               m_objective_consts;
+        vector<expr_ref_vector>        m_objective_assignments;
+
         struct f_target {
             theory_var m_target;
             numeral    m_new_distance;
@@ -189,6 +197,8 @@ namespace smt {
         void del_atoms(unsigned old_size);
         void del_vars(unsigned old_num_vars);
         void init_model();
+        bool internalize_objective(expr * n, rational const& m, rational& r, objective_term & objective);
+        expr_ref mk_ineq(theory_var v, inf_rational const& val, bool is_strict);
 #ifdef Z3DEBUG
         bool check_vector_sizes() const;
         bool check_matrix() const;
@@ -229,6 +239,8 @@ namespace smt {
         virtual void flush_eh();
         virtual void reset_eh();
         
+        bool dump_lemmas() const { return m_params.m_arith_dump_lemmas; }
+
         virtual bool validate_eq_in_model(theory_var v1, theory_var v2, bool is_true) const;
 
         virtual void display(std::ostream & out) const;
@@ -249,6 +261,18 @@ namespace smt {
         virtual void init_model(model_generator & m);
         virtual model_value_proc * mk_value(enode * n, model_generator & mg);
 
+        // -----------------------------------
+        //
+        // Optimization
+        //
+        // -----------------------------------
+
+        virtual inf_eps_rational<inf_rational> maximize(theory_var v, expr_ref& blocker, bool& has_shared);
+        virtual inf_eps_rational<inf_rational> value(theory_var v);
+        virtual theory_var add_objective(app* term);
+        virtual expr_ref mk_gt(theory_var v, inf_rational const& val);
+        virtual expr_ref mk_ge(filter_model_converter& fm, theory_var v, inf_rational const& val);
+        
         // -----------------------------------
         //
         // Main

src/smt/theory_dense_diff_logic_def.h

@@ -23,6 +23,8 @@ Revision History:
 #include"theory_dense_diff_logic.h"
 #include"ast_pp.h"
 #include"smt_model_generator.h"
+#include"simplex.h"
+#include"simplex_def.h"
 
 namespace smt {
 
@@ -120,6 +122,7 @@ namespace smt {
         if (!m_non_diff_logic_exprs) {
             TRACE("non_diff_logic", tout << "found non diff logic expression:\n" << mk_pp(n, get_manager()) << "\n";);
             get_context().push_trail(value_trail<context, bool>(m_non_diff_logic_exprs));
+        IF_VERBOSE(0, verbose_stream() << "(smt.diff_logic: non-diff logic expression " << mk_pp(n, get_manager()) << ")\n";); 
             m_non_diff_logic_exprs = true;
         }
     }
@@ -587,6 +590,11 @@ namespace smt {
             context & ctx = get_context();
             region & r    = ctx.get_region();
             ctx.set_conflict(ctx.mk_justification(theory_conflict_justification(get_id(), r, antecedents.size(), antecedents.c_ptr())));
+
+            if (dump_lemmas()) {
+                ctx.display_lemma_as_smt_problem(antecedents.size(), antecedents.c_ptr(), false_literal, "");
+            }
+
             return;
         }
         
@@ -818,6 +826,273 @@ namespace smt {
         return alloc(expr_wrapper_proc, m_factory->mk_value(num, is_int(v)));
     }
 
+    // TBD: code is common to both sparse and dense difference logic solvers.
+    template<typename Ext>
+    bool theory_dense_diff_logic<Ext>::internalize_objective(expr * n, rational const& m, rational& q, objective_term & objective) {
+
+        // Compile term into objective_term format
+        rational r;
+        expr* x, *y;
+        if (m_autil.is_numeral(n, r)) {
+            q += r;
+        }
+        else if (m_autil.is_add(n)) {
+            for (unsigned i = 0; i < to_app(n)->get_num_args(); ++i) {
+                if (!internalize_objective(to_app(n)->get_arg(i), m, q, objective)) {
+                    return false;
+                }
+            }
+        }
+        else if (m_autil.is_mul(n, x, y) && m_autil.is_numeral(x, r)) {
+            return internalize_objective(y, m*r, q, objective);
+        }
+        else if (m_autil.is_mul(n, y, x) && m_autil.is_numeral(x, r)) {
+            return internalize_objective(y, m*r, q, objective);
+        }
+        else if (!is_app(n)) {
+            return false;
+        }
+        else if (to_app(n)->get_family_id() == m_autil.get_family_id()) {
+            return false;
+        }
+        else {
+            context& ctx = get_context();
+            enode * e = 0;
+            theory_var v = 0;
+            if (ctx.e_internalized(n)) {
+                e = ctx.get_enode(to_app(n));                
+            }
+            else {
+                e = ctx.mk_enode(to_app(n), false, false, true);
+            }            
+            v = e->get_th_var(get_id());
+            if (v == null_theory_var) {
+                v = mk_var(e);                
+            }
+
+            objective.push_back(std::make_pair(v, m));
+        }
+        return true;
+    }
+
+    template<typename Ext>
+    inf_eps_rational<inf_rational> theory_dense_diff_logic<Ext>::value(theory_var v) {
+        objective_term const& objective = m_objectives[v];   
+        inf_eps r = inf_eps(m_objective_consts[v]);
+        for (unsigned i = 0; i < objective.size(); ++i) {
+            numeral n = m_assignment[v];
+            rational r1 = n.get_rational().to_rational();
+            rational r2 = n.get_infinitesimal().to_rational();
+            r += objective[i].second * inf_eps(rational(0), inf_rational(r1, r2));
+        }
+        return r;
+    }
+
+    template<typename Ext>
+    inf_eps_rational<inf_rational> theory_dense_diff_logic<Ext>::maximize(theory_var v, expr_ref& blocker, bool& has_shared) {
+        typedef simplex::simplex<simplex::mpq_ext> Simplex;
+        Simplex S;
+        ast_manager& m = get_manager();
+        objective_term const& objective = m_objectives[v];
+        has_shared = false;
+        
+        IF_VERBOSE(1,
+                   for (unsigned i = 0; i < objective.size(); ++i) {
+                       verbose_stream() << objective[i].second 
+                                        << " * v" << objective[i].first << " ";
+                   }
+                   verbose_stream() << " + " << m_objective_consts[v] << "\n";);
+
+        unsigned num_nodes = get_num_vars();
+        unsigned num_edges = m_edges.size();
+        S.ensure_var(num_nodes + num_edges + m_objectives.size());
+        for (unsigned i = 0; i < num_nodes; ++i) {
+            numeral const& a = m_assignment[i];
+            rational fin = a.get_rational().to_rational();
+            rational inf = a.get_infinitesimal().to_rational();
+            mpq_inf q(fin.to_mpq(), inf.to_mpq());
+            S.set_value(i, q);
+        }
+        for (unsigned i = 0; i < num_nodes; ++i) {
+            enode * n = get_enode(i);
+            if (m_autil.is_zero(n->get_owner())) {
+                S.set_lower(v, mpq_inf(mpq(0), mpq(0)));
+                S.set_upper(v, mpq_inf(mpq(0), mpq(0)));
+                break;
+            }
+        }
+        svector<unsigned> vars;
+        unsynch_mpq_manager mgr;
+        scoped_mpq_vector coeffs(mgr);
+        coeffs.push_back(mpq(1));
+        coeffs.push_back(mpq(-1));
+        coeffs.push_back(mpq(-1));
+        vars.resize(3);
+        for (unsigned i = 0; i < num_edges; ++i) {
+            edge const& e = m_edges[i];
+            if (e.m_source == null_theory_var || e.m_target == null_theory_var) {
+                continue;
+            }
+            unsigned base_var = num_nodes + i;
+            vars[0] = e.m_target;
+            vars[1] = e.m_source;
+            vars[2] = base_var;
+            S.add_row(base_var, 3, vars.c_ptr(), coeffs.c_ptr());
+            // t - s <= w
+            // t - s - b = 0, b >= w
+            numeral const& w = e.m_offset;
+            rational fin = w.get_rational().to_rational();
+            rational inf = w.get_infinitesimal().to_rational();
+            mpq_inf q(fin.to_mpq(),inf.to_mpq());
+            S.set_upper(base_var, q);            
+        }
+        unsigned w = num_nodes + num_edges + v;
+
+        // add objective function as row.
+        coeffs.reset();
+        vars.reset();
+        for (unsigned i = 0; i < objective.size(); ++i) {
+            coeffs.push_back(objective[i].second.to_mpq());
+            vars.push_back(objective[i].first);
+        }
+        coeffs.push_back(mpq(1));
+        vars.push_back(w);
+        Simplex::row row = S.add_row(w, vars.size(), vars.c_ptr(), coeffs.c_ptr());
+        
+        TRACE("opt", S.display(tout); display(tout););
+        
+        // optimize    
+        lbool is_sat = S.make_feasible();
+        if (is_sat == l_undef) {
+            blocker = m.mk_false();
+            return inf_eps::infinity();        
+        }
+        TRACE("opt", S.display(tout); );    
+        SASSERT(is_sat != l_false);
+        lbool is_fin = S.minimize(w);
+        
+        switch (is_fin) {
+        case l_true: {
+            simplex::mpq_ext::eps_numeral const& val = S.get_value(w);
+            inf_rational r(-rational(val.first), -rational(val.second));
+            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);
+            core.reset();
+            for (; it != end; ++it) {
+                unsigned v = it->m_var;
+                if (num_nodes <= v && v < num_nodes + num_edges) {
+                    unsigned edge_id = v - num_nodes;
+                    literal lit = m_edges[edge_id].m_justification;
+                    get_context().literal2expr(lit, tmp);
+                    core.push_back(tmp);
+                }
+            }
+            for (unsigned i = 0; i < num_nodes; ++i) {
+                mpq_inf const& val = S.get_value(i);
+                rational q(val.first), eps(val.second);
+                numeral  a(q);
+                m_assignment[i] = a;
+                // TBD: if epsilon is != 0, then adjust a by some small fraction.
+            }
+            blocker = mk_gt(v, r);
+            IF_VERBOSE(10, verbose_stream() << blocker << "\n";);
+            return inf_eps(rational(0), r);
+        }
+        default:
+            TRACE("opt", tout << "unbounded\n"; );        
+            blocker = m.mk_false();
+            return inf_eps::infinity();        
+        }
+    }
+
+    template<typename Ext>
+    theory_var theory_dense_diff_logic<Ext>::add_objective(app* term) {
+        objective_term objective;
+        theory_var result = m_objectives.size();
+        rational q(1), r(0);
+        expr_ref_vector vr(get_manager());
+        if (!is_linear(get_manager(), term)) {
+            result = null_theory_var;
+        }
+        else if (internalize_objective(term, q, r, objective)) {
+            m_objectives.push_back(objective);
+            m_objective_consts.push_back(r);
+            m_objective_assignments.push_back(vr);
+        }
+        else {
+            result = null_theory_var;
+        }
+        return result;
+    }
+
+    template<typename Ext>
+    expr_ref theory_dense_diff_logic<Ext>::mk_gt(theory_var v, inf_rational const& val) {
+        return mk_ineq(v, val, true);
+    }
+
+    template<typename Ext>
+    expr_ref theory_dense_diff_logic<Ext>::mk_ge(
+        filter_model_converter& fm, theory_var v, inf_rational const& val) {
+        return mk_ineq(v, val, false);
+    }
+
+    template<typename Ext>
+    expr_ref theory_dense_diff_logic<Ext>::mk_ineq(theory_var v, inf_rational const& val, bool is_strict) {
+        ast_manager& m = get_manager();
+        objective_term const& t = m_objectives[v];
+        expr_ref e(m), f(m), f2(m);
+        if (t.size() == 1 && t[0].second.is_one()) {
+            f = get_enode(t[0].first)->get_owner();
+        }
+        else if (t.size() == 1 && t[0].second.is_minus_one()) {
+            f = m_autil.mk_uminus(get_enode(t[0].first)->get_owner());
+        }
+        else if (t.size() == 2 && t[0].second.is_one() && t[1].second.is_minus_one()) {
+            f = get_enode(t[0].first)->get_owner();
+            f2 = get_enode(t[1].first)->get_owner();
+            f = m_autil.mk_sub(f, f2); 
+        }
+        else if (t.size() == 2 && t[1].second.is_one() && t[0].second.is_minus_one()) {
+            f = get_enode(t[1].first)->get_owner();
+            f2 = get_enode(t[0].first)->get_owner();
+            f = m_autil.mk_sub(f, f2);
+        }
+        else {
+            // 
+            expr_ref_vector const& core = m_objective_assignments[v];
+            f = m.mk_and(core.size(), core.c_ptr());
+            if (is_strict) {
+                f = m.mk_not(f);
+            }
+            TRACE("arith", tout << "block: " << f << "\n";);
+            return f;
+        }
+        
+        e = m_autil.mk_numeral(val.get_rational(), m.get_sort(f));
+        
+        if (val.get_infinitesimal().is_neg()) {
+            if (is_strict) {
+                f = m_autil.mk_ge(f, e);
+            }
+            else {
+                expr_ref_vector const& core = m_objective_assignments[v];
+                f = m.mk_and(core.size(), core.c_ptr());
+            }
+        }
+        else {
+            if (is_strict) {
+                f = m_autil.mk_gt(f, e);
+            }
+            else {
+                f = m_autil.mk_ge(f, e);
+            }
+        }
+        return f;
+    }
+
 };
 
 #endif /* _THEORY_DENSE_DIFF_LOGIC_DEF_H_ */

src/smt/theory_diff_logic.cpp

@@ -21,6 +21,7 @@ Revision History:
 
 #include"rational.h"
 #include"theory_diff_logic_def.h"
+#include"sparse_matrix_def.h"
 
 namespace smt {
 
@@ -29,4 +30,10 @@ template class theory_diff_logic<sidl_ext>;
 template class theory_diff_logic<rdl_ext>;
 template class theory_diff_logic<srdl_ext>;
 
+
+};
+
+namespace simplex {
+template class simplex<mpq_ext>;
+template class sparse_matrix<mpq_ext>;
 };

src/smt/theory_diff_logic.h

@@ -37,7 +37,9 @@ Revision History:
 #include"smt_model_generator.h"
 #include"numeral_factory.h"
 #include"smt_clause.h"
-
+#include"theory_opt.h"
+#include"simplex.h"
+#include"simplex_def.h"
 
 // The DL theory can represent term such as n + k, where n is an enode and k is a numeral.
 namespace smt {
@@ -59,9 +61,11 @@ namespace smt {
     };
 
     template<typename Ext>
-    class theory_diff_logic : public theory, private Ext {
+    class theory_diff_logic : public theory, public theory_opt, private Ext {
         
         typedef typename Ext::numeral numeral;
+        typedef simplex::simplex<simplex::mpq_ext> Simplex;
+        typedef inf_eps_rational<inf_rational> inf_eps;
 
         class atom {
             bool_var  m_bvar;
@@ -158,18 +162,21 @@ namespace smt {
             unsigned      m_asserted_atoms_lim;
             unsigned      m_asserted_qhead_old;
         };
+        typedef dl_graph<GExt> Graph;
 
         smt_params &                   m_params;
         arith_util                     m_util;
         arith_eq_adapter               m_arith_eq_adapter;
         theory_diff_logic_statistics   m_stats;
-        dl_graph<GExt>                 m_graph;
-        theory_var                     m_zero_int; // cache the variable representing the zero variable.
-        theory_var                     m_zero_real; // cache the variable representing the zero variable.
+        Graph                          m_graph;
+        theory_var                     m_zero; // cache the variable representing the zero variable.
         int_vector                     m_scc_id;                  // Cheap equality propagation
         eq_prop_info_set               m_eq_prop_info_set;        // set of existing equality prop infos
         ptr_vector<eq_prop_info>       m_eq_prop_infos;
 
+        app_ref_vector                 m_terms;
+        svector<bool>                  m_signs;
+
         ptr_vector<atom>               m_atoms;
         ptr_vector<atom>               m_asserted_atoms;   // set of asserted atoms
         unsigned                       m_asserted_qhead;   
@@ -184,7 +191,16 @@ namespace smt {
 
         arith_factory *                m_factory;
         rational                       m_delta;
-        nc_functor                     m_nc_functor;        
+        nc_functor                     m_nc_functor;   
+
+        // For optimization purpose
+        typedef vector <std::pair<theory_var, rational> > objective_term;
+        vector<objective_term>         m_objectives;
+        vector<rational>               m_objective_consts;
+        vector<expr_ref_vector>        m_objective_assignments;
+        vector<Simplex::row>           m_objective_rows;
+        Simplex                        m_S;
+        unsigned                       m_num_simplex_edges;
 
         // Set a conflict due to a negative cycle.
         void set_neg_cycle_conflict();
@@ -210,8 +226,8 @@ namespace smt {
             m_params(params),
             m_util(m),
             m_arith_eq_adapter(*this, params, m_util),
-            m_zero_int(null_theory_var),
-            m_zero_real(null_theory_var),
+            m_zero(null_theory_var),
+            m_terms(m),
             m_asserted_qhead(0),
             m_num_core_conflicts(0),
             m_num_propagation_calls(0),
@@ -219,7 +235,8 @@ namespace smt {
             m_is_lia(true),
             m_non_diff_logic_exprs(false),
             m_factory(0),
-            m_nc_functor(*this) {
+            m_nc_functor(*this),
+            m_num_simplex_edges(0) {
         }            
 
         virtual ~theory_diff_logic() {
@@ -296,12 +313,33 @@ namespace smt {
         
         virtual void collect_statistics(::statistics & st) const;
 
-    private:        
+        
+        // -----------------------------------
+        //
+        // Optimization
+        //
+        // -----------------------------------
+
+        virtual inf_eps maximize(theory_var v, expr_ref& blocker, bool& has_shared);
+        virtual inf_eps value(theory_var v);
+        virtual theory_var add_objective(app* term);
+        virtual expr_ref mk_gt(theory_var v, inf_rational const& val);
+        virtual expr_ref mk_ge(filter_model_converter& fm, theory_var v, inf_rational const& val);
+
+        bool internalize_objective(expr * n, rational const& m, rational& r, objective_term & objective);
+
+    private:     
+
+        expr_ref mk_ineq(theory_var v, inf_rational const& val, bool is_strict);
 
         virtual void new_eq_eh(theory_var v1, theory_var v2, justification& j);
 
         virtual void new_diseq_eh(theory_var v1, theory_var v2, justification& j);
 
+        bool decompose_linear(app_ref_vector& args, svector<bool>& signs);
+
+        bool is_sign(expr* n, bool& sign);
+
         bool is_negative(app* n, app*& m);
 
         void del_atoms(unsigned old_size);
@@ -336,18 +374,26 @@ namespace smt {
 
         void get_implied_bound_antecedents(edge_id bridge_edge, edge_id subsumed_edge, conflict_resolution & cr);
 
-        theory_var get_zero(sort* s) const { return m_util.is_int(s)?m_zero_int:m_zero_real; }
-
-        theory_var get_zero(expr* e) const { return get_zero(get_manager().get_sort(e)); }
+        theory_var get_zero() const { return m_zero; }
 
         void inc_conflicts();
 
+        // Optimization:
+        // convert variables, edges and objectives to simplex.
+        unsigned node2simplex(unsigned v);
+        unsigned edge2simplex(unsigned e);
+        unsigned obj2simplex(unsigned v);
+        unsigned num_simplex_vars();
+        bool is_simplex_edge(unsigned e);
+        unsigned simplex2edge(unsigned e);
+        void update_simplex(Simplex& S);
     };
 
     struct idl_ext {
         // TODO: It doesn't need to be a rational, but a bignum integer.
         static const bool m_int_theory = true;
         typedef rational numeral;
+        typedef rational fin_numeral;
         numeral     m_epsilon;
         idl_ext() : m_epsilon(1) {}
     };
@@ -356,6 +402,7 @@ namespace smt {
         // TODO: It doesn't need to be a rational, but a bignum integer.
         static const bool m_int_theory = true;
         typedef s_integer numeral;
+        typedef s_integer fin_numeral;
         numeral m_epsilon;
         sidl_ext() : m_epsilon(1) {}
     };
@@ -363,18 +410,19 @@ namespace smt {
     struct rdl_ext {
         static const bool m_int_theory = false;
         typedef inf_int_rational numeral;
-        numeral      m_epsilon;
+        typedef rational fin_numeral;
+        numeral      m_epsilon;        
         rdl_ext() : m_epsilon(rational(), true) {}
     };
 
     struct srdl_ext {
         static const bool m_int_theory = false;
         typedef inf_s_integer numeral;
+        typedef s_integer fin_numeral;
         numeral m_epsilon;
         srdl_ext() : m_epsilon(s_integer(0),true) {}
     };
 
-
     typedef theory_diff_logic<idl_ext>  theory_idl;
     typedef theory_diff_logic<sidl_ext> theory_fidl;
     typedef theory_diff_logic<rdl_ext> theory_rdl;

src/smt/theory_diff_logic_def.h

@@ -28,6 +28,8 @@ Revision History:
 #include"ast_pp.h"
 #include"warning.h"
 #include"smt_model_generator.h"
+#include"model_implicant.h"
+
 
 using namespace smt;
 
@@ -69,12 +71,7 @@ void theory_diff_logic<Ext>::init(context * ctx) {
     zero = m_util.mk_numeral(rational(0), true);
     e = ctx->mk_enode(zero, false, false, true);
     SASSERT(!is_attached_to_var(e));
-    m_zero_int = mk_var(e);
-
-    zero = m_util.mk_numeral(rational(0), false);
-    e = ctx->mk_enode(zero, false, false, true);
-    SASSERT(!is_attached_to_var(e));
-    m_zero_real = mk_var(e);
+    m_zero = mk_var(e);
 }
 
 
@@ -82,8 +79,10 @@ template<typename Ext>
 bool theory_diff_logic<Ext>::internalize_term(app * term) {
     bool result = null_theory_var != mk_term(term);
     CTRACE("arith", !result, tout << "Did not internalize " << mk_pp(term, get_manager()) << "\n";);
-    TRACE("non_diff_logic", tout << "Terms may not be internalized\n";);
-    found_non_diff_logic_expr(term);
+    if (!result) {
+        TRACE("non_diff_logic", tout << "Terms may not be internalized\n";);
+        found_non_diff_logic_expr(term);
+    }
     return result;
 }
 
@@ -159,6 +158,7 @@ template<typename Ext>
 void theory_diff_logic<Ext>::found_non_diff_logic_expr(expr * n) {
     if (!m_non_diff_logic_exprs) {
         TRACE("non_diff_logic", tout << "found non diff logic expression:\n" << mk_pp(n, get_manager()) << "\n";);
+        IF_VERBOSE(0, verbose_stream() << "(smt.diff_logic: non-diff logic expression " << mk_pp(n, get_manager()) << ")\n";); 
         get_context().push_trail(value_trail<context, bool>(m_non_diff_logic_exprs));
         m_non_diff_logic_exprs = true;
     }
@@ -177,7 +177,6 @@ bool theory_diff_logic<Ext>::internalize_atom(app * n, bool gate_ctx) {
     bool is_ge = m_util.is_ge(n);
     bool_var bv;
     rational kr;
-    app * x, *y, *z;
     theory_var source, target; // target - source <= k
     app * lhs = to_app(n->get_arg(0));
     app * rhs = to_app(n->get_arg(1));
@@ -191,25 +190,26 @@ bool theory_diff_logic<Ext>::internalize_atom(app * n, bool gate_ctx) {
     }
     numeral k(kr);
 
-    bool is_add = m_util.is_add(lhs) && lhs->get_num_args() == 2;
-
-    if (is_add) {
-        x = to_app(lhs->get_arg(0));
-        y = to_app(lhs->get_arg(1));
+    m_terms.reset();
+    m_signs.reset();
+    m_terms.push_back(lhs);
+    m_signs.push_back(true);
+    if (!decompose_linear(m_terms, m_signs)) {
+        found_non_diff_logic_expr(n);        
+        return false;
     }
-
-    if (is_add && is_negative(x, z)) {  
-        target = mk_var(y);
-        source = mk_var(z);
-    }
-    else if (is_add && is_negative(y, z)) {
-        target = mk_var(x);
-        source = mk_var(z);
+    if (m_terms.size() == 2 && m_signs[0] != m_signs[1]) {
+        target = mk_var(m_terms[0].get());
+        source = mk_var(m_terms[1].get());
+        if (!m_signs[0]) {
+            std::swap(target, source);
+        }
     }
     else {
         target = mk_var(lhs);
-        source = get_zero(lhs);
+        source = get_zero();
     }
+
     if (is_ge) {
         std::swap(target, source);
         k.neg();
@@ -273,6 +273,8 @@ bool theory_diff_logic<Ext>::internalize_atom(app * n, bool gate_ctx) {
 template<typename Ext>
 void theory_diff_logic<Ext>::internalize_eq_eh(app * atom, bool_var v) {
     context & ctx  = get_context();
+    ast_manager& m = get_manager();
+    TRACE("arith", tout << mk_pp(atom, m) << "\n";);
     app * lhs      = to_app(atom->get_arg(0));
     app * rhs      = to_app(atom->get_arg(1));
     app * s;
@@ -339,7 +341,13 @@ void theory_diff_logic<Ext>::pop_scope_eh(unsigned num_scopes) {
     m_asserted_atoms.shrink(s.m_asserted_atoms_lim);
     m_asserted_qhead = s.m_asserted_qhead_old;
     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) {
+        m_S.reset();
+        m_num_simplex_edges = 0;
+        m_objective_rows.reset();
+    }
     theory::pop_scope_eh(num_scopes);
 }
 
@@ -353,7 +361,7 @@ final_check_status theory_diff_logic<Ext>::final_check_eh() {
 
     TRACE("arith_final", display(tout); );
     // either will already be zero (as we don't do mixed constraints).
-    m_graph.set_to_zero(m_zero_int, m_zero_real);
+    m_graph.set_to_zero(m_zero);
     SASSERT(is_consistent());
     if (m_non_diff_logic_exprs) {
         return FC_GIVEUP; 
@@ -378,6 +386,70 @@ void theory_diff_logic<Ext>::del_atoms(unsigned old_size) {
 }
 
 
+template<typename Ext>
+bool theory_diff_logic<Ext>::decompose_linear(app_ref_vector& terms, svector<bool>& signs) {
+    for (unsigned i = 0; i < terms.size(); ++i) {
+        app* n = terms[i].get();
+        if (m_util.is_add(n)) {
+            expr* arg = n->get_arg(0);
+            if (!is_app(arg)) return false;
+            terms[i] = to_app(arg);
+            for (unsigned j = 1; j < n->get_num_args(); ++j) {
+                arg = n->get_arg(j);
+                if (!is_app(arg)) return false;
+                terms.push_back(to_app(arg));
+                signs.push_back(signs[i]);
+            }
+            --i;
+            continue;
+        }
+        expr* x, *y;
+        bool sign;
+        if (m_util.is_mul(n, x, y)) {
+            if (is_sign(x, sign) && is_app(y)) {
+                terms[i] = to_app(y);
+                signs[i] = (signs[i] == sign);
+                --i;
+            }
+            else if (is_sign(y, sign) && is_app(x)) {
+                terms[i] = to_app(x);
+                signs[i] = (signs[i] == sign);
+                --i;
+            }
+            continue;
+        }
+        if (m_util.is_uminus(n, x) && is_app(x)) {
+            terms[i] = to_app(x);
+            signs[i] = !signs[i];
+            --i;
+            continue;
+        }
+    }
+    return true;
+}
+
+template<typename Ext>
+bool theory_diff_logic<Ext>::is_sign(expr* n, bool& sign) {
+    rational r;
+    expr* x;
+    if (m_util.is_numeral(n, r)) {
+        if (r.is_one()) {
+            sign = true;
+            return true;
+        }
+        if (r.is_minus_one()) {
+            sign = false;
+            return true;
+        }
+    }
+    else if (m_util.is_uminus(n, x)) {
+        if (is_sign(x, sign)) {
+            sign = !sign;
+            return true;
+        }
+    }
+    return false;
+}
 
 template<typename Ext>
 bool theory_diff_logic<Ext>::is_negative(app* n, app*& m) { 
@@ -474,7 +546,7 @@ void theory_diff_logic<Ext>::propagate_core() {
 template<typename Ext>
 bool theory_diff_logic<Ext>::propagate_atom(atom* a) {
     context& ctx = get_context();
-    TRACE("arith", a->display(*this, tout); );
+    TRACE("arith", a->display(*this, tout); tout << "\n";);
     if (ctx.inconsistent()) {
         return false;
     }
@@ -650,6 +722,7 @@ theory_var theory_diff_logic<Ext>::mk_term(app* n) {
     app* a, *offset;
     theory_var source, target;
     enode* e;
+    context& ctx = get_context();
 
     TRACE("arith", tout << mk_pp(n, get_manager()) << "\n";);
 
@@ -660,6 +733,13 @@ theory_var theory_diff_logic<Ext>::mk_term(app* n) {
     else if (is_offset(n, a, offset, r)) {
         // n = a + k
         source = mk_var(a);
+        for (unsigned i = 0; i < n->get_num_args(); ++i) {
+            expr* arg = n->get_arg(i);
+            std::cout << "internalize: " << mk_pp(arg, get_manager()) << " " << ctx.e_internalized(arg) << "\n";
+            if (!ctx.e_internalized(arg)) {
+                ctx.internalize(arg, false);
+            }
+        }
         e = get_context().mk_enode(n, false, false, true);
         target = mk_var(e);
         numeral k(r);
@@ -698,7 +778,7 @@ theory_var theory_diff_logic<Ext>::mk_num(app* n, rational const& r) {
     enode* e = 0;
     context& ctx = get_context();
     if (r.is_zero()) {
-        v = get_zero(n);
+        v = get_zero();
     }
     else if (ctx.e_internalized(n)) {
         e = ctx.get_enode(n);
@@ -706,7 +786,8 @@ theory_var theory_diff_logic<Ext>::mk_num(app* n, rational const& r) {
         SASSERT(v != null_theory_var);
     }
     else {
-        theory_var zero = get_zero(n);
+        theory_var zero = get_zero();
+        SASSERT(n->get_num_args() == 0);
         e = ctx.mk_enode(n, false, false, true);
         v = mk_var(e);
         // internalizer is marking enodes as interpreted whenever the associated ast is a value and a constant.
@@ -763,8 +844,7 @@ void theory_diff_logic<Ext>::reset_eh() {
         dealloc(m_atoms[i]);
     }
     m_graph            .reset();
-    m_zero_int          = null_theory_var;
-    m_zero_real         = null_theory_var;
+    m_zero              = null_theory_var;
     m_atoms            .reset();
     m_asserted_atoms   .reset();
     m_stats            .reset();
@@ -775,6 +855,9 @@ void theory_diff_logic<Ext>::reset_eh() {
     m_agility               = 0.5;
     m_is_lia                = true;
     m_non_diff_logic_exprs  = false;
+    m_objectives      .reset();
+    m_objective_consts.reset();
+    m_objective_assignments.reset();
     theory::reset_eh();
 }
 
@@ -993,5 +1076,313 @@ void theory_diff_logic<Ext>::get_implied_bound_antecedents(edge_id bridge_edge,
     m_graph.explain_subsumed_lazy(bridge_edge, subsumed_edge, f);
 }
 
+template<typename Ext>
+unsigned theory_diff_logic<Ext>::node2simplex(unsigned v) {
+    return m_objectives.size() + 2*v + 1;
+}
+template<typename Ext>
+unsigned theory_diff_logic<Ext>::edge2simplex(unsigned e) {
+    return m_objectives.size() + 2*e;
+}
+template<typename Ext>
+unsigned theory_diff_logic<Ext>::obj2simplex(unsigned e) {
+    return e;
+}
+
+template<typename Ext>
+unsigned theory_diff_logic<Ext>::num_simplex_vars() {
+    return m_objectives.size() + std::max(2*m_graph.get_num_edges(),2*m_graph.get_num_nodes()+1);
+}
+
+template<typename Ext>
+bool theory_diff_logic<Ext>::is_simplex_edge(unsigned e) {
+    if (e < m_objectives.size()) return false;
+    e -= m_objectives.size();
+    return (0 == (e & 0x1));
+}
+
+template<typename Ext> 
+unsigned theory_diff_logic<Ext>::simplex2edge(unsigned e) {
+    SASSERT(is_simplex_edge(e));
+    return (e - m_objectives.size())/2;
+}
+
+template<typename Ext> 
+void theory_diff_logic<Ext>::update_simplex(Simplex& S) {
+    unsigned num_nodes = m_graph.get_num_nodes();
+    vector<dl_edge<GExt> > const& es = m_graph.get_all_edges();
+    S.ensure_var(num_simplex_vars());
+    for (unsigned i = 0; i < num_nodes; ++i) {
+        numeral const& a = m_graph.get_assignment(i);
+        rational fin = a.get_rational().to_rational();
+        rational inf = a.get_infinitesimal().to_rational();
+        mpq_inf q(fin.to_mpq(), inf.to_mpq());
+        S.set_value(node2simplex(i), q);
+    }
+    S.set_lower(node2simplex(get_zero()), mpq_inf(mpq(0), mpq(0)));
+    S.set_upper(node2simplex(get_zero()), mpq_inf(mpq(0), mpq(0)));
+    svector<unsigned> vars;
+    unsynch_mpq_manager mgr;
+    scoped_mpq_vector coeffs(mgr);
+    coeffs.push_back(mpq(1));
+    coeffs.push_back(mpq(-1));
+    coeffs.push_back(mpq(-1));
+    vars.resize(3);
+    for (unsigned i = m_num_simplex_edges; i < es.size(); ++i) {
+        //    t - s <= w 
+        // =>
+        //    t - s - b = 0, b >= w
+        dl_edge<GExt> const& e = es[i];
+        unsigned base_var = edge2simplex(i);
+        vars[0] = node2simplex(e.get_target());
+        vars[1] = node2simplex(e.get_source());
+        vars[2] = base_var;
+        S.add_row(base_var, 3, vars.c_ptr(), coeffs.c_ptr());        
+    }
+    m_num_simplex_edges = es.size();
+    for (unsigned i = 0; i < es.size(); ++i) {
+        dl_edge<GExt> const& e = es[i];
+        unsigned base_var = edge2simplex(i);
+        if (e.is_enabled()) {
+            numeral const& w = e.get_weight();
+            rational fin = w.get_rational().to_rational();
+            rational inf = w.get_infinitesimal().to_rational();
+            mpq_inf q(fin.to_mpq(),inf.to_mpq());
+            S.set_upper(base_var, q);
+        }
+        else {
+            S.unset_upper(base_var);
+        }
+    }
+    for (unsigned v = m_objective_rows.size(); v < m_objectives.size(); ++v) {
+        unsigned w = obj2simplex(v);
+        objective_term const& objective = m_objectives[v];
+
+        // add objective function as row.
+        coeffs.reset();
+        vars.reset();
+        for (unsigned i = 0; i < objective.size(); ++i) {
+            coeffs.push_back(objective[i].second.to_mpq());
+            vars.push_back(node2simplex(objective[i].first));
+        }
+        coeffs.push_back(mpq(1));
+        vars.push_back(w);
+        Simplex::row row = S.add_row(w, vars.size(), vars.c_ptr(), coeffs.c_ptr());
+        m_objective_rows.push_back(row);
+    }
+}
+
+template<typename Ext>
+typename theory_diff_logic<Ext>::inf_eps theory_diff_logic<Ext>::value(theory_var v) {
+     objective_term const& objective = m_objectives[v];   
+     inf_eps r = inf_eps(m_objective_consts[v]);
+     for (unsigned i = 0; i < objective.size(); ++i) {
+         numeral n = m_graph.get_assignment(v);
+         rational r1 = n.get_rational().to_rational();
+         rational r2 = n.get_infinitesimal().to_rational();
+         r += objective[i].second * inf_eps(rational(0), inf_rational(r1, r2));
+     }
+     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) {
+    
+    has_shared = false;
+    Simplex& S = m_S;
+    ast_manager& m = get_manager();
+
+    update_simplex(S);
+    objective_term const& objective = m_objectives[v];
+
+    TRACE("arith",
+          for (unsigned i = 0; i < objective.size(); ++i) {
+              tout << "Coefficient " << objective[i].second 
+                   << " of theory_var " << objective[i].first << "\n";
+          }
+          tout << "Free coefficient " << m_objective_consts[v] << "\n";
+          );
+
+    TRACE("opt", S.display(tout); display(tout););
+    
+    // optimize    
+    lbool is_sat = S.make_feasible();
+    if (is_sat == l_undef) {
+        blocker = m.mk_false();
+        return inf_eps::infinity();        
+    }
+    TRACE("opt", S.display(tout); );    
+    SASSERT(is_sat != l_false);
+    unsigned w = obj2simplex(v);
+    lbool is_fin = S.minimize(w);
+    switch (is_fin) {
+    case l_true: {
+        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);
+        core.reset();
+        for (; it != end; ++it) {
+            unsigned v = it->m_var;
+            if (is_simplex_edge(v)) {
+                unsigned edge_id = simplex2edge(v);
+                literal lit = m_graph.get_explanation(edge_id);
+                get_context().literal2expr(lit, tmp);
+                core.push_back(tmp);
+            }
+        }
+        compute_delta();
+        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);
+            m_graph.set_assignment(i, numeral(r));
+        }
+        blocker = mk_gt(v, r);
+        return inf_eps(rational(0), r + m_objective_consts[v]);
+    }
+    default:
+        TRACE("opt", tout << "unbounded\n"; );        
+        blocker = m.mk_false();
+        return inf_eps::infinity();        
+    }
+}
+
+template<typename Ext>
+theory_var theory_diff_logic<Ext>::add_objective(app* term) {
+    objective_term objective;
+    theory_var result = m_objectives.size();
+    rational q(1), r(0);
+    expr_ref_vector vr(get_manager());
+    if (!is_linear(get_manager(), term)) {
+        result = null_theory_var;
+    }
+    else if (internalize_objective(term, q, r, objective)) {
+        m_objectives.push_back(objective);
+        m_objective_consts.push_back(r);
+        m_objective_assignments.push_back(vr);
+    }
+    else {
+        result = null_theory_var;
+    }
+    return result; 
+}
+
+template<typename Ext>
+expr_ref theory_diff_logic<Ext>::mk_ineq(theory_var v, inf_rational const& val, bool is_strict) {
+    ast_manager& m = get_manager();
+    objective_term const& t = m_objectives[v];
+    expr_ref e(m), f(m), f2(m);
+    if (t.size() == 1 && t[0].second.is_one()) {
+        f = get_enode(t[0].first)->get_owner();
+    }
+    else if (t.size() == 1 && t[0].second.is_minus_one()) {
+        f = m_util.mk_uminus(get_enode(t[0].first)->get_owner());
+    }
+    else if (t.size() == 2 && t[0].second.is_one() && t[1].second.is_minus_one()) {
+        f = get_enode(t[0].first)->get_owner();
+        f2 = get_enode(t[1].first)->get_owner();
+        f = m_util.mk_sub(f, f2); 
+    }
+    else if (t.size() == 2 && t[1].second.is_one() && t[0].second.is_minus_one()) {
+        f = get_enode(t[1].first)->get_owner();
+        f2 = get_enode(t[0].first)->get_owner();
+        f = m_util.mk_sub(f, f2);
+    }
+    else {
+        // 
+        expr_ref_vector const& core = m_objective_assignments[v];
+        f = m.mk_and(core.size(), core.c_ptr());
+        if (is_strict) {
+            f = m.mk_not(f);
+        }
+        TRACE("arith", tout << "block: " << f << "\n";);
+        return f;
+    }
+
+    inf_rational new_val = val; // - inf_rational(m_objective_consts[v]);
+    e = m_util.mk_numeral(new_val.get_rational(), m.get_sort(f));
+    
+    if (new_val.get_infinitesimal().is_neg()) {
+        if (is_strict) {
+            f = m_util.mk_ge(f, e);
+        }
+        else {
+            expr_ref_vector const& core = m_objective_assignments[v];
+            f = m.mk_and(core.size(), core.c_ptr());            
+        }
+    }
+    else {
+        if (is_strict) {
+            f = m_util.mk_gt(f, e);
+        }
+        else {
+            f = m_util.mk_ge(f, e);
+        }
+    }
+    return f;
+}
+
+template<typename Ext>
+expr_ref theory_diff_logic<Ext>::mk_gt(theory_var v, inf_rational const& val) {
+    return mk_ineq(v, val, true);
+}
+
+template<typename Ext>
+expr_ref theory_diff_logic<Ext>::mk_ge(filter_model_converter& fm, theory_var v, inf_rational const& val) {
+    return mk_ineq(v, val, false);
+}
+
+#if 0
+    context & ctx = get_context();
+    model_ref mdl;
+    ctx.get_model(mdl);
+    ptr_vector<expr> formulas(ctx.get_num_asserted_formulas(), ctx.get_asserted_formulas());
+    ast_manager& m = get_manager();
+    model_implicant impl_extractor(m);
+    expr_ref_vector implicants = impl_extractor.minimize_literals(formulas, mdl);
+    return m.mk_and(o, m.mk_not(m.mk_and(implicants.size(), implicants.c_ptr())));
+#endif
+
+template<typename Ext>
+bool theory_diff_logic<Ext>::internalize_objective(expr * n, rational const& m, rational& q, objective_term & objective) {
+
+    // Compile term into objective_term format
+    rational r;
+    expr* x, *y;
+    if (m_util.is_numeral(n, r)) {
+        q += r;
+    }
+    else if (m_util.is_add(n)) {
+        for (unsigned i = 0; i < to_app(n)->get_num_args(); ++i) {
+            if (!internalize_objective(to_app(n)->get_arg(i), m, q, objective)) {
+                return false;
+            }
+        }
+    }
+    else if (m_util.is_mul(n, x, y) && m_util.is_numeral(x, r)) {
+        return internalize_objective(y, m*r, q, objective);
+    }
+    else if (m_util.is_mul(n, y, x) && m_util.is_numeral(x, r)) {
+        return internalize_objective(y, m*r, q, objective);
+    }
+    else if (!is_app(n)) {
+        return false;
+    }
+    else if (to_app(n)->get_family_id() == m_util.get_family_id()) {
+        return false;
+    }
+    else {
+        theory_var v = mk_var(to_app(n));
+        objective.push_back(std::make_pair(v, m));
+    }
+    return true;
+}
+
 #endif /* _THEORY_DIFF_LOGIC_DEF_H_ */
 

src/smt/theory_opt.cpp

@@ -0,0 +1,84 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    theory_opt.cpp
+
+Abstract:
+
+    Interface utilities used by optimization providing 
+    theory solvers.
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2013-10-18
+
+Notes:
+
+--*/
+
+#include "arith_decl_plugin.h"
+#include "smt_types.h"
+#include "smt_theory.h"
+#include "theory_opt.h"
+
+namespace smt {
+    bool theory_opt::is_linear(ast_manager& m, expr* term) {
+        arith_util a(m);
+        ptr_vector<expr> todo;
+        ast_mark mark;
+        todo.push_back(term);
+        expr* t1, *t2;
+        while (!todo.empty()) {
+            term = todo.back();
+            todo.pop_back();
+            if (mark.is_marked(term)) {
+                continue;
+            }
+            mark.mark(term, true);
+            if (!is_app(term)) {
+                return false;
+            }
+            app* t = to_app(term);
+            if (t->get_family_id() != a.get_family_id()) {
+                // done
+            }
+            else if (a.is_add(t) || a.is_to_real(t) || a.is_to_int(t) ||
+                a.is_uminus(t) || a.is_numeral(t) || a.is_sub(t)) {
+                todo.append(t->get_num_args(), t->get_args());
+            }
+            else if (a.is_mul(t, t1, t2)) {
+                if (is_numeral(a, t1)) {
+                    todo.push_back(t2);
+                }
+                else if (is_numeral(a, t2)) {
+                    todo.push_back(t1);
+                }
+                else {
+                    return false;
+                }
+            }
+            else {
+                return false;
+            } 
+        }
+        return true;
+    }
+
+    
+    bool theory_opt::is_numeral(arith_util& a, expr* term) {
+        while (true) {
+            if (a.is_uminus(term) || a.is_to_real(term) || a.is_to_int(term)) {
+                term = to_app(term)->get_arg(0);
+            }
+            else if (a.is_numeral(term)) {
+                return true;
+            }
+            else {
+                return false;
+            }
+        }
+    }
+
+};

src/smt/theory_opt.h

@@ -0,0 +1,42 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    theory_opt.h
+
+Abstract:
+
+    Interface utilities used by optimization providing 
+    theory solvers.
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2013-10-18
+
+Notes:
+
+--*/
+
+#include "inf_rational.h"
+#include "inf_eps_rational.h"
+#include "arith_decl_plugin.h"
+
+#ifndef _THEORY_OPT_H_
+#define _THEORY_OPT_H_
+
+class filter_model_converter;
+namespace smt {
+    class theory_opt {
+    public:
+        typedef inf_eps_rational<inf_rational> inf_eps;
+        virtual inf_eps value(theory_var) = 0;
+        virtual inf_eps maximize(theory_var v, expr_ref& blocker, bool& has_shared) = 0; 
+        virtual theory_var add_objective(app* term) = 0;
+        virtual expr_ref mk_ge(filter_model_converter& fm, theory_var v, inf_eps const& val) { UNREACHABLE(); return expr_ref(*((ast_manager*)0)); }
+        bool is_linear(ast_manager& m, expr* term);
+        bool is_numeral(arith_util& a, expr* term);
+    };
+}
+
+#endif

src/smt/theory_pb.h

@@ -0,0 +1,322 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    theory_pb.h
+
+Abstract:
+
+    Cardinality theory plugin.
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2013-11-05
+
+Notes:
+
+    This custom theory handles cardinality constraints
+    It performs unit propagation and switches to creating
+    sorting circuits if it keeps having to propagate (create new clauses).
+--*/
+
+#include "smt_theory.h"
+#include "pb_decl_plugin.h"
+#include "smt_clause.h"
+#include "theory_pb_params.h"
+#include "simplex.h"
+
+namespace smt {
+    class theory_pb : public theory {
+
+        struct psort_expr;
+        class  pb_justification;
+        class  pb_model_value_proc;
+        class  unwatch_ge;
+        class  rewatch_vars;
+        class  negate_ineq;
+        class  remove_var;
+        class  undo_bound;
+        typedef rational numeral;
+        typedef simplex::simplex<simplex::mpz_ext> simplex;
+        typedef simplex::row row;
+        typedef simplex::row_iterator row_iterator;
+        typedef unsynch_mpq_inf_manager eps_manager;
+        typedef _scoped_numeral<eps_manager> scoped_eps_numeral;
+
+        struct arg_t : public vector<std::pair<literal, numeral> > {
+            numeral         m_k;        // invariants: m_k > 0, coeffs[i] > 0
+
+            unsigned get_hash() const;
+            bool operator==(arg_t const& other) const;
+
+            numeral const& k() const { return m_k; }
+
+            struct hash {
+                unsigned operator()(arg_t const& i) const { return i.get_hash(); }
+            };
+            struct eq {
+                bool operator()(arg_t const& a, arg_t const& b) const { 
+                    return a == b;
+                }
+            };
+            struct child_hash {
+                unsigned operator()(arg_t const& args, unsigned idx) const {
+                    return args[idx].first.hash() ^ args[idx].second.hash();
+                }
+            };
+            struct kind_hash {
+                unsigned operator()(arg_t const& args) const {
+                    return args.size();
+                }
+            };   
+
+            void remove_negations();         
+
+            void negate();
+
+            lbool normalize(bool is_eq);
+
+            void  prune(bool is_eq);
+            
+            literal lit(unsigned i) const { 
+                return (*this)[i].first; 
+            }
+
+            numeral const & coeff(unsigned i) const { return (*this)[i].second; }
+
+            std::ostream& display(context& ctx, std::ostream& out, bool values = false) const;
+
+            app_ref to_expr(bool is_eq, context& ctx, ast_manager& m);
+
+            bool well_formed() const;
+        };
+
+        struct stats {
+            unsigned m_num_conflicts;
+            unsigned m_num_propagations;
+            unsigned m_num_predicates;
+            unsigned m_num_compiles;
+            unsigned m_num_compiled_vars;
+            unsigned m_num_compiled_clauses;
+            void reset() { memset(this, 0, sizeof(*this)); }
+            stats() { reset(); }
+        };
+
+
+        struct ineq {
+            unsynch_mpz_manager& m_mpz;    // mpz manager.
+            literal         m_lit;      // literal repesenting predicate
+            bool            m_is_eq;    // is this an = or >=.
+            arg_t           m_args[2];  // encode args[0]*coeffs[0]+...+args[n-1]*coeffs[n-1] >= k();
+            // Watch the first few positions until the sum satisfies:
+            // sum coeffs[i] >= m_lower + max_watch            
+            scoped_mpz      m_max_watch;    // maximal coefficient.
+            unsigned        m_watch_sz;     // number of literals being watched.
+            scoped_mpz      m_watch_sum;    // maximal sum of watch literals.
+            // Watch infrastructure for = and unassigned >=:
+            unsigned        m_nfixed;       // number of variables that are fixed.
+            scoped_mpz      m_max_sum;      // maximal possible sum.
+            scoped_mpz      m_min_sum;      // minimal possible sum.
+            unsigned        m_num_propagations;
+            unsigned        m_compilation_threshold;
+            lbool           m_compiled;
+            
+            ineq(unsynch_mpz_manager& m, literal l, bool is_eq) : 
+                m_mpz(m), m_lit(l), m_is_eq(is_eq), 
+                m_max_watch(m), m_watch_sum(m), 
+                m_max_sum(m), m_min_sum(m) {
+                reset();
+            }
+
+            arg_t const& args() const { return m_args[m_lit.sign()]; }
+            arg_t& args() { return m_args[m_lit.sign()]; }
+
+            literal lit() const { return m_lit; }
+            numeral const & k() const { return args().m_k; }
+            mpz const & mpz_k() const { return k().to_mpq().numerator(); }
+
+            literal lit(unsigned i) const { return args()[i].first; }
+            numeral const & coeff(unsigned i) const { return args()[i].second; }
+            class mpz const& ncoeff(unsigned i) const { return coeff(i).to_mpq().numerator(); }
+
+            unsigned size() const { return args().size(); }
+
+            scoped_mpz const& watch_sum() const { return m_watch_sum; }
+            scoped_mpz const& max_watch() const { return m_max_watch; }
+            void set_max_watch(mpz const& n) { m_max_watch = n; }
+            unsigned watch_size() const { return m_watch_sz; }
+
+            // variable watch infrastructure
+            scoped_mpz const& min_sum() const { return m_min_sum; }
+            scoped_mpz const& max_sum() const { return m_max_sum; }
+            unsigned nfixed() const { return m_nfixed; }
+            bool vwatch_initialized() const { return !m_mpz.is_zero(max_sum()); }
+            void vwatch_reset() { m_min_sum.reset(); m_max_sum.reset(); m_nfixed = 0; }
+
+            unsigned find_lit(bool_var v, unsigned begin, unsigned end) {
+                while (lit(begin).var() != v) {
+                    ++begin;
+                    SASSERT(begin < end);
+                }
+                return begin;
+            }
+
+            void reset();
+
+            void negate();
+
+            lbool normalize();
+
+            void unique();
+
+            void prune();
+
+            void post_prune();
+
+            app_ref to_expr(context& ctx, ast_manager& m);
+
+            bool is_eq() const { return m_is_eq; }
+            bool is_ge() const { return !m_is_eq; }
+
+        };
+
+        struct row_info {
+            unsigned     m_slack;   // slack variable in simplex tableau
+            numeral      m_bound;   // bound
+            arg_t        m_rep;     // representative
+            row_info(theory_var slack, numeral const& b, arg_t const& r):
+                m_slack(slack), m_bound(b), m_rep(r) {}
+            row_info(): m_slack(0) {}
+        };
+
+        typedef ptr_vector<ineq> watch_list;
+        typedef map<arg_t, bool_var, arg_t::hash, arg_t::eq> arg_map;
+
+        theory_pb_params         m_params;        
+        u_map<watch_list*>       m_lwatch;      // per literal.
+        u_map<watch_list*>       m_vwatch;      // per variable.
+        u_map<ineq*>             m_ineqs;       // per inequality.
+        arg_map                  m_ineq_rep;       // Simplex: representative inequality
+        u_map<row_info>          m_ineq_row_info;  // Simplex: row information per variable
+        uint_set                 m_vars;           // Simplex: 0-1 variables.
+        simplex                  m_simplex;        // Simplex: tableau
+        literal_vector           m_explain_lower;  // Simplex: explanations for lower bounds
+        literal_vector           m_explain_upper;  // Simplex: explanations for upper bounds
+        unsynch_mpq_inf_manager  m_mpq_inf_mgr;    // Simplex: manage inf_mpq numerals
+        mutable unsynch_mpz_manager      m_mpz_mgr;        // Simplex: manager mpz numerals
+        unsigned_vector          m_ineqs_trail;
+        unsigned_vector          m_ineqs_lim;
+        literal_vector           m_literals;    // temporary vector
+        pb_util                  m_util;
+        stats                    m_stats;
+        ptr_vector<ineq>         m_to_compile;  // inequalities to compile.
+        unsigned                 m_conflict_frequency;
+        bool                     m_learn_complements;
+        bool                     m_enable_compilation;
+        bool                     m_enable_simplex;
+        rational                 m_max_compiled_coeff;
+
+        // internalize_atom:
+        literal compile_arg(expr* arg);
+        void add_watch(ineq& c, unsigned index);
+        void del_watch(watch_list& watch, unsigned index, ineq& c, unsigned ineq_index);
+        void init_watch_literal(ineq& c);
+        void init_watch_var(ineq& c);
+        void clear_watch(ineq& c);
+        void watch_literal(literal lit, ineq* c);
+        void watch_var(bool_var v, ineq* c);
+        void unwatch_literal(literal w, ineq* c);
+        void unwatch_var(bool_var v, ineq* c);
+        void remove(ptr_vector<ineq>& ineqs, ineq* c);
+        bool assign_watch_ge(bool_var v, bool is_true, watch_list& watch, unsigned index);
+        void assign_watch(bool_var v, bool is_true, ineq& c);
+        void assign_ineq(ineq& c, bool is_true);
+        void assign_eq(ineq& c, bool is_true);
+
+        // simplex:
+        literal set_explain(literal_vector& explains, unsigned var, literal expl);
+        bool update_bound(bool_var v, literal explain, bool is_lower, mpq_inf const& bound);
+        bool check_feasible();
+
+        std::ostream& display(std::ostream& out, ineq const& c, bool values = false) const;
+        std::ostream& display(std::ostream& out, arg_t const& c, bool values = false) const;
+        virtual void display(std::ostream& out) const;
+        void display_resolved_lemma(std::ostream& out) const;
+
+        void add_clause(ineq& c, literal_vector const& lits);
+        void add_assign(ineq& c, literal_vector const& lits, literal l);
+        literal_vector& get_lits();
+
+        literal_vector& get_all_literals(ineq& c, bool negate);
+        literal_vector& get_helpful_literals(ineq& c, bool negate);
+        literal_vector& get_unhelpful_literals(ineq& c, bool negate);
+
+        //
+        // Utilities to compile cardinality 
+        // constraints into a sorting network.
+        //
+        void compile_ineq(ineq& c);
+        void inc_propagations(ineq& c);
+        unsigned get_compilation_threshold(ineq& c);
+
+        //
+        // Conflict resolution, cutting plane derivation.
+        // 
+        unsigned          m_num_marks;
+        unsigned          m_conflict_lvl;
+        arg_t             m_lemma;
+        literal_vector    m_ineq_literals;
+        svector<bool_var> m_marked;
+
+        // bool_var |-> index into m_lemma
+        unsigned_vector   m_conseq_index;
+        static const unsigned null_index;
+        bool is_marked(bool_var v) const;
+        void set_mark(bool_var v, unsigned idx);
+        void unset_mark(bool_var v);
+        void unset_marks();
+
+        bool resolve_conflict(ineq& c);
+        void process_antecedent(literal l, numeral coeff);
+        void process_ineq(ineq& c, literal conseq, numeral coeff);
+        void remove_from_lemma(unsigned idx);
+        bool is_proof_justification(justification const& j) const;
+
+        void hoist_maximal_values();
+
+        void validate_final_check();
+        void validate_final_check(ineq& c);
+        void validate_assign(ineq const& c, literal_vector const& lits, literal l) const;
+        void validate_watch(ineq const& c) const;
+
+        bool proofs_enabled() const { return get_manager().proofs_enabled(); }
+        justification* justify(literal l1, literal l2);
+        justification* justify(literal_vector const& lits);
+
+    public:
+        theory_pb(ast_manager& m, theory_pb_params& p);
+        
+        virtual ~theory_pb();
+
+        virtual theory * mk_fresh(context * new_ctx);
+        virtual bool internalize_atom(app * atom, bool gate_ctx);
+        virtual bool internalize_term(app * term) { UNREACHABLE(); return false; }
+        virtual void new_eq_eh(theory_var v1, theory_var v2);
+        virtual void new_diseq_eh(theory_var v1, theory_var v2) { }
+        virtual bool use_diseqs() const { return false; }
+        virtual bool build_models() const { return false; }
+        virtual final_check_status final_check_eh();
+        virtual void reset_eh();
+        virtual void assign_eh(bool_var v, bool is_true);
+        virtual void init_search_eh();
+        virtual void push_scope_eh();
+        virtual void pop_scope_eh(unsigned num_scopes);
+        virtual void restart_eh();
+        virtual void collect_statistics(::statistics & st) const;
+        virtual model_value_proc * mk_value(enode * n, model_generator & mg);
+        virtual void init_model(model_generator & m);        
+
+        static literal assert_ge(context& ctx, unsigned k, unsigned n, literal const* xs);
+    };
+};

src/smt/theory_utvpi.h

@@ -138,8 +138,7 @@ namespace smt {
         smt_params              m_params;
         arith_util              a;
         arith_eq_adapter        m_arith_eq_adapter;
-        th_var                  m_zero_int; // cache the variable representing the zero variable.
-        th_var                  m_zero_real; // cache the variable representing the zero variable.
+        th_var                  m_zero; //cache the variable representing the zero variable.
 
         dl_graph<GExt>          m_graph;
         nc_functor              m_nc_functor;
@@ -304,7 +303,7 @@ namespace smt {
 
         void new_eq_or_diseq(bool is_eq, th_var v1, th_var v2, justification& eq_just);
         
-        th_var get_zero(sort* s) const { return a.is_int(s)?m_zero_int:m_zero_real; }
+        th_var get_zero(sort* s) const { return m_zero; }
 
         th_var get_zero(expr* e) const { return get_zero(get_manager().get_sort(e)); }
 

src/smt/theory_utvpi_def.h

@@ -62,8 +62,7 @@ namespace smt {
             theory(m.mk_family_id("arith")),
             a(m),
             m_arith_eq_adapter(*this, m_params, a),
-            m_zero_int(null_theory_var),
-            m_zero_real(null_theory_var),
+            m_zero(null_theory_var),
             m_nc_functor(*this),
             m_asserted_qhead(0),
             m_agility(0.5),
@@ -134,8 +133,7 @@ namespace smt {
     template<typename Ext>
     void theory_utvpi<Ext>::reset_eh() {
         m_graph            .reset();
-        m_zero_int          = null_theory_var;
-        m_zero_real         = null_theory_var;
+        m_zero              = null_theory_var;
         m_atoms            .reset();
         m_asserted_atoms   .reset();
         m_stats            .reset();
@@ -261,8 +259,7 @@ namespace smt {
     template<typename Ext>
     void theory_utvpi<Ext>::init(context* ctx) {
         theory::init(ctx);
-        m_zero_int  = mk_var(ctx->mk_enode(a.mk_numeral(rational(0), true), false, false, true));
-        m_zero_real = mk_var(ctx->mk_enode(a.mk_numeral(rational(0), false), false, false, true));
+        m_zero  = mk_var(ctx->mk_enode(a.mk_numeral(rational(0), true), false, false, true));
     }
 
     /**
@@ -556,7 +553,7 @@ namespace smt {
         theory_var v = null_theory_var;
         context& ctx = get_context();
         if (r.is_zero()) {            
-            v = a.is_int(n)?m_zero_int:m_zero_real;
+            v = m_zero;
         }
         else if (ctx.e_internalized(n)) {
             enode* e = ctx.get_enode(n);
@@ -778,9 +775,7 @@ namespace smt {
         m_factory = alloc(arith_factory, get_manager());
         m.register_factory(m_factory);
         enforce_parity();
-        m_graph.set_to_zero(to_var(m_zero_int), to_var(m_zero_real));
-        m_graph.set_to_zero(neg(to_var(m_zero_int)), neg(to_var(m_zero_real)));
-        m_graph.set_to_zero(to_var(m_zero_int), neg(to_var(m_zero_int)));
+        m_graph.set_to_zero(to_var(m_zero), neg(to_var(m_zero)));
         compute_delta();   
         DEBUG_CODE(validate_model(););
     }

src/smt/theory_wmaxsat.cpp

@@ -0,0 +1,311 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    theory_wmaxsat.h
+
+Abstract:
+
+    Weighted Max-SAT theory plugin.
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2013-11-05
+
+Notes:
+
+--*/
+
+#include <typeinfo>
+#include "smt_context.h"
+#include "ast_pp.h"
+#include "theory_wmaxsat.h"
+
+namespace smt {
+
+theory_wmaxsat::theory_wmaxsat(ast_manager& m, filter_model_converter& mc):
+    theory(m.mk_family_id("weighted_maxsat")),
+    m_mc(mc),
+    m_vars(m),
+    m_fmls(m),
+    m_zweights(m_mpz),
+    m_old_values(m_mpz),
+    m_zcost(m_mpz),
+    m_zmin_cost(m_mpz),
+    m_found_optimal(false),
+    m_propagate(false),
+    m_normalize(false)
+{}
+
+theory_wmaxsat::~theory_wmaxsat() { 
+    m_old_values.reset();
+}
+
+/**
+   \brief return the complement of variables that are currently assigned.
+*/
+void theory_wmaxsat::get_assignment(svector<bool>& result) {
+    result.reset();
+    
+    if (!m_found_optimal) {
+        for (unsigned i = 0; i < m_vars.size(); ++i) {
+            result.push_back(false);
+        }
+    }
+    else {
+        std::sort(m_cost_save.begin(), m_cost_save.end());
+        for (unsigned i = 0, j = 0; i < m_vars.size(); ++i) {
+            if (j < m_cost_save.size() && m_cost_save[j] == i) {
+                result.push_back(false);
+                ++j;
+            }
+            else {
+                result.push_back(true);
+            }
+        }
+    }
+    TRACE("opt",
+          tout << "cost save: ";
+          for (unsigned i = 0; i < m_cost_save.size(); ++i) {
+              tout << m_cost_save[i] << " ";
+          }
+          tout << "\nvars: ";
+          for (unsigned i = 0; i < m_vars.size(); ++i) {
+              tout << mk_pp(m_vars[i].get(), get_manager()) << " ";
+          }
+          tout << "\nassignment: ";
+          for (unsigned i = 0; i < result.size(); ++i) {
+              tout << result[i] << " ";
+          }
+          tout << "\n";);                  
+    
+}
+
+void theory_wmaxsat::init_search_eh() {
+    m_propagate = true;
+}
+
+bool_var theory_wmaxsat::assert_weighted(expr* fml, rational const& w) {
+    context & ctx = get_context();
+    ast_manager& m = get_manager();
+    app_ref var(m), wfml(m);
+    var = m.mk_fresh_const("w", m.mk_bool_sort());
+    m_mc.insert(var->get_decl());
+    wfml = m.mk_or(var, fml);
+    ctx.assert_expr(wfml);
+    m_rweights.push_back(w);
+    m_vars.push_back(var);
+    m_fmls.push_back(fml);
+    m_assigned.push_back(false);
+    m_rmin_cost += w;
+    m_normalize = true;
+    return register_var(var, true);
+}
+
+bool_var theory_wmaxsat::register_var(app* var, bool attach) {
+    context & ctx = get_context();
+    bool_var bv;
+    SASSERT(!ctx.e_internalized(var));
+    enode* x = ctx.mk_enode(var, false, true, true);
+    if (ctx.b_internalized(var)) {
+        bv = ctx.get_bool_var(var);
+    }
+    else {
+        bv = ctx.mk_bool_var(var);
+    }
+    ctx.set_enode_flag(bv, true);
+    if (attach) {
+        ctx.set_var_theory(bv, get_id());                    
+        theory_var v = mk_var(x);
+        ctx.attach_th_var(x, this, v);
+        m_bool2var.insert(bv, v);
+        SASSERT(v == m_var2bool.size());
+        m_var2bool.push_back(bv);
+        SASSERT(ctx.bool_var2enode(bv));
+    }
+    return bv;
+}
+
+rational const& theory_wmaxsat::get_min_cost() { 
+    unsynch_mpq_manager mgr;
+    scoped_mpq q(mgr);
+    mgr.set(q, m_zmin_cost, m_den.to_mpq().numerator());
+    m_rmin_cost = rational(q);
+    return m_rmin_cost; 
+}
+
+void theory_wmaxsat::assign_eh(bool_var v, bool is_true) {
+    TRACE("opt", tout << "Assign " << mk_pp(m_vars[m_bool2var[v]].get(), get_manager()) << " " << is_true << "\n";);
+    if (is_true) {
+        if (m_normalize) normalize();
+        context& ctx = get_context();
+        theory_var tv = m_bool2var[v];
+        if (m_assigned[tv]) return;
+        scoped_mpz w(m_mpz);
+        w = m_zweights[tv];
+        ctx.push_trail(numeral_trail(m_zcost, m_old_values));
+        ctx.push_trail(push_back_vector<context, svector<theory_var> >(m_costs));
+        ctx.push_trail(value_trail<context, bool>(m_assigned[tv]));
+        m_zcost += w;
+        m_costs.push_back(tv);
+        m_assigned[tv] = true;
+        if (m_zcost > m_zmin_cost) {
+            block();
+        }
+    }
+}
+
+final_check_status theory_wmaxsat::final_check_eh() {
+    if (m_normalize) normalize();
+    return FC_DONE;
+}
+
+
+void theory_wmaxsat::reset_eh() {
+    theory::reset_eh();
+    reset_local();
+}
+
+void theory_wmaxsat::reset_local() {
+    m_vars.reset();
+    m_fmls.reset();
+    m_rweights.reset();
+    m_rmin_cost.reset();
+    m_rcost.reset();
+    m_zweights.reset();
+    m_zcost.reset();
+    m_zmin_cost.reset();
+    m_cost_save.reset();
+    m_bool2var.reset();
+    m_var2bool.reset();
+    m_propagate = false;
+    m_found_optimal = false;
+    m_assigned.reset();
+}
+
+
+void theory_wmaxsat::propagate() {
+    context& ctx = get_context();
+    for (unsigned i = 0; m_propagate && i < m_vars.size(); ++i) {
+        bool_var bv = m_var2bool[i];
+        lbool asgn = ctx.get_assignment(bv);
+        if (asgn == l_true) {
+            assign_eh(bv, true);
+        }
+    }
+    m_propagate = false;
+}       
+
+bool theory_wmaxsat::is_optimal() const {
+    return !m_found_optimal || m_zcost < m_zmin_cost;
+}
+
+expr_ref theory_wmaxsat::mk_block() {
+    ++m_stats.m_num_blocks;
+    ast_manager& m = get_manager();
+    expr_ref_vector disj(m);
+    compare_cost compare_cost(*this);
+    svector<theory_var> costs(m_costs);
+    std::sort(costs.begin(), costs.end(), compare_cost);
+    scoped_mpz weight(m_mpz);
+    m_mpz.reset(weight);
+    for (unsigned i = 0; i < costs.size() && m_mpz.lt(weight, m_zmin_cost); ++i) {
+        weight += m_zweights[costs[i]];
+        disj.push_back(m.mk_not(m_vars[costs[i]].get()));
+    }
+    if (is_optimal()) {
+        unsynch_mpq_manager mgr;
+        scoped_mpq q(mgr);
+        mgr.set(q, m_zmin_cost, m_den.to_mpq().numerator());
+        rational rw = rational(q);
+        m_zmin_cost = weight;
+        m_found_optimal = true;
+        m_cost_save.reset();
+        m_cost_save.append(m_costs);
+        TRACE("opt",
+              tout << "costs: ";
+              for (unsigned i = 0; i < m_costs.size(); ++i) {
+                  tout << mk_pp(get_enode(m_costs[i])->get_owner(), get_manager()) << " ";
+              }
+              tout << "\n";
+              get_context().display(tout);                      
+              );
+    }
+    expr_ref result(m.mk_or(disj.size(), disj.c_ptr()), m);
+    TRACE("opt",
+          tout << result << " weight: " << weight << "\n";
+          tout << "cost: " << m_zcost << " min-cost: " << m_zmin_cost << "\n";);
+    return result;
+}
+
+expr_ref theory_wmaxsat::mk_optimal_block(svector<bool_var> const& ws, rational const& weight) {
+    ast_manager& m = get_manager();
+    expr_ref_vector disj(m);
+    rational new_w = weight*m_den;
+    m_zmin_cost = new_w.to_mpq().numerator();
+    m_cost_save.reset();
+    for (unsigned i = 0; i < ws.size(); ++i) {
+        bool_var bv = ws[i];
+        theory_var v = m_bool2var[bv];
+        m_cost_save.push_back(v);
+        disj.push_back(m.mk_not(m_vars[v].get()));
+    }
+    expr_ref result(m.mk_or(disj.size(), disj.c_ptr()), m);
+    return result;
+}
+
+void theory_wmaxsat::block() {
+    if (m_vars.empty()) {
+        return;
+    }
+    ++m_stats.m_num_blocks;
+    ast_manager& m = get_manager();
+    context& ctx = get_context();
+    literal_vector lits;
+    compare_cost compare_cost(*this);
+    svector<theory_var> costs(m_costs);
+    std::sort(costs.begin(), costs.end(), compare_cost);
+    
+    scoped_mpz weight(m_mpz);
+    m_mpz.reset(weight);
+    for (unsigned i = 0; i < costs.size() && weight < m_zmin_cost; ++i) {
+        weight += m_zweights[costs[i]];
+        lits.push_back(~literal(m_var2bool[costs[i]]));
+    }
+    TRACE("opt",
+          tout << "block: ";
+          for (unsigned i = 0; i < lits.size(); ++i) {
+              expr_ref tmp(m);
+              ctx.literal2expr(lits[i], tmp);
+              tout << tmp << " ";
+          }
+          tout << "\n";
+          );
+    
+    ctx.mk_th_axiom(get_id(), lits.size(), lits.c_ptr());
+}                
+
+
+void theory_wmaxsat::normalize() {
+    m_den = rational::one();
+    for (unsigned i = 0; i < m_rweights.size(); ++i) {
+        m_den = lcm(m_den, denominator(m_rweights[i]));
+    }
+    m_den = lcm(m_den, denominator(m_rmin_cost));
+    SASSERT(!m_den.is_zero());
+    m_zweights.reset();
+    for (unsigned i = 0; i < m_rweights.size(); ++i) {
+        rational r = m_rweights[i]*m_den;
+        SASSERT(r.is_int());
+        mpq const& q = r.to_mpq();
+        m_zweights.push_back(q.numerator());
+    }
+    rational r = m_rcost* m_den;
+    m_zcost = r.to_mpq().numerator();
+    r = m_rmin_cost * m_den;
+    m_zmin_cost = r.to_mpq().numerator();
+    m_normalize = false;
+}
+
+};

src/smt/theory_wmaxsat.h

@@ -0,0 +1,128 @@
+/*++
+Copyright (c) 2013 Microsoft Corporation
+
+Module Name:
+
+    theory_wmaxsat.h
+
+Abstract:
+
+    Weighted Max-SAT theory plugin.
+
+Author:
+
+    Nikolaj Bjorner (nbjorner) 2013-11-05
+
+Notes:
+
+--*/
+
+#ifndef _THEORY_WMAXSAT_H_
+#define _THEORY_WMAXSAT_H_
+
+#include "smt_theory.h"
+#include "smt_clause.h"
+#include "filter_model_converter.h"
+
+namespace smt {
+    class theory_wmaxsat : public theory {
+        struct stats {
+            unsigned m_num_blocks;
+            void reset() { memset(this, 0, sizeof(*this)); }
+            stats() { reset(); }
+        };        
+        filter_model_converter&     m_mc;
+        mutable unsynch_mpz_manager m_mpz;
+        app_ref_vector           m_vars;        // Auxiliary variables per soft clause
+        expr_ref_vector          m_fmls;        // Formulas per soft clause
+        vector<rational>         m_rweights;    // weights of theory variables.
+        scoped_mpz_vector        m_zweights;
+        scoped_mpz_vector        m_old_values;
+        svector<theory_var>      m_costs;       // set of asserted theory variables
+        svector<theory_var>      m_cost_save;   // set of asserted theory variables
+        rational                 m_rcost;       // current sum of asserted costs
+        rational                 m_rmin_cost;   // current maximal cost assignment.
+        scoped_mpz               m_zcost;       // current sum of asserted costs
+        scoped_mpz               m_zmin_cost;   // current maximal cost assignment.
+        bool                     m_found_optimal; 
+        u_map<theory_var>        m_bool2var;    // bool_var -> theory_var
+        svector<bool_var>        m_var2bool;    // theory_var -> bool_var
+        bool                     m_propagate;  
+        bool                     m_normalize; 
+        rational                 m_den;         // lcm of denominators for rational weights.
+        svector<bool>            m_assigned;
+        stats                    m_stats;
+    public:
+        theory_wmaxsat(ast_manager& m, filter_model_converter& mc);
+        virtual ~theory_wmaxsat();
+        void get_assignment(svector<bool>& result);
+        virtual void init_search_eh();
+        bool_var assert_weighted(expr* fml, rational const& w);
+        bool_var register_var(app* var, bool attach);
+        rational const& get_min_cost();
+        class numeral_trail : public trail<context> {
+            typedef scoped_mpz T;
+            T & m_value;
+            scoped_mpz_vector&  m_old_values;            
+        public:
+            numeral_trail(T & value, scoped_mpz_vector& old):
+                m_value(value),
+                m_old_values(old) {
+                old.push_back(value);
+            }
+            
+            virtual ~numeral_trail() {
+            }
+            
+            virtual void undo(context & ctx) {
+                m_value = m_old_values.back();
+                m_old_values.shrink(m_old_values.size() - 1);
+            }
+        };
+        virtual void assign_eh(bool_var v, bool is_true);
+        virtual final_check_status final_check_eh();
+        virtual bool use_diseqs() const { 
+            return false;
+        }
+        virtual bool build_models() const { 
+            return false;
+        }
+        void reset_local();
+        virtual void reset_eh();
+        virtual theory * mk_fresh(context * new_ctx) { return 0; }
+        virtual bool internalize_atom(app * atom, bool gate_ctx) { return false; }
+        virtual bool internalize_term(app * term) { return false; }
+        virtual void new_eq_eh(theory_var v1, theory_var v2) { }
+        virtual void new_diseq_eh(theory_var v1, theory_var v2) { }
+
+        virtual void collect_statistics(::statistics & st) const {
+            st.update("wmaxsat num blocks", m_stats.m_num_blocks);
+        }
+        virtual bool can_propagate() {
+            return m_propagate;
+        }
+
+        virtual void propagate();
+        bool is_optimal() const;
+        expr_ref mk_block();
+
+        expr_ref mk_optimal_block(svector<bool_var> const& ws, rational const& weight);
+    private:
+
+        void block();
+        void normalize();
+       
+        class compare_cost {
+            theory_wmaxsat& m_th;
+        public:
+            compare_cost(theory_wmaxsat& t):m_th(t) {}
+            bool operator() (theory_var v, theory_var w) const { 
+                return m_th.m_mpz.gt(m_th.m_zweights[v], m_th.m_zweights[w]); 
+            }
+        };
+
+
+    };
+};
+
+#endif

src/smt/union_find.h

@@ -1,160 +0,0 @@
-/*++
-Copyright (c) 2006 Microsoft Corporation
-
-Module Name:
-
-    union_find.h
-
-Abstract:
-
-    <abstract>
-
-Author:
-
-    Leonardo de Moura (leonardo) 2008-05-31.
-
-Revision History:
-
---*/
-#ifndef _UNION_FIND_H_
-#define _UNION_FIND_H_
-
-#include "trail.h"
-
-class union_find_default_ctx {
-public:
-    typedef trail_stack<union_find_default_ctx> _trail_stack;
-    union_find_default_ctx() : m_stack(*this) {}
-       
-    void unmerge_eh(unsigned, unsigned) {}
-    void merge_eh(unsigned, unsigned, unsigned, unsigned) {}
-    void after_merge_eh(unsigned, unsigned, unsigned, unsigned) {}
-
-    _trail_stack& get_trail_stack() { return m_stack; }
-
-private:
-    _trail_stack m_stack;
-};
-
-template<typename Ctx = union_find_default_ctx>
-class union_find {
-    Ctx &                         m_ctx;
-    trail_stack<Ctx> &            m_trail_stack;
-    svector<unsigned>             m_find;
-    svector<unsigned>             m_size;
-    svector<unsigned>             m_next;
-    
-    class mk_var_trail;
-    friend class mk_var_trail;
-
-    class mk_var_trail : public trail<Ctx> {
-        union_find & m_owner;
-    public:
-        mk_var_trail(union_find & o):m_owner(o) {}
-        virtual ~mk_var_trail() {}
-        virtual void undo(Ctx & ctx) {
-            m_owner.m_find.pop_back();
-            m_owner.m_size.pop_back();
-            m_owner.m_next.pop_back();
-        }
-    };
-
-    mk_var_trail                  m_mk_var_trail;
-
-    class merge_trail;
-    friend class merge_trail;
-
-    class merge_trail : public trail<Ctx> {
-        union_find & m_owner;
-        unsigned     m_r1;
-    public:
-        merge_trail(union_find & o, unsigned r1):m_owner(o), m_r1(r1) {}
-        virtual ~merge_trail() {}
-        virtual void undo(Ctx & ctx) { m_owner.unmerge(m_r1); }
-    };
-
-    void unmerge(unsigned r1) {
-        unsigned r2 = m_find[r1];
-        TRACE("union_find", tout << "unmerging " << r1 << " " << r2 << "\n";);
-        SASSERT(find(r2) == r2);
-        m_size[r2] -= m_size[r1];
-        m_find[r1]  = r1;
-        std::swap(m_next[r1], m_next[r2]);
-        m_ctx.unmerge_eh(r2, r1);
-        CASSERT("union_find", check_invariant());
-    }
-
-public:    
-    union_find(Ctx & ctx):m_ctx(ctx), m_trail_stack(ctx.get_trail_stack()), m_mk_var_trail(*this) {}
-
-    unsigned mk_var() {
-        unsigned r = m_find.size();
-        m_find.push_back(r);
-        m_size.push_back(1);
-        m_next.push_back(r);
-        m_trail_stack.push_ptr(&m_mk_var_trail);
-        return r;
-    }
-
-    unsigned get_num_vars() const { return m_find.size(); }
-
-    unsigned find(unsigned v) const {
-        while (true) {
-            unsigned new_v = m_find[v];
-            if (new_v == v)
-                return v;
-            v = new_v;
-        }
-    }
-
-    unsigned next(unsigned v) const { return m_next[v]; }
-
-    bool is_root(unsigned v) const { return m_find[v] == v; }
-
-    void merge(unsigned v1, unsigned v2) {
-        unsigned r1 = find(v1);
-        unsigned r2 = find(v2);
-        TRACE("union_find", tout << "merging " << r1 << " " << r2 << "\n";);
-        if (r1 == r2)
-            return;
-        if (m_size[r1] > m_size[r2])
-            std::swap(r1, r2);
-        m_ctx.merge_eh(r2, r1, v2, v1);
-        m_find[r1] = r2;
-        m_size[r2] += m_size[r1];
-        std::swap(m_next[r1], m_next[r2]);
-        m_trail_stack.push(merge_trail(*this, r1));
-        m_ctx.after_merge_eh(r2, r1, v2, v1);
-        CASSERT("union_find", check_invariant());
-    }
-
-    void display(std::ostream & out) const {
-        unsigned num = get_num_vars(); 
-        for (unsigned v = 0; v < num; v++) {
-            out << "v" << v << " --> v" << m_find[v] << "\n";
-        }
-    }
-
-#ifdef Z3DEBUG
-    bool check_invariant() const {
-        unsigned num = get_num_vars(); 
-        for (unsigned v = 0; v < num; v++) {
-            if (is_root(v)) {
-                unsigned curr = v;
-                unsigned sz   = 0;
-                do {
-                    SASSERT(find(curr) == v);
-                    sz++;
-                    curr = next(curr);
-                }
-                while (curr != v);
-                SASSERT(m_size[v] == sz);
-            }
-        }
-        return true;
-    }
-#endif
-};
-
-#endif /* _UNION_FIND_H_ */
-