From eaf52f4c3265d7dcc85354f697ef1cac684a57bf Mon Sep 17 00:00:00 2001
From: Nikolaj Bjorner <nbjorner@microsoft.com>
Date: Tue, 18 Oct 2022 14:52:30 -0700
Subject: [PATCH] fix infinite loop issue in def::operator+, other issues
 remain though for #6404

The example from #6404 results in an incorrect result. It uses integer division on private variables where MBQI support is new and not tested for substitutions.
---
 src/math/simplex/model_based_opt.cpp | 34 ++++++++++++++++++++++++----
 src/math/simplex/model_based_opt.h   |  1 +
 2 files changed, 31 insertions(+), 4 deletions(-)

diff --git a/src/math/simplex/model_based_opt.cpp b/src/math/simplex/model_based_opt.cpp
index b90851582..7ab246c6e 100644
--- a/src/math/simplex/model_based_opt.cpp
+++ b/src/math/simplex/model_based_opt.cpp
@@ -98,11 +98,13 @@ namespace opt {
             }
             else if (v1 < v2) {
                 vs.push_back(vs1[i]);
-                vs.back().m_coeff *= c1;                 
+                vs.back().m_coeff *= c1;    
+                ++i;
             }
             else {
                 vs.push_back(vs2[j]);
-                vs.back().m_coeff *= c2;                 
+                vs.back().m_coeff *= c2;   
+                ++j;
             }
         }
         result.m_div = c1*m_div;
@@ -111,6 +113,26 @@ namespace opt {
         return result;
     }
 
+    model_based_opt::def model_based_opt::def::substitute(unsigned v, def const& other) const {
+        def result;
+        vector<var> const& vs1 = m_vars;
+        vector<var> const& vs2 = other.m_vars;
+        rational coeff(0);
+        for (auto const& [id, c] : vs1) {
+            if (id == v) {
+                coeff = c;
+                break;
+            }
+        }
+        if (coeff == 0) {
+            return *this;
+        }
+
+        NOT_IMPLEMENTED_YET();       
+        result.normalize();
+        return result;
+    }
+
     model_based_opt::def model_based_opt::def::operator/(rational const& r) const {
         def result(*this);
         result.m_div *= r;
@@ -1418,13 +1440,17 @@ namespace opt {
             project(v, false);
 
         // project internal variables.
-
-        def y_def = project(y, compute_def);
         def z_def = project(z, compute_def);
+        def y_def = project(y, compute_def); // may depend on z
+
 
         if (compute_def) {
+
             result = (y_def * K) + z_def;
             m_var2value[x] = eval(result);
+            TRACE("opt", tout << y << " := " << y_def << "\n";
+                         tout << z << " := " << z_def << "\n";
+                         tout << x << " := " << result << "\n");
         }
         TRACE("opt", display(tout << "solve_div done v" << x << "\n"));
         return result;
diff --git a/src/math/simplex/model_based_opt.h b/src/math/simplex/model_based_opt.h
index 5ba9bd619..8150f945f 100644
--- a/src/math/simplex/model_based_opt.h
+++ b/src/math/simplex/model_based_opt.h
@@ -86,6 +86,7 @@ namespace opt {
             def operator/(rational const& n) const;
             def operator*(rational const& n) const;
             def operator+(rational const& n) const;
+            def substitute(unsigned v, def const& other) const;
             void normalize();
         };