factor out model-based-opt code

Signed-off-by: Nikolaj Bjorner <nbjorner@microsoft.com>

contrib/cmake/src/test/CMakeLists.txt

@@ -61,6 +61,7 @@ add_executable(test-z3
   "${CMAKE_CURRENT_BINARY_DIR}/mem_initializer.cpp"
   memory.cpp
   model2expr.cpp
+  model_based_opt.cpp
   model_evaluator.cpp
   model_retrieval.cpp
   mpbq.cpp

src/math/simplex/model_based_opt.cpp

@@ -23,6 +23,13 @@ Revision History:
 namespace opt {
     
 
+    model_based_opt::model_based_opt():
+        m_objective_id(0)
+    {
+        m_rows.push_back(row());
+    }
+
+        
     bool model_based_opt::invariant() {
         // variables in each row are sorted.
         for (unsigned i = 0; i < m_rows.size(); ++i) {
@@ -30,12 +37,13 @@ namespace opt {
                 return false;
             }
         }
-        return invariant(m_objective);
+        return true;
     }
 
     bool model_based_opt::invariant(row const& r) {
         rational val = r.m_coeff;
         vector<var> const& vars = r.m_vars;
+        SASSERT(!vars.empty());
         for (unsigned i = 0; i < vars.size(); ++i) {
             var const& v = vars[i];
             SASSERT(i + 1 == vars.size() || v.m_id < vars[i+1].m_id);
@@ -72,12 +80,10 @@ namespace opt {
     // then replace a3*x + t3 by t3/a3 - t2/a2 <= 0
     // 
     bound_type model_based_opt::maximize(rational& value) {
-        // tbd
         SASSERT(invariant());
-        vector<var> & vars = m_objective.m_vars;
         unsigned_vector other;
-        while (!vars.empty()) {
+        while (!objective().m_vars.empty()) {
-            var const& v = vars.back();
+            var const& v = objective().m_vars.back();
             unsigned x = v.m_id;
             rational const& coeff = v.m_coeff;
             rational const& x_val = m_var2value[x];
@@ -90,23 +96,25 @@ namespace opt {
                     resolve(other[i], bound_coeff, bound_index, x);
                 }
                 // coeff*x + objective -> coeff*(bound) + objective
-                // tbd:
                 multiply(coeff/bound_coeff, bound_index);
-                //add(m_objective_id, bound_index);
+                SASSERT(invariant(m_rows[bound_index]));
+                objective().m_vars.back().m_coeff.reset();
+                add(m_objective_id, bound_index);
+                SASSERT(invariant(objective()));
                 m_rows[bound_index].m_alive = false;
             }
             else {
                 return unbounded;
             }
         }
-        value = m_objective.m_coeff;
+        value = objective().m_coeff;
-        switch (m_objective.m_type) {
+        if (objective().m_type == t_lt) {
-        case t_lt: return strict;
+            return strict;
-        case t_le: return non_strict;
-        case t_eq: return non_strict;
         }
+        else {
             return non_strict;
         }
+    }
 
     bool model_based_opt::find_bound(unsigned x, unsigned& bound_index, unsigned_vector& other, bool is_pos) {
         bound_index = UINT_MAX;
@@ -143,6 +151,7 @@ namespace opt {
         
     rational model_based_opt::get_coefficient(unsigned row_id, unsigned var_id) {
         row const& r = m_rows[row_id];
+        SASSERT(!r.m_vars.empty());
         unsigned lo = 0, hi = r.m_vars.size();
         while (lo < hi) {
             unsigned mid = lo + (hi - lo)/2;
@@ -228,9 +237,11 @@ namespace opt {
             else if (i == r1.m_vars.size()) {
                 for (; j < r2.m_vars.size(); ++j) {
                     m_new_vars.push_back(r2.m_vars[j]);
+                    if (row_id1 != m_objective_id) {
                         m_var2row_ids[r2.m_vars[j].m_id].push_back(row_id1);
                     }
                 }
+            }
             else {
                 unsigned v1 = r1.m_vars[i].m_id;
                 unsigned v2 = r2.m_vars[j].m_id;
@@ -249,7 +260,9 @@ namespace opt {
                 }
                 else {
                     m_new_vars.push_back(r2.m_vars[j]);
+                    if (row_id1 != m_objective_id) {
                         m_var2row_ids[r2.m_vars[j].m_id].push_back(row_id1);
+                    }
                     ++j;                        
                 }
             }                
@@ -260,6 +273,7 @@ namespace opt {
         if (r2.m_type == t_lt) {
             r1.m_type = t_lt;
         }
+        SASSERT(invariant(r1));
     }
     
     void model_based_opt::display(std::ostream& out) const {
@@ -291,16 +305,41 @@ namespace opt {
     }
 
     unsigned model_based_opt::add_var(rational const& value) {
-        NOT_IMPLEMENTED_YET();
+        unsigned v = m_var2value.size();
-        return 0;
+        m_var2value.push_back(value);
+        m_var2row_ids.push_back(unsigned_vector());
+        return v;
     }
 
-    void model_based_opt::add_constraint(vector<var> const& coeffs, rational const& c, ineq_type r) {
+    void model_based_opt::set_row(row& r, vector<var> const& coeffs, rational const& c, ineq_type rel) {
-        NOT_IMPLEMENTED_YET();
+        rational val(c);
+        SASSERT(r.m_vars.empty());
+        r.m_vars.append(coeffs.size(), coeffs.c_ptr());
+        std::sort(r.m_vars.begin(), r.m_vars.end(), var::compare());
+        for (unsigned i = 0; i < coeffs.size(); ++i) {
+            val += m_var2value[coeffs[i].m_id] * coeffs[i].m_coeff;
+        }
+        r.m_alive = true;
+        r.m_coeff = c;
+        r.m_value = val;
+        r.m_type = rel;
+        SASSERT(invariant(r));
+    }
+
+    void model_based_opt::add_constraint(vector<var> const& coeffs, rational const& c, ineq_type rel) {
+        rational val(c);
+        row r0;
+        unsigned row_id = m_rows.size();
+        m_rows.push_back(r0);
+        row& r = m_rows.back();
+        set_row(r, coeffs, c, rel);
+        for (unsigned i = 0; i < coeffs.size(); ++i) {
+            m_var2row_ids[coeffs[i].m_id].push_back(row_id); 
+        }
     }
 
     void model_based_opt::set_objective(vector<var> const& coeffs, rational const& c) {
-        NOT_IMPLEMENTED_YET();
+        set_row(objective(), coeffs, c, t_le);
     }
 
 }

src/math/simplex/model_based_opt.h

@@ -44,24 +44,33 @@ namespace opt {
             unsigned m_id;
             rational m_coeff;
             var(unsigned id, rational const& c): m_id(id), m_coeff(c) {}
+            struct compare {
+                bool operator()(var x, var y) {
+                    return x.m_id < y.m_id;
+                }
+            };
         };
     private:
         struct row {
+            row(): m_type(t_le), m_value(0), m_alive(false) {}
             vector<var> m_vars;         // variables with coefficients
             rational    m_coeff;        // constant in inequality
             ineq_type   m_type;         // inequality type
             rational    m_value;        // value of m_vars + m_coeff under interpretation of m_var2value.
             bool        m_alive;        // rows can be marked dead if they have been processed.
         };
+
         vector<row>             m_rows;
+        unsigned                m_objective_id;
         vector<unsigned_vector> m_var2row_ids;
         vector<rational>        m_var2value;
-        row                     m_objective;
         vector<var>             m_new_vars;
         
         bool invariant();
         bool invariant(row const& r);
 
+        row& objective() { return m_rows[0]; }
+        
 
         bool find_bound(unsigned x, unsigned& bound_index, unsigned_vector& other, bool is_pos);
         
@@ -73,8 +82,12 @@ namespace opt {
 
         void add(unsigned row_id1, unsigned row_id2);
 
+        void set_row(row& r, vector<var> const& coeffs, rational const& c, ineq_type rel);
+
     public:
 
+        model_based_opt();
+
         // add a fresh variable with value 'value'.
         unsigned add_var(rational const& value);
 

src/test/main.cpp

@@ -186,6 +186,7 @@ int main(int argc, char ** argv) {
     TST(smt_context);
     TST(theory_dl);
     TST(model_retrieval);
+    TST(model_based_opt);
     TST(factor_rewriter);
     TST(smt2print_parse);
     TST(substitution);

src/test/model_based_opt.cpp

@@ -0,0 +1,49 @@
+#include "model_based_opt.h"
+
+static void test1() {
+    opt::model_based_opt mbo;
+    typedef opt::model_based_opt::var var;
+    vector<var> vars;
+    unsigned x = mbo.add_var(rational(2));
+    unsigned y = mbo.add_var(rational(3));
+    unsigned z = mbo.add_var(rational(4));
+    unsigned u = mbo.add_var(rational(5));
+
+    vars.reset();
+    vars.push_back(var(x, rational(1)));
+    vars.push_back(var(y, rational(-1)));
+    mbo.add_constraint(vars, rational(0), opt::t_le);
+
+    vars.reset();
+    vars.push_back(var(x, rational(1)));
+    vars.push_back(var(z, rational(-1)));
+    mbo.add_constraint(vars, rational(0), opt::t_le);
+
+    vars.reset();
+    vars.push_back(var(y, rational(1)));
+    vars.push_back(var(u, rational(-1)));
+    mbo.add_constraint(vars, rational(0), opt::t_le);
+
+    vars.reset();
+    vars.push_back(var(z, rational(1)));
+    vars.push_back(var(u, rational(-1)));
+    mbo.add_constraint(vars, rational(-1), opt::t_le);
+
+    vars.reset();
+    vars.push_back(var(u, rational(1)));
+    mbo.add_constraint(vars, rational(4), opt::t_le);
+
+    vars.reset();
+    vars.push_back(var(x, rational(2)));
+    mbo.set_objective(vars, rational(0));
+
+    rational value;
+    opt::bound_type bound = mbo.maximize(value);
+    
+    std::cout << bound << ": " << value << "\n";
+    
+}
+
+void tst_model_based_opt() {
+    test1();
+}