1 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
2 /*                                                                       */
3 /*    This file is part of the HiGHS linear optimization suite           */
4 /*                                                                       */
5 /*    Written and engineered 2008-2020 at the University of Edinburgh    */
6 /*                                                                       */
7 /*    Available as open-source under the MIT License                     */
8 /*                                                                       */
9 /* * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * */
10 /**@file lp_data/HSimplex.h
11  * @brief
12  * @author Julian Hall, Ivet Galabova, Qi Huangfu and Michael Feldmeier
13  */
14 #ifndef SIMPLEX_HSIMPLEX_H_
15 #define SIMPLEX_HSIMPLEX_H_
16 
17 #include "lp_data/HighsModelObject.h"
18 
19 enum class LpAction {
20   DUALISE = 0,
21   PERMUTE,
22   SCALE,
23   NEW_COSTS,
24   NEW_BOUNDS,
25   NEW_BASIS,
26   NEW_COLS,
27   NEW_ROWS,
28   DEL_COLS,
29   DEL_ROWS,
30   DEL_ROWS_BASIS_OK,
31   SCALED_COL,
32   SCALED_ROW,
33   BACKTRACKING
34 };
35 
36 void setSimplexOptions(
37     HighsModelObject& highs_model_object  //!< Model object in which simplex
38                                           //!< options are to be set
39 );
40 
41 int initialiseSimplexLpBasisAndFactor(HighsModelObject& highs_model_object,
42                                       const bool only_from_known_basis = false);
43 
44 HighsStatus transition(HighsModelObject& highs_model_object  //!< Model object
45 );
46 
47 void setNonbasicFlag(const HighsLp& simplex_lp, vector<int>& nonbasicFlag,
48                      const HighsBasisStatus* col_status = NULL,
49                      const HighsBasisStatus* row_status = NULL);
50 
51 void setNonbasicMove(const HighsLp& simplex_lp, const HighsScale& scale,
52                      const bool have_highs_basis, const HighsBasis& basis,
53                      const bool have_highs_solution,
54                      const HighsSolution& solution,
55                      SimplexBasis& simplex_basis);
56 
57 void initialiseNonbasicWorkValue(const HighsLp& simplex_lp,
58                                  const SimplexBasis& simplex_basis,
59                                  HighsSimplexInfo& simplex_info);
60 
61 bool basisConditionOk(HighsModelObject& highs_model_object,
62                       const double tolerance);
63 
64 // Methods not requiring HighsModelObject
65 
66 bool dual_infeasible(const double value, const double lower, const double upper,
67                      const double dual, const double value_tolerance,
68                      const double dual_tolerance);
69 
70 void appendNonbasicColsToBasis(HighsLp& lp, HighsBasis& basis, int XnumNewCol);
71 void appendNonbasicColsToBasis(HighsLp& lp, SimplexBasis& basis,
72                                int XnumNewCol);
73 
74 void appendBasicRowsToBasis(HighsLp& lp, HighsBasis& basis, int XnumNewRow);
75 void appendBasicRowsToBasis(HighsLp& lp, SimplexBasis& basis, int XnumNewRow);
76 
77 void reportBasis(const HighsOptions options, const HighsLp& lp,
78                  const HighsBasis& basis);
79 void reportBasis(const HighsOptions options, const HighsLp& lp,
80                  const SimplexBasis& simplex_basis);
81 
82 void computeDualObjectiveValue(HighsModelObject& highs_model_object,
83                                int phase = 2);
84 
85 void computePrimalObjectiveValue(HighsModelObject& highs_model_object);
86 
87 int setSourceOutFmBd(const HighsModelObject& highs_model_object,
88                      const int column_out);
89 
90 #ifdef HiGHSDEV
91 void getPrimalValue(const HighsModelObject& highs_model_object,
92                     vector<double>& primal_value);
93 void analysePrimalObjectiveValue(const HighsModelObject& highs_model_object);
94 #endif
95 
96 void initialiseSimplexLpDefinition(HighsModelObject& highs_model);
97 void initialiseSimplexLpRandomVectors(HighsModelObject& highs_model);
98 void extendSimplexLpRandomVectors(HighsModelObject& highs_model_object,
99                                   const int num_new_col, const int num_new_row);
100 
101 // SCALE:
102 
103 void scaleHighsModelInit(HighsModelObject& highs_model);
104 
105 void scaleCosts(HighsModelObject& highs_model);
106 
107 void scaleFactorRanges(HighsModelObject& highs_model_object,
108                        double& min_col_scale, double& max_col_scale,
109                        double& min_row_scale, double& max_row_scale);
110 
111 void scaleSimplexLp(HighsModelObject& highs_model);
112 bool equilibrationScaleMatrix(HighsModelObject& highs_model);
113 bool maxValueScaleMatrix(HighsModelObject& highs_model);
114 
115 HighsStatus deleteScale(const HighsOptions& options, vector<double>& scale,
116                         const HighsIndexCollection& index_collection);
117 // PERMUTE:
118 
119 void permuteSimplexLp(HighsModelObject& highs_model);
120 
121 #ifdef HiGHSDEV
122 // Only used to analyse the row and column status after Crash
123 void initialise_basic_index(HighsModelObject& highs_model_object);
124 #endif
125 
126 void allocateWorkAndBaseArrays(HighsModelObject& highs_model_object);
127 
128 void initialiseValueAndNonbasicMove(HighsModelObject& highs_model_object);
129 
130 void initialisePhase2ColBound(HighsModelObject& highs_model_object);
131 
132 void initialisePhase2RowBound(HighsModelObject& highs_model_object);
133 
134 void initialiseBound(HighsModelObject& highs_model_object, int phase = 2);
135 
136 void initialisePhase2ColCost(HighsModelObject& highs_model_object);
137 
138 void initialisePhase2RowCost(HighsModelObject& highs_model_object);
139 
140 void initialiseCost(HighsModelObject& highs_model_object, int perturb = 0);
141 
142 #ifdef HiGHSDEV
143 void reportSimplexProfiling(HighsModelObject& highs_model_object);
144 
145 #endif
146 void setRunQuiet(HighsModelObject& highs_model_object);
147 /**
148  * @brief Get the Hager condition number estimate for the basis matrix of a
149  * model
150  */
151 double computeBasisCondition(const HighsModelObject& highs_model_object);
152 
153 void flip_bound(HighsModelObject& highs_model_object, int iCol);
154 
155 void simplexHandleRankDeficiency(HighsModelObject& highs_model_object);
156 
157 int computeFactor(HighsModelObject& highs_model_object);
158 
159 // Compute the primal values (in baseValue) and set the lower and upper bounds
160 // of basic variables
161 void computePrimal(HighsModelObject& highs_model_object);
162 
163 void computeSimplexInfeasible(HighsModelObject& highs_model_object);
164 void computeSimplexPrimalInfeasible(HighsModelObject& highs_model_object);
165 void computeSimplexDualInfeasible(HighsModelObject& highs_model_object);
166 
167 void computeDualInfeasibleWithFlips(HighsModelObject& highs_model_object);
168 
169 void computeSimplexLpDualInfeasible(HighsModelObject& highs_model_object);
170 
171 void copySimplexInfeasible(HighsModelObject& highs_model_object);
172 void copySimplexDualInfeasible(HighsModelObject& highs_model_object);
173 void copySimplexPrimalInfeasible(HighsModelObject& highs_model_object);
174 
175 void choosePriceTechnique(const int price_strategy, const double row_ep_density,
176                           bool& use_col_price, bool& use_row_price_w_switch);
177 
178 void computeTableauRowFromPiP(HighsModelObject& highs_model_object,
179                               const HVector& row_ep, HVector& row_ap);
180 
181 void computeDual(HighsModelObject& highs_model_object);
182 
183 void correctDual(HighsModelObject& highs_model_object,
184                  int* free_infeasibility_count);
185 void correctDual(HighsModelObject& highs_model_object);
186 
187 // Record the shift in the cost of a particular column
188 void shift_cost(HighsModelObject& highs_model_object, int iCol, double amount);
189 
190 // Undo the shift in the cost of a particular column
191 void shift_back(HighsModelObject& highs_model_object, int iCol);
192 
193 // The major model updates. Factor calls factor.update; Matrix
194 // calls matrix.update; updatePivots does everything---and is
195 // called from the likes of HDual::updatePivots
196 void update_factor(HighsModelObject& highs_model_object, HVector* column,
197                    HVector* row_ep, int* iRow, int* hint);
198 
199 void update_pivots(HighsModelObject& highs_model_object, int columnIn,
200                    int rowOut, int sourceOut);
201 
202 void update_matrix(HighsModelObject& highs_model_object, int columnIn,
203                    int columnOut);
204 
205 bool reinvertOnNumericalTrouble(const std::string method_name,
206                                 HighsModelObject& highs_model_object,
207                                 double& numerical_trouble_measure,
208                                 const double alpha_from_col,
209                                 const double alpha_from_row,
210                                 const double numerical_trouble_tolerance);
211 
212 // Analyse the unscaled solution from a Simplex basic solution to get
213 // suggested feasibility tolerances for resolving the scaled LP
214 // This sets highs_model_object.unscaled_solution_params_
215 HighsStatus getNewInfeasibilityTolerancesFromSimplexBasicSolution(
216     const HighsModelObject& highs_model_object,
217     HighsSolutionParams& get_unscaled_solution_params,
218     double& new_scaled_primal_feasibility_tolerance,
219     double& new_scaled_dual_feasibility_tolerance);
220 
221 HighsStatus getInfeasibilitiesAndNewTolerances(
222     const HighsOptions& options, const HighsLp& lp, const HighsScale& scale,
223     const SimplexBasis& basis, const HighsSimplexInfo& simplex_info,
224     const HighsModelStatus scaled_model_status,
225     const HighsSolutionParams& unscaled_solution_params,
226     const HighsSolutionParams& scaled_solution_params,
227     HighsSolutionParams& get_unscaled_solution_params,
228     HighsSolutionParams& get_scaled_solution_params,
229     double& new_scaled_primal_feasibility_tolerance,
230     double& new_scaled_dual_feasibility_tolerance);
231 
232 void logRebuild(HighsModelObject& highs_model_object, const bool primal,
233                 const int solve_phase);
234 
235 void reportSimplexLpStatus(
236     HighsSimplexLpStatus&
237         simplex_lp_status,    // !< Status of simplex LP to be reported
238     const char* message = ""  // !< Message to be written in report
239 );
240 
241 void invalidateSimplexLpBasisArtifacts(
242     HighsSimplexLpStatus&
243         simplex_lp_status  // !< Status of simplex LP whose
244                            // basis artifacts are to be invalidated
245 );
246 
247 void invalidateSimplexLpBasis(
248     HighsSimplexLpStatus& simplex_lp_status  // !< Status of simplex LP whose
249                                              // basis is to be invalidated
250 );
251 
252 void invalidateSimplexLp(
253     HighsSimplexLpStatus&
254         simplex_lp_status  // !< Status of simplex LP to be invalidated
255 );
256 
257 void updateSimplexLpStatus(
258     HighsSimplexLpStatus&
259         simplex_lp_status,  // !< Status of simplex LP to be updated
260     LpAction action         // !< Action prompting update
261 );
262 
263 bool isBasisRightSize(const HighsLp& lp, const SimplexBasis& basis);
264 #endif  // SIMPLEX_HSIMPLEX_H_
265