1 // Copyright 2010-2021 Google LLC
2 // Licensed under the Apache License, Version 2.0 (the "License");
3 // you may not use this file except in compliance with the License.
4 // You may obtain a copy of the License at
5 //
6 //     http://www.apache.org/licenses/LICENSE-2.0
7 //
8 // Unless required by applicable law or agreed to in writing, software
9 // distributed under the License is distributed on an "AS IS" BASIS,
10 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
11 // See the License for the specific language governing permissions and
12 // limitations under the License.
13 
14 // Common types and constants used by the Linear Programming solver.
15 
16 #ifndef OR_TOOLS_LP_DATA_LP_TYPES_H_
17 #define OR_TOOLS_LP_DATA_LP_TYPES_H_
18 
19 #include <cmath>
20 #include <cstdint>
21 #include <limits>
22 
23 #include "ortools/base/basictypes.h"
DeferredMessage()24 #include "ortools/base/int_type.h"
25 #include "ortools/base/logging.h"
26 #include "ortools/base/strong_vector.h"
27 #include "ortools/util/bitset.h"
28 
29 // We use typedefs as much as possible to later permit the usage of
30 // types such as quad-doubles or rationals.
31 
32 namespace operations_research {
33 namespace glop {
34 
35 // This type is defined to avoid cast issues during index conversions,
36 // e.g. converting ColIndex into RowIndex.
37 // All types should use 'Index' instead of int32_t.
38 typedef int32_t Index;
39 
40 // ColIndex is the type for integers representing column/variable indices.
41 // int32s are enough for handling even the largest problems.
42 DEFINE_INT_TYPE(ColIndex, Index);
43 
44 // RowIndex is the type for integers representing row/constraint indices.
45 // int32s are enough for handling even the largest problems.
46 DEFINE_INT_TYPE(RowIndex, Index);
47 
48 // Get the ColIndex corresponding to the column # row.
49 inline ColIndex RowToColIndex(RowIndex row) { return ColIndex(row.value()); }
50 
51 // Get the RowIndex corresponding to the row # col.
52 inline RowIndex ColToRowIndex(ColIndex col) { return RowIndex(col.value()); }
53 
54 // Get the integer index corresponding to the col.
55 inline Index ColToIntIndex(ColIndex col) { return col.value(); }
56 
57 // Get the integer index corresponding to the row.
58 inline Index RowToIntIndex(RowIndex row) { return row.value(); }
59 
60 // EntryIndex is the type for integers representing entry indices.
61 // An entry in a sparse matrix is a pair (row, value) for a given known column.
62 // See classes SparseColumn and SparseMatrix.
63 #if defined(__ANDROID__)
64 DEFINE_INT_TYPE(EntryIndex, int32_t);
65 #else
66 DEFINE_INT_TYPE(EntryIndex, int64_t);
67 #endif
68 
69 static inline double ToDouble(double f) { return f; }
70 
71 static inline double ToDouble(long double f) { return static_cast<double>(f); }
72 
73 // The type Fractional denotes the type of numbers on which the computations are
74 // performed. This is defined as double here, but it could as well be float,
75 // DoubleDouble, QuadDouble, or infinite-precision rationals.
76 // Floating-point representations are binary fractional numbers, thus the name.
77 // (See http://en.wikipedia.org/wiki/Fraction_(mathematics) .)
78 typedef double Fractional;
79 
80 // Range max for type Fractional. DBL_MAX for double for example.
81 const double kRangeMax = std::numeric_limits<double>::max();
82 
83 // Infinity for type Fractional.
84 const double kInfinity = std::numeric_limits<double>::infinity();
85 
86 // Epsilon for type Fractional, i.e. the smallest e such that 1.0 + e != 1.0 .
87 const double kEpsilon = std::numeric_limits<double>::epsilon();
88 
89 // Returns true if the given value is finite, that means for a double:
90 // not a NaN and not +/- infinity.
91 inline bool IsFinite(Fractional value) {
92   return value > -kInfinity && value < kInfinity;
93 }
94 
95 // Constants to represent invalid row or column index.
96 // It is important that their values be the same because during transposition,
97 // one needs to be converted into the other.
98 const RowIndex kInvalidRow(-1);
99 const ColIndex kInvalidCol(-1);
100 
101 // Different statuses for a given problem.
102 enum class ProblemStatus : int8_t {
103   // The problem has been solved to optimality. Both the primal and dual have
104   // a feasible solution.
105   OPTIMAL,
106 
107   // The problem has been proven primal-infeasible. Note that the problem is not
108   // necessarily DUAL_UNBOUNDED (See Chvatal p.60). The solver does not have a
109   // dual unbounded ray in this case.
110   PRIMAL_INFEASIBLE,
111 
112   // The problem has been proven dual-infeasible. Note that the problem is not
113   // necessarily PRIMAL_UNBOUNDED (See Chvatal p.60). The solver does
114   // note have a primal unbounded ray in this case,
115   DUAL_INFEASIBLE,
116 
117   // The problem is either INFEASIBLE or UNBOUNDED (this applies to both the
118   // primal and dual algorithms). This status is only returned by the presolve
119   // step and means that a primal or dual unbounded ray was found during
120   // presolve. Note that because some presolve techniques assume that a feasible
121   // solution exists to simplify the problem further, it is difficult to
122   // distinguish between infeasibility and unboundedness.
123   //
124   // If a client needs to distinguish, it is possible to run the primal
125   // algorithm on the same problem with a 0 objective function to know if the
126   // problem was PRIMAL_INFEASIBLE.
127   INFEASIBLE_OR_UNBOUNDED,
128 
129   // The problem has been proven feasible and unbounded. That means that the
130   // problem is DUAL_INFEASIBLE and that the solver has a primal unbounded ray.
131   PRIMAL_UNBOUNDED,
132 
133   // The problem has been proven dual-feasible and dual-unbounded. That means
134   // the problem is PRIMAL_INFEASIBLE and that the solver has a dual unbounded
135   // ray to prove it.
136   DUAL_UNBOUNDED,
137 
138   // All the statuses below correspond to a case where the solver was
139   // interrupted. This can happen because of a timeout, an iteration limit or an
140   // error.
141 
142   // The solver didn't had a chance to prove anything.
143   INIT,
144 
145   // The problem has been proven primal-feasible but may still be
146   // PRIMAL_UNBOUNDED.
147   PRIMAL_FEASIBLE,
148 
149   // The problem has been proven dual-feasible, but may still be DUAL_UNBOUNDED.
150   // That means that if the primal is feasible, then it has a finite optimal
151   // solution.
152   DUAL_FEASIBLE,
153 
154   // An error occurred during the solving process.
155   ABNORMAL,
156 
157   // The input problem was invalid (see LinearProgram.IsValid()).
158   INVALID_PROBLEM,
159 
160   // The problem was solved to a feasible status, but the solution checker found
161   // the primal and/or dual infeasibilities too important for the specified
162   // parameters.
163   IMPRECISE,
164 };
165 
166 // Returns the string representation of the ProblemStatus enum.
167 std::string GetProblemStatusString(ProblemStatus problem_status);
168 
169 inline std::ostream& operator<<(std::ostream& os, ProblemStatus status) {
170   os << GetProblemStatusString(status);
171   return os;
172 }
173 
174 // Different types of variables.
175 enum class VariableType : int8_t {
176   UNCONSTRAINED,
177   LOWER_BOUNDED,
178   UPPER_BOUNDED,
179   UPPER_AND_LOWER_BOUNDED,
180   FIXED_VARIABLE
181 };
182 
183 // Returns the string representation of the VariableType enum.
184 std::string GetVariableTypeString(VariableType variable_type);
185 
186 inline std::ostream& operator<<(std::ostream& os, VariableType type) {
187   os << GetVariableTypeString(type);
188   return os;
189 }
190 
191 // Different variables statuses.
192 // If a solution is OPTIMAL or FEASIBLE, then all the properties described here
193 // should be satisfied. These properties should also be true during the
194 // execution of the revised simplex algorithm, except that because of
195 // bound-shifting, the variable may not be at their exact bounds until the
196 // shifts are removed.
197 enum class VariableStatus : int8_t {
198   // The basic status is special and takes precedence over all the other
199   // statuses. It means that the variable is part of the basis.
200   BASIC,
201   // Only possible status of a FIXED_VARIABLE not in the basis. The variable
202   // value should be exactly equal to its bounds (which are the same).
203   FIXED_VALUE,
204   // Only possible statuses of a non-basic variable which is not UNCONSTRAINED
205   // or FIXED. The variable value should be at its exact specified bound (which
206   // must be finite).
207   AT_LOWER_BOUND,
208   AT_UPPER_BOUND,
209   // Only possible status of an UNCONSTRAINED non-basic variable.
210   // Its value should be zero.
211   //
212   // Note that during crossover, this status is relaxed, and any variable that
213   // can currently move in both directions can be marked as free.
214   FREE,
215 };
216 
217 // Returns the string representation of the VariableStatus enum.
218 std::string GetVariableStatusString(VariableStatus status);
219 
220 inline std::ostream& operator<<(std::ostream& os, VariableStatus status) {
221   os << GetVariableStatusString(status);
222   return os;
223 }
224 
225 // Different constraints statuses.
226 // The meaning is the same for the constraint activity relative to its bounds as
227 // it is for a variable value relative to its bounds. Actually, this is the
228 // VariableStatus of the slack variable associated to a constraint modulo a
229 // change of sign. The difference is that because of precision error, a
230 // constraint activity cannot exactly be equal to one of its bounds or to zero.
231 enum class ConstraintStatus : int8_t {
232   BASIC,
233   FIXED_VALUE,
234   AT_LOWER_BOUND,
235   AT_UPPER_BOUND,
236   FREE,
237 };
238 
239 // Returns the string representation of the ConstraintStatus enum.
240 std::string GetConstraintStatusString(ConstraintStatus status);
241 
242 inline std::ostream& operator<<(std::ostream& os, ConstraintStatus status) {
243   os << GetConstraintStatusString(status);
244   return os;
245 }
246 
247 // Returns the ConstraintStatus corresponding to a given VariableStatus.
248 ConstraintStatus VariableToConstraintStatus(VariableStatus status);
249 
250 // Wrapper around an ITIVector to allow (and enforce) creation/resize/assign
251 // to use the index type for the size.
252 //
253 // TODO(user): This should probably move into ITIVector, but note that this
254 // version is more strict and does not allow any other size types.
255 template <typename IntType, typename T>
256 class StrictITIVector : public absl::StrongVector<IntType, T> {
257  public:
258   typedef IntType IndexType;
259   typedef absl::StrongVector<IntType, T> ParentType;
260 // This allows for brace initialization, which is really useful in tests.
261 // It is not 'explicit' by design, so one can do vector = {...};
262 #if !defined(__ANDROID__) && (!defined(_MSC_VER) || (_MSC_VER >= 1800))
263   StrictITIVector(std::initializer_list<T> init_list)  // NOLINT
264       : ParentType(init_list.begin(), init_list.end()) {}
265 #endif
266   StrictITIVector() : ParentType() {}
267   explicit StrictITIVector(IntType size) : ParentType(size.value()) {}
268   StrictITIVector(IntType size, const T& v) : ParentType(size.value(), v) {}
269   template <typename InputIteratorType>
270   StrictITIVector(InputIteratorType first, InputIteratorType last)
271       : ParentType(first, last) {}
272 
273   void resize(IntType size) { ParentType::resize(size.value()); }
274   void resize(IntType size, const T& v) { ParentType::resize(size.value(), v); }
275 
276   void reserve(IntType size) { ParentType::reserve(size.value()); }
277 
278   void assign(IntType size, const T& v) { ParentType::assign(size.value(), v); }
279 
280   IntType size() const { return IntType(ParentType::size()); }
281 
282   IntType capacity() const { return IntType(ParentType::capacity()); }
283 
284   // Since calls to resize() must use a default value, we introduce a new
285   // function for convenience to reduce the size of a vector.
286   void resize_down(IntType size) {
287     DCHECK_GE(size, IntType(0));
288     DCHECK_LE(size, IntType(ParentType::size()));
289     ParentType::resize(size.value());
290   }
291 
292   // This function can be up to 4 times faster than calling assign(size, 0).
293   // Note that it only works with StrictITIVector of basic types.
294   void AssignToZero(IntType size) {
295     resize(size, 0);
296     memset(ParentType::data(), 0, size.value() * sizeof(T));
297   }
298 };
299 
300 // Row-vector types. Row-vector types are indexed by a column index.
301 
302 // Row of fractional values.
303 typedef StrictITIVector<ColIndex, Fractional> DenseRow;
304 
305 // Row of booleans.
306 typedef StrictITIVector<ColIndex, bool> DenseBooleanRow;
307 
308 // Row of column indices. Used to represent mappings between columns.
309 typedef StrictITIVector<ColIndex, ColIndex> ColMapping;
310 
311 // Vector of row or column indices. Useful to list the non-zero positions.
312 typedef std::vector<ColIndex> ColIndexVector;
313 typedef std::vector<RowIndex> RowIndexVector;
314 
315 // Row of row indices.
316 // Useful for knowing which row corresponds to a particular column in the basis,
317 // or for storing the number of rows for a given column.
318 typedef StrictITIVector<ColIndex, RowIndex> ColToRowMapping;
319 
320 // Row of variable types.
321 typedef StrictITIVector<ColIndex, VariableType> VariableTypeRow;
322 
323 // Row of variable statuses.
324 typedef StrictITIVector<ColIndex, VariableStatus> VariableStatusRow;
325 
326 // Row of bits.
327 typedef Bitset64<ColIndex> DenseBitRow;
328 
329 // Column-vector types. Column-vector types are indexed by a row index.
330 
331 // Column of fractional values.
332 typedef StrictITIVector<RowIndex, Fractional> DenseColumn;
333 
334 // Column of booleans.
335 typedef StrictITIVector<RowIndex, bool> DenseBooleanColumn;
336 
337 // Column of bits.
338 typedef Bitset64<RowIndex> DenseBitColumn;
339 
340 // Column of row indices. Used to represent mappings between rows.
341 typedef StrictITIVector<RowIndex, RowIndex> RowMapping;
342 
343 // Column of column indices.
344 // Used to represent which column corresponds to a particular row in the basis,
345 // or for storing the number of columns for a given row.
346 typedef StrictITIVector<RowIndex, ColIndex> RowToColMapping;
347 
348 // Column of constraints (slack variables) statuses.
349 typedef StrictITIVector<RowIndex, ConstraintStatus> ConstraintStatusColumn;
350 
351 // --------------------------------------------------------
352 // VectorIterator
353 // --------------------------------------------------------
354 
355 // An iterator over the elements of a sparse data structure that stores the
356 // elements in arrays for indices and coefficients. The iterator is
357 // built as a wrapper over a sparse vector entry class; the concrete entry class
358 // is provided through the template argument EntryType.
359 template <typename EntryType>
360 class VectorIterator : EntryType {
361  public:
362   using Index = typename EntryType::Index;
363   using Entry = EntryType;
364 
365   VectorIterator(const Index* indices, const Fractional* coefficients,
366                  EntryIndex i)
367       : EntryType(indices, coefficients, i) {}
368 
369   void operator++() { ++this->i_; }
370   bool operator!=(const VectorIterator& other) const {
371     // This operator is intended for use in natural range iteration ONLY.
372     // Therefore, we prefer to use '<' so that a buggy range iteration which
373     // start point is *after* its end point stops immediately, instead of
374     // iterating 2^(number of bits of EntryIndex) times.
375     return this->i_ < other.i_;
376   }
377   const Entry& operator*() const { return *this; }
378 };
379 
380 // This is used during the deterministic time computation to convert a given
381 // number of floating-point operations to something in the same order of
382 // magnitude as a second (on a 2014 desktop).
383 static inline double DeterministicTimeForFpOperations(int64_t n) {
384   const double kConversionFactor = 2e-9;
385   return kConversionFactor * static_cast<double>(n);
386 }
387 
388 }  // namespace glop
389 }  // namespace operations_research
390 
391 #endif  // OR_TOOLS_LP_DATA_LP_TYPES_H_
392