Now applying user scaling to incumbent solution before/after solve

jajhall · jajhall · commit a8a840d3f524 · 2025-09-22T13:12:52.000+01:00
diff --git a/check/TestUserScale.cpp b/check/TestUserScale.cpp
@@ -7,71 +7,51 @@
 const bool dev_run = true;  // false;
 const double inf = kHighsInf;
 
-void checkModelScaling(const HighsInt user_bound_scale,
-                       const HighsInt user_cost_scale,
-                       const HighsModel& unscaled_model,
-                       const HighsModel& scaled_model);
-
-void checkLpScaling(const HighsInt user_bound_scale,
-                    const HighsInt user_cost_scale, const HighsLp& unscaled_lp,
-                    const HighsLp& scaled_lp);
-
 bool doubleEqual0(const double v0, const double v1) {
-  return std::fabs(v0 - v1)/std::max(1.0, std::fabs(v0)) < 1e-12; 
+  double rel_difference = std::fabs(v0 - v1) / std::max(1.0, std::fabs(v0));
+  bool ok = rel_difference < 1e-12;
+  if (dev_run && !ok)
+    printf("UserScaleDoubleEqual: %g and %g have relative difference = %g\n",
+           v0, v1, rel_difference);
+  return ok;
 }
 
 TEST_CASE("user-scale-after-run", "[highs_user_scale]") {
-  const std::string mip_model = "flugpl";//"rgn";//
+  const std::string mip_model = "flugpl";  //"rgn";//
   std::string model = "avgas";
   Highs highs;
   const HighsInfo& info = highs.getInfo();
   highs.setOptionValue("output_flag", dev_run);
-  //    REQUIRE(highs.setOptionValue("presolve", kHighsOffString) == HighsStatus::kOk);
-  HighsInt num_k = 1;
+  //    REQUIRE(highs.setOptionValue("presolve", kHighsOffString) ==
+  //    HighsStatus::kOk);
+  HighsInt num_k = 2;
   if (num_k == 1) model = mip_model;
   for (HighsInt k = 0; k < num_k; k++) {
     std::string filename =
-      std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps";
+        std::string(HIGHS_DIR) + "/check/instances/" + model + ".mps";
     highs.readModel(filename);
     HighsLp unscaled_lp = highs.getLp();
     const bool is_lp = !unscaled_lp.isMip();
-    highs.writeModel("unscaled.mps");
+    //    highs.writeModel("unscaled.mps");
 
     highs.run();
     double unscaled_objective = highs.getInfo().objective_function_value;
 
     HighsInt user_bound_scale = 1;
     double user_bound_scale_value = std::pow(2, user_bound_scale);
-    
-    REQUIRE(highs.setOptionValue("user_bound_scale", user_bound_scale) == HighsStatus::kOk);
-
-    HighsLp bound_scaled_lp = highs.getLp();
-    checkLpScaling(user_bound_scale, 0, unscaled_lp, bound_scaled_lp);
+    REQUIRE(highs.setOptionValue("user_bound_scale", user_bound_scale) ==
+            HighsStatus::kOk);
 
     HighsInt user_cost_scale = 4;
     double user_cost_scale_value = std::pow(2, user_cost_scale);
-    REQUIRE(highs.setOptionValue("user_cost_scale", user_cost_scale) == HighsStatus::kOk);
-
-    HighsLp scaled_lp = highs.getLp();
-    checkLpScaling(user_bound_scale, user_cost_scale, unscaled_lp, scaled_lp);
-
-    REQUIRE(highs.getModelStatus() == HighsModelStatus::kNotset);
-    REQUIRE(info.primal_solution_status == kSolutionStatusNone);
-    REQUIRE(info.dual_solution_status == kSolutionStatusNone);
-    REQUIRE(info.num_primal_infeasibilities == kHighsIllegalInfeasibilityCount);
-    REQUIRE(info.max_primal_infeasibility == kHighsIllegalInfeasibilityMeasure);
-    REQUIRE(info.sum_primal_infeasibilities == kHighsIllegalInfeasibilityMeasure);
-    REQUIRE(info.num_dual_infeasibilities == kHighsIllegalInfeasibilityCount);
-    REQUIRE(info.max_dual_infeasibility == kHighsIllegalInfeasibilityMeasure);
-    REQUIRE(info.sum_dual_infeasibilities == kHighsIllegalInfeasibilityMeasure);
-
-    highs.writeModel("scaled0.mps");
-    highs.passModel(scaled_lp);
+    REQUIRE(highs.setOptionValue("user_cost_scale", user_cost_scale) ==
+            HighsStatus::kOk);
+
     highs.run();
-    
+
     REQUIRE(highs.getModelStatus() == HighsModelStatus::kOptimal);
     REQUIRE(doubleEqual0(highs.getInfo().objective_function_value,
-			 user_bound_scale_value*user_cost_scale_value*unscaled_objective));
+                         unscaled_objective));
 
     model = mip_model;
     REQUIRE(highs.setOptionValue("user_bound_scale", 0) == HighsStatus::kOk);
@@ -81,38 +61,6 @@ TEST_CASE("user-scale-after-run", "[highs_user_scale]") {
   highs.resetGlobalScheduler(true);
 }
 
-TEST_CASE("user-scale-after-load", "[highs_user_scale]") {
-  std::string filename =
-      std::string(HIGHS_DIR) + "/check/instances/adlittle.mps";
-  Highs highs;
-  const HighsInfo& info = highs.getInfo();
-  highs.setOptionValue("output_flag", dev_run);
-
-  for (HighsInt k = 0; k < 2; k++) {
-    highs.readModel(filename);
-    HighsLp unscaled_lp = highs.getLp();
-
-    HighsInt user_bound_scale = 10;
-    REQUIRE(highs.setOptionValue("user_bound_scale", user_bound_scale) == HighsStatus::kOk);
-
-    HighsInt user_cost_scale = 30;
-    REQUIRE(highs.setOptionValue("user_cost_scale", user_cost_scale) == HighsStatus::kOk);
-
-    REQUIRE(highs.readModel(filename) == HighsStatus::kOk);
-    HighsLp scaled_lp = highs.getLp();
-
-    checkLpScaling(user_bound_scale, user_cost_scale, unscaled_lp, scaled_lp);
-
-    filename =
-      std::string(HIGHS_DIR) + "/check/instances/flugpl.mps";
-    REQUIRE(highs.setOptionValue("user_bound_scale", 0) == HighsStatus::kOk);
-    REQUIRE(highs.setOptionValue("user_cost_scale", 0) == HighsStatus::kOk);
-    
-  }
-
-  highs.resetGlobalScheduler(true);
-}
-
 TEST_CASE("user-small-cost-scale", "[highs_user_scale]") {
   Highs highs;
   const HighsInfo& info = highs.getInfo();
@@ -151,53 +99,3 @@ TEST_CASE("user-small-cost-scale", "[highs_user_scale]") {
 
   highs.resetGlobalScheduler(true);
 }
-
-void checkModelScaling(const HighsInt user_bound_scale,
-                       const HighsInt user_cost_scale,
-                       const HighsModel& unscaled_model,
-                       const HighsModel& scaled_model) {
-  checkLpScaling(user_bound_scale, user_cost_scale, unscaled_model.lp_,
-                 scaled_model.lp_);
-}
-
-void checkLpScaling(const HighsInt user_bound_scale,
-                    const HighsInt user_cost_scale, const HighsLp& unscaled_lp,
-                    const HighsLp& scaled_lp) {
-  const double user_bound_scale_value = std::pow(2, user_bound_scale);
-  const double user_cost_scale_value = std::pow(2, user_cost_scale);
-  REQUIRE(unscaled_lp.num_col_ == scaled_lp.num_col_);
-  REQUIRE(unscaled_lp.num_row_ == scaled_lp.num_row_);
-  const bool has_integrality = scaled_lp.integrality_.size() > 0;
-  const HighsSparseMatrix& unscaled_matrix = unscaled_lp.a_matrix_;
-  REQUIRE(unscaled_matrix.isColwise());
-  for (HighsInt iCol = 0; iCol < unscaled_lp.num_col_; iCol++) {
-    bool continuous = !has_integrality || scaled_lp.integrality_[iCol] == HighsVarType::kContinuous;
-    double value = unscaled_lp.col_cost_[iCol] * user_cost_scale_value;
-    if (!continuous) value *= user_bound_scale_value;
-    REQUIRE(scaled_lp.col_cost_[iCol] == value);
-    if (unscaled_lp.col_lower_[iCol] > -inf) {
-      value = unscaled_lp.col_lower_[iCol];
-      if (continuous) value *= user_bound_scale_value;
-      REQUIRE(scaled_lp.col_lower_[iCol] == value);
-    }
-    if (unscaled_lp.col_upper_[iCol] < inf) {
-      value = unscaled_lp.col_upper_[iCol];
-      if (continuous) value *= user_bound_scale_value;
-      REQUIRE(scaled_lp.col_upper_[iCol] == value);
-    }
-    for (HighsInt iEl = unscaled_matrix.start_[iCol]; iEl < unscaled_matrix.start_[iCol+1]; iEl++) {
-      value = unscaled_matrix.value_[iEl];
-      if (!continuous) value *= user_bound_scale_value;
-      REQUIRE(scaled_lp.a_matrix_.value_[iEl] == value);
-    }
-  }
-  for (HighsInt iRow = 0; iRow < unscaled_lp.num_row_; iRow++) {
-    if (unscaled_lp.row_lower_[iRow] > -inf)
-      REQUIRE(scaled_lp.row_lower_[iRow] ==
-              unscaled_lp.row_lower_[iRow] * user_bound_scale_value);
-    if (unscaled_lp.row_upper_[iRow] < inf)
-      REQUIRE(scaled_lp.row_upper_[iRow] ==
-              unscaled_lp.row_upper_[iRow] * user_bound_scale_value);
-  }
-}
-
diff --git a/highs/Highs.h b/highs/Highs.h
@@ -1709,7 +1709,8 @@ class Highs {
   HighsStatus optionChangeAction() const { return HighsStatus::kOk; };
 
   HighsStatus userScaleModel(HighsUserScaleData& data);
-  HighsStatus userScaleSolution(HighsUserScaleData& data, bool update_kkt = false);
+  HighsStatus userScaleSolution(HighsUserScaleData& data,
+                                bool update_kkt = false);
   HighsStatus computeIllConditioning(HighsIllConditioning& ill_conditioning,
                                      const bool constraint,
                                      const HighsInt method,
diff --git a/highs/lp_data/Highs.cpp b/highs/lp_data/Highs.cpp
@@ -946,10 +946,10 @@ HighsStatus Highs::run() {
   HighsUserScaleData user_scale_data;
   initialiseUserScaleData(this->options_, user_scale_data);
   const bool user_scaling =
-    user_scale_data.user_cost_scale ||
-    user_scale_data.user_bound_scale;
+      user_scale_data.user_cost_scale || user_scale_data.user_bound_scale;
   if (user_scaling) {
-    if (this->userScaleModel(user_scale_data) == HighsStatus::kError) return HighsStatus::kError;
+    if (this->userScaleModel(user_scale_data) == HighsStatus::kError)
+      return HighsStatus::kError;
     this->userScaleSolution(user_scale_data);
     // Zero the user scale values to prevent further scaling
     this->options_.user_cost_scale = 0;
@@ -977,13 +977,14 @@ HighsStatus Highs::run() {
   if (user_scaling) {
     // Unscale the incumbent model and solution
     //
-    // Flip the scaling sign 
+    // Flip the scaling sign
     user_scale_data.user_cost_scale *= -1;
     user_scale_data.user_bound_scale *= -1;
     HighsStatus unscale_status = this->userScaleModel(user_scale_data);
     if (unscale_status == HighsStatus::kError) {
-      highsLogUser(this->options_.log_options, HighsLogType::kError,
-		   "Unexpected error removing user scaling from the incumbent model\n");
+      highsLogUser(
+          this->options_.log_options, HighsLogType::kError,
+          "Unexpected error removing user scaling from the incumbent model\n");
       assert(unscale_status != HighsStatus::kError);
     }
     const bool update_kkt = true;
@@ -992,6 +993,10 @@ HighsStatus Highs::run() {
     // negated to undo user scaling
     this->options_.user_cost_scale = -user_scale_data.user_cost_scale;
     this->options_.user_bound_scale = -user_scale_data.user_bound_scale;
+    highsLogUser(this->options_.log_options, HighsLogType::kInfo,
+                 "After solving the user-scaled model the unscaled solution "
+                 "has objective value %.12g\n",
+                 this->info_.objective_function_value);
   }
   return status;
 }
diff --git a/highs/lp_data/HighsInterface.cpp b/highs/lp_data/HighsInterface.cpp
@@ -2994,7 +2994,8 @@ HighsStatus Highs::userScaleModel(HighsUserScaleData& data) {
   return return_status;
 }
 
-HighsStatus Highs::userScaleSolution(HighsUserScaleData& data, bool update_kkt) {
+HighsStatus Highs::userScaleSolution(HighsUserScaleData& data,
+                                     bool update_kkt) {
   HighsStatus return_status = HighsStatus::kOk;
   if (!data.user_cost_scale && !data.user_bound_scale) return HighsStatus::kOk;
   double cost_scale_value = std::pow(2, data.user_cost_scale);
@@ -3004,23 +3005,25 @@ HighsStatus Highs::userScaleSolution(HighsUserScaleData& data, bool update_kkt)
   if (info_.primal_solution_status != kSolutionStatusNone) {
     if (data.user_bound_scale) {
       for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++) {
-	if (has_integrality && lp.integrality_[iCol] != HighsVarType::kContinuous) continue;
-	this->solution_.col_value[iCol] *= bound_scale_value;
+        if (has_integrality &&
+            lp.integrality_[iCol] != HighsVarType::kContinuous)
+          continue;
+        this->solution_.col_value[iCol] *= bound_scale_value;
       }
-      for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++) 
-	this->solution_.row_value[iRow] *= bound_scale_value;
+      for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++)
+        this->solution_.row_value[iRow] *= bound_scale_value;
     }
   }
   if (info_.dual_solution_status != kSolutionStatusNone) {
     if (data.user_cost_scale) {
-      for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++) 
-	this->solution_.col_dual[iCol] *= cost_scale_value;
-      for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++) 
-	this->solution_.row_dual[iRow] *= cost_scale_value;
+      for (HighsInt iCol = 0; iCol < lp.num_col_; iCol++)
+        this->solution_.col_dual[iCol] *= cost_scale_value;
+      for (HighsInt iRow = 0; iRow < lp.num_row_; iRow++)
+        this->solution_.row_dual[iRow] *= cost_scale_value;
     }
   }
   if (!update_kkt) return return_status;
-  info_.objective_function_value *= (bound_scale_value*cost_scale_value);
+  info_.objective_function_value *= (bound_scale_value * cost_scale_value);
   getKktFailures(options_, model_, solution_, basis_, info_);
   return return_status;
 }
diff --git a/highs/lp_data/HighsLpUtils.cpp b/highs/lp_data/HighsLpUtils.cpp
@@ -769,8 +769,7 @@ HighsStatus cleanBounds(const HighsOptions& options, HighsLp& lp) {
   return HighsStatus::kOk;
 }
 
-HighsStatus userScaleLp(HighsLp& lp,
-                        HighsUserScaleData& data,
+HighsStatus userScaleLp(HighsLp& lp, HighsUserScaleData& data,
                         const HighsLogOptions& log_options) {
   userScaleLp(lp, data, false);
   HighsStatus return_status = userScaleStatus(log_options, data);
@@ -779,10 +778,10 @@ HighsStatus userScaleLp(HighsLp& lp,
   return return_status;
 }
 
-void userScaleLp(HighsLp& lp, HighsUserScaleData& data,
-                 const bool apply) {
+void userScaleLp(HighsLp& lp, HighsUserScaleData& data, const bool apply) {
   userScaleCosts(lp.integrality_, lp.col_cost_, data, apply);
-  userScaleColBounds(lp.integrality_, lp.col_lower_, lp.col_upper_, data, apply);
+  userScaleColBounds(lp.integrality_, lp.col_lower_, lp.col_upper_, data,
+                     apply);
   userScaleMatrix(lp.integrality_, lp.a_matrix_, data, apply);
   userScaleRowBounds(lp.row_lower_, lp.row_upper_, data, apply);
 }
@@ -1734,10 +1733,10 @@ HighsStatus changeLpIntegrality(HighsLp& lp,
   HighsInt lp_col;
   HighsInt usr_col = -1;
 
- // If changing integrality for a problem without an integrality
- // vector (ie an LP), have to create it for the incumbent columns -
- // which are naturally continuous
-  if (lp.integrality_.size() == 0) 
+  // If changing integrality for a problem without an integrality
+  // vector (ie an LP), have to create it for the incumbent columns -
+  // which are naturally continuous
+  if (lp.integrality_.size() == 0)
     lp.integrality_.assign(lp.num_col_, HighsVarType::kContinuous);
 
   assert(HighsInt(lp.integrality_.size()) == lp.num_col_);
@@ -3688,8 +3687,9 @@ void HighsUserScaleData::initialise(const HighsInt& user_cost_scale_,
 }
 
 bool HighsUserScaleData::scaleError(std::string& message) const {
-  if (this->num_infinite_costs + this->num_infinite_hessian_values + this->num_infinite_col_bounds +
-          this->num_infinite_row_bounds + this->num_large_matrix_values ==
+  if (this->num_infinite_costs + this->num_infinite_hessian_values +
+          this->num_infinite_col_bounds + this->num_infinite_row_bounds +
+          this->num_large_matrix_values ==
       0)
     return false;
   assert(this->user_cost_scale != 0 || this->user_bound_scale != 0);
@@ -3710,19 +3710,18 @@ bool HighsUserScaleData::scaleError(std::string& message) const {
   }
   if (this->num_infinite_hessian_values) {
     if (this->num_infinite_costs) {
-      if (this->num_infinite_col_bounds ||
-	  this->num_infinite_row_bounds) {
+      if (this->num_infinite_col_bounds || this->num_infinite_row_bounds) {
         ss << ",";
       } else {
         ss << " and";
       }
     }
-    ss << " " << this->num_infinite_hessian_values << " infinite Hessian values";
+    ss << " " << this->num_infinite_hessian_values
+       << " infinite Hessian values";
     if (this->num_infinite_hessian_values > 1) ss << "s";
   }
   if (this->num_infinite_col_bounds) {
-    if (this->num_infinite_costs ||
-	this->num_infinite_hessian_values) {
+    if (this->num_infinite_costs || this->num_infinite_hessian_values) {
       if (this->num_infinite_row_bounds) {
         ss << ",";
       } else {
@@ -3733,17 +3732,15 @@ bool HighsUserScaleData::scaleError(std::string& message) const {
     if (this->num_infinite_col_bounds > 1) ss << "s";
   }
   if (this->num_infinite_row_bounds) {
-    if (this->num_infinite_costs ||
-	this->num_infinite_hessian_values ||
-	this->num_infinite_col_bounds) ss << " and";
+    if (this->num_infinite_costs || this->num_infinite_hessian_values ||
+        this->num_infinite_col_bounds)
+      ss << " and";
     ss << " " << this->num_infinite_row_bounds << " infinite row bound";
     if (this->num_infinite_row_bounds > 1) ss << "s";
   }
   if (this->num_large_matrix_values) {
-    if (this->num_infinite_costs +
-	this->num_infinite_hessian_values +
-	this->num_infinite_col_bounds +
-	this->num_infinite_row_bounds >
+    if (this->num_infinite_costs + this->num_infinite_hessian_values +
+            this->num_infinite_col_bounds + this->num_infinite_row_bounds >
         0)
       ss << ", and";
     ss << " " << this->num_large_matrix_values << " large matrix value";
diff --git a/highs/lp_data/HighsLpUtils.h b/highs/lp_data/HighsLpUtils.h
@@ -61,12 +61,10 @@ HighsStatus assessBounds(const HighsOptions& options, const char* type,
 
 HighsStatus cleanBounds(const HighsOptions& options, HighsLp& lp);
 
-HighsStatus userScaleLp(HighsLp& lp,
-                        HighsUserScaleData& data,
+HighsStatus userScaleLp(HighsLp& lp, HighsUserScaleData& data,
                         const HighsLogOptions& log_options);
 
-void userScaleLp(HighsLp& lp,
-		 HighsUserScaleData& data,
+void userScaleLp(HighsLp& lp, HighsUserScaleData& data,
                  const bool apply = true);
 
 void userScaleCosts(const vector<HighsVarType>& integrality,
diff --git a/highs/lp_data/HighsModelUtils.cpp b/highs/lp_data/HighsModelUtils.cpp
diff --git a/highs/lp_data/HighsSolve.h b/highs/lp_data/HighsSolve.h
diff --git a/highs/model/HighsHessianUtils.cpp b/highs/model/HighsHessianUtils.cpp
