Improve repairSolution handling of numerical issues

Fixed highspy unit test that had assumed unique solution

Author: mathgeekcoder

highs/lp_data/HighsCallback.cpp

@@ -9,6 +9,7 @@
  * @brief
  */
 #include "lp_data/HighsCallback.h"
+#include "util/HighsIntegers.h"
 
 #include <algorithm>
 #include <cassert>
@@ -254,10 +255,30 @@ HighsStatus HighsCallbackInput::repairSolution() {
     clone.setOptionValue("output_flag", false);
     clone.passModel(highs->getModel());
 
+    HighsVarType vtype = HighsVarType::kContinuous;
+    double tolerance = highs->getOptions().mip_feasibility_tolerance;
+
     // fix the variables
     for (HighsInt i = 0; i < user_solution.size(); i++) {
       if (user_solution[i] != kHighsUndefined) {
-        clone.changeColBounds(i, user_solution[i], user_solution[i]);
+        double value = user_solution[i];
+        highs->getColIntegrality(i, vtype);
+
+        if (vtype == HighsVarType::kInteger) {
+          // check if the provided solution value is integerfeasible
+          if (!HighsIntegers::isIntegral(value, tolerance)) {
+            highsLogUser(highs->getOptions().log_options, HighsLogType::kError,
+                         "repairSolution: User solution (index %d) is outside "
+                         "integrality tolerance (value %f)\n",
+                         i, value);
+            return HighsStatus::kError;
+          } else {
+            // since the variable is integral, round it to avoid numerical issues
+            value = std::round(value);
+          }
+        }
+
+        clone.changeColBounds(i, value, value);
       }
     }
 

tests/test_highspy.py

@@ -1454,13 +1454,13 @@ def partial_solution(e):
 
             # get every 2nd element from 2d numpy sol array
             self.assertEqual(e.data_in.setSolution(x[:, ::2], sol[:, ::2]), highspy.HighsStatus.kOk)
-            self.assertTrue((e.data_in.user_solution[0:4] == [0,highspy.kHighsUndefined,0,highspy.kHighsUndefined]).all())
+            self.assertEqual(list(e.data_in.user_solution[0:4]), [sol[0,0],highspy.kHighsUndefined,sol[0,2],highspy.kHighsUndefined])
             self.assertEqual(e.data_in.repairSolution(), highspy.HighsStatus.kOk)
-            self.assertEqual(list(e.data_in.user_solution[0:8]), [0., 0., 0., 0., 0., 1., 0., 0.])
+            self.assertEqual(list(e.data_in.user_solution[0:8]), list(sol[0,0:8]))
 
             def try_change_ptr(e):
                 e.data_in.user_solution = [0] * (N*N)
-        
+
             self.assertRaises(Exception, lambda: try_change_ptr(e))
 
             # modify directly
@@ -1475,6 +1475,10 @@ def try_change_ptr(e):
             self.assertEqual(list(e.data_in.user_solution[0:8]), [0., 0., 0., 0., 0., 0., 1., 0.])
             self.assertEqual(e.data_in.user_has_solution, True)
 
+            # set partial solution with fractional value
+            self.assertEqual(e.data_in.setSolution([0., 0.5]), highspy.HighsStatus.kOk)
+            self.assertEqual(e.data_in.repairSolution(), highspy.HighsStatus.kError)
+
         # verify partial solution
         h.cbMipUserSolution += partial_solution