OR-Tools  8.2
zero_half_cuts.cc
Go to the documentation of this file.
1// Copyright 2010-2018 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
15
16namespace operations_research {
17namespace sat {
18
20 rows_.clear();
21 shifted_lp_values_.clear();
22 bound_parity_.clear();
23 col_to_rows_.clear();
24 col_to_rows_.resize(size);
25 tmp_marked_.resize(size);
26}
27
29 const std::vector<double>& lp_values,
30 const std::vector<IntegerValue>& lower_bounds,
31 const std::vector<IntegerValue>& upper_bounds) {
32 Reset(lp_values.size());
33
34 // Shift all variables to their closest bound.
35 lp_values_ = lp_values;
36 for (int i = 0; i < lp_values.size(); ++i) {
37 const double lb_dist = lp_values[i] - ToDouble(lower_bounds[i]);
38 const double ub_dist = ToDouble(upper_bounds[i]) - lp_values_[i];
39 if (lb_dist < ub_dist) {
40 shifted_lp_values_.push_back(lb_dist);
41 bound_parity_.push_back(lower_bounds[i].value() & 1);
42 } else {
43 shifted_lp_values_.push_back(ub_dist);
44 bound_parity_.push_back(upper_bounds[i].value() & 1);
45 }
46 }
47}
48
50 // No point pushing an all zero row with a zero rhs.
51 if (binary_row.cols.empty() && !binary_row.rhs_parity) return;
52 for (const int col : binary_row.cols) {
53 col_to_rows_[col].push_back(rows_.size());
54 }
55 rows_.push_back(binary_row);
56}
57
59 const glop::RowIndex row,
60 const std::vector<std::pair<glop::ColIndex, IntegerValue>>& terms,
61 IntegerValue lb, IntegerValue ub) {
62 if (terms.size() > kMaxInputConstraintSize) return;
63
64 double activity = 0.0;
65 IntegerValue magnitude(0);
66 CombinationOfRows binary_row;
67 int rhs_adjust = 0;
68 for (const auto& term : terms) {
69 const int col = term.first.value();
70 activity += ToDouble(term.second) * lp_values_[col];
71 magnitude = std::max(magnitude, IntTypeAbs(term.second));
72
73 // Only consider odd coefficient.
74 if ((term.second.value() & 1) == 0) continue;
75
76 // Ignore column in the binary matrix if its lp value is almost zero.
77 if (shifted_lp_values_[col] > 1e-2) {
78 binary_row.cols.push_back(col);
79 }
80
81 // Because we work on the shifted variable, the rhs needs to be updated.
82 rhs_adjust ^= bound_parity_[col];
83 }
84
85 // We ignore constraint with large coefficient, since there is little chance
86 // to cancel them and because of that the efficacity of a generated cut will
87 // be limited.
88 if (magnitude > kMaxInputConstraintMagnitude) return;
89
90 // TODO(user): experiment with the best value. probably only tight rows are
91 // best? and we could use the basis status rather than recomputing the
92 // activity for that.
93 //
94 // TODO(user): Avoid adding duplicates and just randomly pick one. Note
95 // that we should also remove duplicate in a generic way.
96 const double tighteness_threshold = 1e-2;
97 if (ToDouble(ub) - activity < tighteness_threshold) {
98 binary_row.multipliers = {{row, IntegerValue(1)}};
99 binary_row.slack = ToDouble(ub) - activity;
100 binary_row.rhs_parity = (ub.value() & 1) ^ rhs_adjust;
101 AddBinaryRow(binary_row);
102 }
103 if (activity - ToDouble(lb) < tighteness_threshold) {
104 binary_row.multipliers = {{row, IntegerValue(-1)}};
105 binary_row.slack = activity - ToDouble(lb);
106 binary_row.rhs_parity = (lb.value() & 1) ^ rhs_adjust;
107 AddBinaryRow(binary_row);
108 }
109}
110
112 std::function<bool(int)> extra_condition, const std::vector<int>& a,
113 std::vector<int>* b) {
114 for (const int v : *b) tmp_marked_[v] = true;
115 for (const int v : a) {
116 if (tmp_marked_[v]) {
117 tmp_marked_[v] = false;
118 } else {
119 tmp_marked_[v] = true;
120
121 // TODO(user): optim by doing that at the end?
122 b->push_back(v);
123 }
124 }
125
126 // Remove position that are not marked, and clear tmp_marked_.
127 int new_size = 0;
128 for (const int v : *b) {
129 if (tmp_marked_[v]) {
130 if (extra_condition(v)) {
131 (*b)[new_size++] = v;
132 }
133 tmp_marked_[v] = false;
134 }
135 }
136 b->resize(new_size);
137}
138
139void ZeroHalfCutHelper::ProcessSingletonColumns() {
140 for (const int singleton_col : singleton_cols_) {
141 if (col_to_rows_[singleton_col].empty()) continue;
142 CHECK_EQ(col_to_rows_[singleton_col].size(), 1);
143 const int row = col_to_rows_[singleton_col][0];
144 int new_size = 0;
145 auto& mutable_cols = rows_[row].cols;
146 for (const int col : mutable_cols) {
147 if (col == singleton_col) continue;
148 mutable_cols[new_size++] = col;
149 }
150 CHECK_LT(new_size, mutable_cols.size());
151 mutable_cols.resize(new_size);
152 col_to_rows_[singleton_col].clear();
153 rows_[row].slack += shifted_lp_values_[singleton_col];
154 }
155 singleton_cols_.clear();
156}
157
158// This is basically one step of a Gaussian elimination with the given pivot.
160 int eliminated_row) {
161 CHECK_LE(rows_[eliminated_row].slack, 1e-6);
162 CHECK(!rows_[eliminated_row].cols.empty());
163
164 // First update the row representation of the matrix.
165 tmp_marked_.resize(std::max(col_to_rows_.size(), rows_.size()));
166 DCHECK(std::all_of(tmp_marked_.begin(), tmp_marked_.end(),
167 [](bool b) { return !b; }));
168 int new_size = 0;
169 for (const int other_row : col_to_rows_[eliminated_col]) {
170 if (other_row == eliminated_row) continue;
171 col_to_rows_[eliminated_col][new_size++] = other_row;
172
173 SymmetricDifference([](int i) { return true; }, rows_[eliminated_row].cols,
174 &rows_[other_row].cols);
175
176 // Update slack & parity.
177 rows_[other_row].rhs_parity ^= rows_[eliminated_row].rhs_parity;
178 rows_[other_row].slack += rows_[eliminated_row].slack;
179
180 // Update the multipliers the same way.
181 {
182 auto& mutable_multipliers = rows_[other_row].multipliers;
183 mutable_multipliers.insert(mutable_multipliers.end(),
184 rows_[eliminated_row].multipliers.begin(),
185 rows_[eliminated_row].multipliers.end());
186 std::sort(mutable_multipliers.begin(), mutable_multipliers.end());
187 int new_size = 0;
188 for (const auto& entry : mutable_multipliers) {
189 if (new_size > 0 && entry == mutable_multipliers[new_size - 1]) {
190 // Cancel both.
191 --new_size;
192 } else {
193 mutable_multipliers[new_size++] = entry;
194 }
195 }
196 mutable_multipliers.resize(new_size);
197 }
198 }
199 col_to_rows_[eliminated_col].resize(new_size);
200
201 // Then update the col representation of the matrix.
202 //
203 // Note that we remove from the col-wise representation any rows with a large
204 // slack.
205 {
206 int new_size = 0;
207 for (const int other_col : rows_[eliminated_row].cols) {
208 if (other_col == eliminated_col) continue;
209 const int old_size = col_to_rows_[other_col].size();
210 rows_[eliminated_row].cols[new_size++] = other_col;
212 [this](int i) { return rows_[i].slack < kSlackThreshold; },
213 col_to_rows_[eliminated_col], &col_to_rows_[other_col]);
214 if (old_size != 1 && col_to_rows_[other_col].size() == 1) {
215 singleton_cols_.push_back(other_col);
216 }
217 }
218 rows_[eliminated_row].cols.resize(new_size);
219 }
220
221 // Clear col.
222 col_to_rows_[eliminated_col].clear();
223 rows_[eliminated_row].slack += shifted_lp_values_[eliminated_col];
224}
225
226std::vector<std::vector<std::pair<glop::RowIndex, IntegerValue>>>
228 std::vector<std::vector<std::pair<glop::RowIndex, IntegerValue>>> result;
229
230 // Initialize singleton_cols_.
231 singleton_cols_.clear();
232 for (int col = 0; col < col_to_rows_.size(); ++col) {
233 if (col_to_rows_[col].size() == 1) singleton_cols_.push_back(col);
234 }
235
236 // Process rows by increasing size, but randomize if same size.
237 std::vector<int> to_process;
238 for (int row = 0; row < rows_.size(); ++row) to_process.push_back(row);
239 std::shuffle(to_process.begin(), to_process.end(), *random);
240 std::stable_sort(to_process.begin(), to_process.end(), [this](int a, int b) {
241 return rows_[a].cols.size() < rows_[b].cols.size();
242 });
243
244 for (const int row : to_process) {
245 ProcessSingletonColumns();
246
247 if (rows_[row].cols.empty()) continue;
248 if (rows_[row].slack > 1e-6) continue;
249 if (rows_[row].multipliers.size() > kMaxAggregationSize) continue;
250
251 // Heuristic: eliminate the variable with highest shifted lp value.
252 int eliminated_col = -1;
253 double max_lp_value = 0.0;
254 for (const int col : rows_[row].cols) {
255 if (shifted_lp_values_[col] > max_lp_value) {
256 max_lp_value = shifted_lp_values_[col];
257 eliminated_col = col;
258 }
259 }
260 if (eliminated_col == -1) continue;
261
262 EliminateVarUsingRow(eliminated_col, row);
263 }
264
265 // As an heuristic, we just try to add zero rows with an odd rhs and a low
266 // enough slack.
267 for (const auto& row : rows_) {
268 if (row.cols.empty() && row.rhs_parity && row.slack < kSlackThreshold) {
269 result.push_back(row.multipliers);
270 }
271 }
272 VLOG(1) << "#candidates: " << result.size() << " / " << rows_.size();
273 return result;
274}
275
276} // namespace sat
277} // namespace operations_research
int64 max
Definition: alldiff_cst.cc:139
#define CHECK(condition)
Definition: base/logging.h:495
#define CHECK_LT(val1, val2)
Definition: base/logging.h:700
#define CHECK_EQ(val1, val2)
Definition: base/logging.h:697
#define DCHECK(condition)
Definition: base/logging.h:884
#define CHECK_LE(val1, val2)
Definition: base/logging.h:699
#define VLOG(verboselevel)
Definition: base/logging.h:978
void AddBinaryRow(const CombinationOfRows &binary_row)
void SymmetricDifference(std::function< bool(int)> extra_condition, const std::vector< int > &a, std::vector< int > *b)
void ProcessVariables(const std::vector< double > &lp_values, const std::vector< IntegerValue > &lower_bounds, const std::vector< IntegerValue > &upper_bounds)
void AddOneConstraint(glop::RowIndex, const std::vector< std::pair< glop::ColIndex, IntegerValue > > &terms, IntegerValue lb, IntegerValue ub)
std::vector< std::vector< std::pair< glop::RowIndex, IntegerValue > > > InterestingCandidates(ModelRandomGenerator *random)
int64 value
ColIndex col
Definition: markowitz.cc:176
RowIndex row
Definition: markowitz.cc:175
IntType IntTypeAbs(IntType t)
Definition: integer.h:77
double ToDouble(IntegerValue value)
Definition: integer.h:69
The vehicle routing library lets one model and solve generic vehicle routing problems ranging from th...
std::vector< double > lower_bounds
std::vector< double > upper_bounds
std::vector< std::pair< glop::RowIndex, IntegerValue > > multipliers