Add option for ordering

filikat · filikat · commit b3df7d74e403 · 2025-11-29T16:45:09.000Z
diff --git a/docs/src/options/definitions.md b/docs/src/options/definitions.md
@@ -434,6 +434,11 @@
 - Type: string
 - Default: "both"
 
+## hipo\_ordering
+- HiPO matrix reordering: "metis", "amd", "rcm" or "choose"
+- Type: string
+- Default: "choose"
+
 ## hipo\_block\_size
 - Block size for dense linear algebra within HiPO
 - Type: integer
diff --git a/highs/ipm/IpxWrapper.cpp b/highs/ipm/IpxWrapper.cpp
@@ -552,6 +552,18 @@ HighsStatus solveLpHipo(const HighsOptions& options, HighsTimer& timer,
     return HighsStatus::kError;
   }
 
+  // Reordering heuristic
+  if (options.hipo_ordering != kHipoMetisString &&
+      options.hipo_ordering != kHipoAmdString &&
+      options.hipo_ordering != kHipoRcmString &&
+      options.hipo_ordering != kHighsChooseString) {
+    highsLogUser(options.log_options, HighsLogType::kError,
+                 "Unknown value of option %s\n", kHipoOrderingString.c_str());
+    model_status = HighsModelStatus::kSolveError;
+    return HighsStatus::kError;
+  }
+  hipo_options.ordering = options.hipo_ordering;
+
   // block size option
   hipo_options.block_size = options.hipo_block_size;
 
diff --git a/highs/ipm/hipo/factorhighs/Analyse.cpp b/highs/ipm/hipo/factorhighs/Analyse.cpp
@@ -13,6 +13,7 @@
 #include "ReturnValues.h"
 #include "ipm/hipo/auxiliary/Auxiliary.h"
 #include "ipm/hipo/auxiliary/Log.h"
+#include "rcm/rcm.h"
 
 namespace hipo {
 
@@ -21,8 +22,8 @@ const Int64 int64_limit = std::numeric_limits<int64_t>::max();
 
 Analyse::Analyse(const std::vector<Int>& rows, const std::vector<Int>& ptr,
                  const std::vector<Int>& signs, Int nb, const Log* log,
-                 DataCollector& data)
-    : log_{log}, data_{data} {
+                 DataCollector& data, const std::string& ordering)
+    : log_{log}, data_{data}, ordering_{ordering} {
   // Input the symmetric matrix to be analysed in CSC format.
   // rows contains the row indices.
   // ptr contains the starting points of each column.
@@ -116,46 +117,74 @@ Int Analyse::getPermutation(bool metis_no2hop) {
     }
   }
 
-  idx_t options[METIS_NOPTIONS];
-  METIS_SetDefaultOptions(options);
-  options[METIS_OPTION_SEED] = kMetisSeed;
+  if (ordering_ == "metis") {
+    // ----------------------------
+    // ----- METIS ----------------
+    // ----------------------------
+    idx_t options[METIS_NOPTIONS];
+    METIS_SetDefaultOptions(options);
+    options[METIS_OPTION_SEED] = kMetisSeed;
 
-  // set logging of Metis depending on debug level
-  options[METIS_OPTION_DBGLVL] = 0;
-  if (log_->debug(2))
-    options[METIS_OPTION_DBGLVL] = METIS_DBG_INFO | METIS_DBG_COARSEN;
+    // set logging of Metis depending on debug level
+    options[METIS_OPTION_DBGLVL] = 0;
+    if (log_->debug(2))
+      options[METIS_OPTION_DBGLVL] = METIS_DBG_INFO | METIS_DBG_COARSEN;
 
-  // set no2hop=1 if the user requested it
-  if (metis_no2hop) options[METIS_OPTION_NO2HOP] = 1;
+    // set no2hop=1 if the user requested it
+    if (metis_no2hop) options[METIS_OPTION_NO2HOP] = 1;
 
-  if (log_) log_->printDevInfo("Running Metis\n");
+    if (log_) log_->printDevInfo("Running Metis\n");
 
-  Int status = METIS_NodeND(&n_, temp_ptr.data(), temp_rows.data(), NULL,
-                            options, perm_.data(), iperm_.data());
+    Int status = METIS_NodeND(&n_, temp_ptr.data(), temp_rows.data(), NULL,
+                              options, perm_.data(), iperm_.data());
 
-  if (log_) log_->printDevInfo("Metis done\n");
-  if (status != METIS_OK) {
-    if (log_) log_->printDevInfo("Error with Metis\n");
-    return kRetMetisError;
+    if (log_) log_->printDevInfo("Metis done\n");
+    if (status != METIS_OK) {
+      if (log_) log_->printDevInfo("Error with Metis\n");
+      return kRetMetisError;
+    }
+
+  } else if (ordering_ == "amd") {
+    // ----------------------------
+    // ------ AMD -----------------
+    // ----------------------------
+    double control[AMD_CONTROL];
+    amd_defaults(control);
+    double info[AMD_INFO];
+
+    if (log_) log_->printDevInfo("Running AMD\n");
+    Int status = amd_order(n_, temp_ptr.data(), temp_rows.data(), perm_.data(),
+                           control, info);
+    if (log_) log_->printDevInfo("AMD done\n");
+
+    if (status != AMD_OK) {
+      if (log_) log_->printDevInfo("Error with AMD\n");
+      return kRetMetisError;
+    }
+    inversePerm(perm_, iperm_);
   }
 
-  /*
-  double control[AMD_CONTROL];
-  amd_defaults(control);
+  else if (ordering_ == "rcm") {
+    // ----------------------------
+    // ------ RCM -----------------
+    // ----------------------------
 
-  double info[AMD_INFO];
+    if (log_) log_->printDevInfo("Running RCM\n");
+    Int status = genrcm(n_, temp_ptr.back(), temp_ptr.data(), temp_rows.data(),
+                        perm_.data());
+    if (log_) log_->printDevInfo("RCM done\n");
 
-  Int status = amd_order(n_, temp_ptr.data(), temp_rows.data(), perm_.data(),
-                       control, info);
+    if (status != 0) {
+      if (log_) log_->printDevInfo("Error with RCM\n");
+      return kRetMetisError;
+    }
+    inversePerm(perm_, iperm_);
 
-  if (status != AMD_OK) {
-  if (log_) log_->printDevInfo("Error with AMD\n");
-  return kRetMetisError;
+  } else {
+    if (log_) log_->printe("Invalid reordering option\n");
+    return kRetMetisError;
   }
 
-  inversePerm(perm_, iperm_);
-  */
-
   return kRetOk;
 }
 
diff --git a/highs/ipm/hipo/factorhighs/Analyse.h b/highs/ipm/hipo/factorhighs/Analyse.h
@@ -85,6 +85,8 @@ class Analyse {
   const Log* log_;
   DataCollector& data_;
 
+  const std::string& ordering_;
+
   // Functions to perform analyse phase
   Int getPermutation(bool metis_no2hop);
   void permute(const std::vector<Int>& iperm);
@@ -110,7 +112,7 @@ class Analyse {
   // Constructor: matrix must be in lower triangular format
   Analyse(const std::vector<Int>& rows, const std::vector<Int>& ptr,
           const std::vector<Int>& signs, Int nb, const Log* log,
-          DataCollector& data);
+          DataCollector& data, const std::string& ordering);
 
   // Run analyse phase and save the result in Symbolic object S
   Int run(Symbolic& S);
diff --git a/highs/ipm/hipo/factorhighs/FactorHiGHS.cpp b/highs/ipm/hipo/factorhighs/FactorHiGHS.cpp
@@ -36,8 +36,9 @@ void FHsolver::setRegularisation(double reg_p, double reg_d) {
 
 Int FHsolver::analyse(Symbolic& S, const std::vector<Int>& rows,
                       const std::vector<Int>& ptr,
-                      const std::vector<Int>& signs) {
-  Analyse an_obj(rows, ptr, signs, nb_, log_, data_);
+                      const std::vector<Int>& signs,
+                      const std::string& ordering) {
+  Analyse an_obj(rows, ptr, signs, nb_, log_, data_, ordering);
   return an_obj.run(S);
 }
 
diff --git a/highs/ipm/hipo/factorhighs/FactorHiGHS.h b/highs/ipm/hipo/factorhighs/FactorHiGHS.h
@@ -83,7 +83,8 @@ class FHsolver {
   // ptr, and store symbolic factorisation in object S.
   // See ReturnValues.h for errors.
   Int analyse(Symbolic& S, const std::vector<Int>& rows,
-              const std::vector<Int>& ptr, const std::vector<Int>& signs);
+              const std::vector<Int>& ptr, const std::vector<Int>& signs,
+              const std::string& ordering = "metis");
 
   // Perform factorise phase of matrix given by rows, ptr, vals, and store
   // numerical factorisation in object N. Matrix is moved into the object, so
diff --git a/highs/ipm/hipo/ipm/FactorHiGHSSolver.cpp b/highs/ipm/hipo/ipm/FactorHiGHSSolver.cpp
@@ -398,8 +398,11 @@ Int FactorHiGHSSolver::analyseAS(Symbolic& S) {
 
   log_.printDevInfo("Performing AS analyse phase\n");
 
+  std::string ordering = options_.ordering;
+  if (ordering == kHighsChooseString) ordering = kHipoMetisString;
+
   clock.start();
-  Int status = FH_.analyse(S, rowsLower, ptrLower, pivot_signs);
+  Int status = FH_.analyse(S, rowsLower, ptrLower, pivot_signs, ordering);
   if (info_) info_->analyse_AS_time = clock.stop();
 
   if (status && log_.debug(2)) {
@@ -437,8 +440,11 @@ Int FactorHiGHSSolver::analyseNE(Symbolic& S, Int64 nz_limit) {
 
   log_.printDevInfo("Performing NE analyse phase\n");
 
+  std::string ordering = options_.ordering;
+  if (ordering == kHighsChooseString) ordering = kHipoMetisString;
+
   clock.start();
-  Int status = FH_.analyse(S, rowsNE_, ptrNE_, pivot_signs);
+  Int status = FH_.analyse(S, rowsNE_, ptrNE_, pivot_signs, ordering);
   if (info_) info_->analyse_NE_time = clock.stop();
 
   if (status && log_.debug(2)) {
diff --git a/highs/ipm/hipo/ipm/Options.h b/highs/ipm/hipo/ipm/Options.h
@@ -40,6 +40,7 @@ struct Options {
   OptionNla nla = kOptionNlaDefault;
   OptionCrossover crossover = kOptionCrossoverDefault;
   OptionParallel parallel = kOptionParallelDefault;
+  std::string ordering = kHighsChooseString;
 
   // Metis options
   bool metis_no2hop = false;
diff --git a/highs/lp_data/HighsOptions.h b/highs/lp_data/HighsOptions.h
@@ -307,6 +307,12 @@ const string kHipoTreeString = "tree";
 const string kHipoNodeString = "node";
 const string kHipoBothString = "both";
 
+// Strings for HiPO matrix reordering
+const string kHipoOrderingString = "hipo_ordering";
+const string kHipoMetisString = "metis";
+const string kHipoAmdString = "amd";
+const string kHipoRcmString = "rcm";
+
 struct HighsOptionsStruct {
   // Run-time options read from the command line
   std::string presolve;
@@ -371,6 +377,7 @@ struct HighsOptionsStruct {
   HighsInt ipm_iteration_limit;
   std::string hipo_system;
   std::string hipo_parallel_type;
+  std::string hipo_ordering;
   HighsInt hipo_block_size;
   bool hipo_metis_no2hop;
 
@@ -543,6 +550,7 @@ struct HighsOptionsStruct {
         ipm_iteration_limit(0),
         hipo_system(""),
         hipo_parallel_type(""),
+        hipo_ordering(""),
         hipo_block_size(0),
         hipo_metis_no2hop(false),
         pdlp_scaling(false),
@@ -1262,6 +1270,13 @@ class HighsOptions : public HighsOptionsStruct {
                                advanced, &hipo_parallel_type, kHipoBothString);
     records.push_back(record_string);
 
+    record_string =
+        new OptionRecordString(kHipoOrderingString,
+                               "HiPO matrix reordering option: \"metis\", "
+                               "\"amd\", \"rcm\" or \"choose\".",
+                               advanced, &hipo_ordering, kHighsChooseString);
+    records.push_back(record_string);
+
     record_int = new OptionRecordInt(
         "hipo_block_size", "Block size for dense linear algebra within HiPO",
         advanced, &hipo_block_size, 0, 128, kHighsIInf);
