Google OR-Tools v9.12
a fast and portable software suite for combinatorial optimization
Loading...
Searching...
No Matches
set_cover_invariant.cc
Go to the documentation of this file.
1// Copyright 2010-2025 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
16#include <algorithm>
17#include <limits>
18#include <tuple>
19#include <vector>
20
21#include "absl/log/check.h"
22#include "absl/types/span.h"
26
27namespace operations_research {
28
30
31// Note: in many of the member functions, variables have "crypterse" names
32// to avoid confusing them with member data. For example mrgnl_impcts is used
33// to avoid confusion with num_free_elements_.
35 DCHECK(model_->ComputeFeasibility());
36 model_->CreateSparseRowView();
37 Clear();
38}
39
41 cost_ = 0.0;
42
43 const BaseInt num_subsets = model_->num_subsets();
44 const BaseInt num_elements = model_->num_elements();
45
46 is_selected_.assign(num_subsets, false);
47 num_free_elements_.assign(num_subsets, 0);
48 num_non_overcovered_elements_.assign(num_subsets, 0);
49 is_redundant_.assign(num_subsets, false);
50
51 const SparseColumnView& columns = model_->columns();
52 for (const SubsetIndex subset : model_->SubsetRange()) {
53 num_free_elements_[subset] = columns[subset].size();
54 num_non_overcovered_elements_[subset] = columns[subset].size();
55 }
56
57 coverage_.assign(num_elements, 0);
58
59 // No need to reserve for trace_ and other vectors as extending with
60 // push_back is fast enough.
61 trace_.clear();
62 newly_removable_subsets_.clear();
63 newly_non_removable_subsets_.clear();
64
65 num_uncovered_elements_ = num_elements;
66 consistency_level_ = CL::kRedundancy;
67}
68
70 CHECK(consistency <= CL::kRedundancy);
71 if (consistency == CL::kInconsistent) {
72 return true;
73 }
74 auto [cst, cvrg] = ComputeCostAndCoverage(is_selected_);
75 CHECK(MathUtil::AlmostEquals(cost_, cst));
76 for (const ElementIndex element : model_->ElementRange()) {
77 CHECK_EQ(cvrg[element], coverage_[element]);
78 }
79 if (consistency == CL::kCostAndCoverage) {
80 return true;
81 }
82 auto [num_uncvrd_elts, num_free_elts] =
83 ComputeNumUncoveredAndFreeElements(coverage_);
84 for (const SubsetIndex subset : model_->SubsetRange()) {
85 CHECK_EQ(num_free_elts[subset], num_free_elements_[subset]);
86 }
87 if (consistency == CL::kFreeAndUncovered) {
88 return true;
89 }
90 auto [num_non_ovrcvrd_elts, is_rdndnt] = ComputeRedundancyInfo(coverage_);
91 for (const SubsetIndex subset : model_->SubsetRange()) {
92 CHECK_EQ(is_rdndnt[subset], is_redundant_[subset]);
93 CHECK_EQ(is_rdndnt[subset], num_non_ovrcvrd_elts[subset] == 0);
94 }
95 if (consistency == CL::kRedundancy) {
96 return true;
97 }
98 LOG(FATAL) << "Consistency level not supported: "
99 << static_cast<int>(consistency);
100 return false;
101}
102
104 is_selected_ = solution;
105 ClearTrace();
107 SubsetIndex subset(0);
108 for (const bool b : solution) {
109 if (b) {
110 trace_.push_back(SetCoverDecision(subset, true));
111 }
112 ++subset;
113 }
114 consistency_level_ = CL::kInconsistent;
115 Recompute(CL::kCostAndCoverage);
116}
117
118bool SetCoverInvariant::NeedToRecompute(ConsistencyLevel cheched_consistency,
119 ConsistencyLevel target_consistency) {
120 return consistency_level_ < cheched_consistency &&
121 cheched_consistency <= target_consistency;
122}
123
125 CHECK(target_consistency >= CL::kCostAndCoverage);
126 CHECK(target_consistency <= CL::kRedundancy);
127 DCHECK(CheckConsistency(consistency_level_));
128 if (NeedToRecompute(CL::kCostAndCoverage, target_consistency)) {
129 std::tie(cost_, coverage_) = ComputeCostAndCoverage(is_selected_);
130 }
131 if (NeedToRecompute(CL::kFreeAndUncovered, target_consistency)) {
132 std::tie(num_uncovered_elements_, num_free_elements_) =
133 ComputeNumUncoveredAndFreeElements(coverage_);
134 }
135 if (NeedToRecompute(CL::kRedundancy, target_consistency)) {
136 std::tie(num_non_overcovered_elements_, is_redundant_) =
137 ComputeRedundancyInfo(coverage_);
138 }
139 consistency_level_ = target_consistency;
140}
141
142// NOTE(user): This piece of code is for reference because it seems to be
143// faster to update the invariant. const BaseInt num_subsets =
144// model_->num_subsets(); is_redundant_.assign(num_subsets, false);
145// num_non_overcovered_elements_.assign(num_subsets, 0);
146// const SparseColumnView& columns = model_->columns();
147// for (const ElementIndex element : model_->ElementRange()) {
148// if (coverage_[element] >= 1) {
149// --num_uncovered_elements_;
150// }
151// }
152// for (const SubsetIndex subset : model_->SubsetRange()) {
153// for (const ElementIndex element : columns[subset]) {
154// if (coverage_[element] <= 1) {
155// ++num_non_overcovered_elements_[subset];
156// }
157// if (coverage_[element] >= 1) {
158// --num_free_elements_[subset];
159// }
160// }
161// is_redundant_[subset] = (num_non_overcovered_elements_[subset] == 0);
162// }
163
164std::tuple<Cost, ElementToIntVector> SetCoverInvariant::ComputeCostAndCoverage(
165 const SubsetBoolVector& choices) const {
166 Cost cst = 0.0;
167 ElementToIntVector cvrg(model_->num_elements(), 0);
168 const SparseColumnView& columns = model_->columns();
169 // Initialize coverage, update cost, and compute the coverage for
170 // all the elements covered by the selected subsets.
171 const SubsetCostVector& subset_costs = model_->subset_costs();
172 SubsetIndex subset(0);
173 for (const bool b : choices) {
174 if (b) {
175 cst += subset_costs[subset];
176 for (const ElementIndex element : columns[subset]) {
177 ++cvrg[element];
178 }
179 }
180 ++subset;
181 }
182 return {cst, cvrg};
183}
184
186 const absl::Span<const SubsetIndex> focus) const {
187 ElementToIntVector coverage(coverage_.size());
188 for (const SubsetIndex subset : focus) {
189 if (is_selected_[subset]) {
190 for (const ElementIndex element : model_->columns()[subset]) {
191 ++coverage[element];
192 }
193 }
194 }
195 return coverage;
196}
197
198std::tuple<BaseInt, SubsetToIntVector>
199SetCoverInvariant::ComputeNumUncoveredAndFreeElements(
200 const ElementToIntVector& cvrg) const {
201 BaseInt num_uncvrd_elts = model_->num_elements();
202
203 const BaseInt num_subsets(model_->num_subsets());
204 SubsetToIntVector num_free_elts(num_subsets, 0);
205
206 const SparseColumnView& columns = model_->columns();
207 // Initialize number of free elements and number of elements covered 0 or 1.
208 for (const SubsetIndex subset : model_->SubsetRange()) {
209 num_free_elts[subset] = columns[subset].size();
210 }
211
212 const SparseRowView& rows = model_->rows();
213 for (const ElementIndex element : model_->ElementRange()) {
214 if (cvrg[element] >= 1) {
215 --num_uncvrd_elts;
216 for (const SubsetIndex subset : rows[element]) {
217 --num_free_elts[subset];
218 }
219 }
220 }
221 return {num_uncvrd_elts, num_free_elts};
222}
223
224std::tuple<SubsetToIntVector, SubsetBoolVector>
225SetCoverInvariant::ComputeRedundancyInfo(const ElementToIntVector& cvrg) const {
226 const BaseInt num_subsets(model_->num_subsets());
227 SubsetToIntVector num_cvrg_le_1_elts(num_subsets, 0);
228 SubsetBoolVector is_rdndnt(num_subsets, false);
229
230 const SparseColumnView& columns = model_->columns();
231 // Initialize number of free elements and number of elements covered 0 or 1.
232 for (const SubsetIndex subset : model_->SubsetRange()) {
233 num_cvrg_le_1_elts[subset] = columns[subset].size();
234 }
235
236 const SparseRowView& rows = model_->rows();
237 for (const ElementIndex element : model_->ElementRange()) {
238 if (cvrg[element] >= 2) {
239 for (const SubsetIndex subset : rows[element]) {
240 --num_cvrg_le_1_elts[subset];
241 if (num_cvrg_le_1_elts[subset] == 0) {
242 is_rdndnt[subset] = true;
243 }
244 }
245 }
246 }
247 return {num_cvrg_le_1_elts, is_rdndnt};
248}
249
251 // As of 2024-05-14, this is as fast as "smarter" alternatives that try to
252 // avoid some memory writes by keeping track of already visited subsets.
253 // We also tried to use is_selected_ as a helper, which slowed down
254 // everything.
255 const SubsetIndex num_subsets(model_->num_subsets());
256 SubsetBoolVector last_value_seen(num_subsets, false);
257 for (BaseInt i = 0; i < trace_.size(); ++i) {
258 const SubsetIndex subset(trace_[i].subset());
259 last_value_seen[subset] = trace_[i].decision();
260 }
261 BaseInt w = 0; // Write index.
262 for (BaseInt i = 0; i < trace_.size(); ++i) {
263 const SubsetIndex subset(trace_[i].subset());
264 if (last_value_seen[subset]) {
265 last_value_seen[subset] = false;
266 trace_[w] = SetCoverDecision(subset, true);
267 ++w;
268 }
269 }
270 trace_.resize(w);
271}
272
273bool SetCoverInvariant::ComputeIsRedundant(SubsetIndex subset) const {
274 if (consistency_level_ >= CL::kRedundancy) {
275 return is_redundant_[subset];
276 }
277 if (is_selected_[subset]) {
278 for (const ElementIndex element : model_->columns()[subset]) {
279 if (coverage_[element] <= 1) { // If deselected, it will be <= 0...
280 return false;
281 }
282 }
283 } else {
284 for (const ElementIndex element : model_->columns()[subset]) {
285 if (coverage_[element] == 0) { // Cannot be removed from the problem.
286 return false;
287 }
288 }
289 }
290 return true;
291}
292
294 BaseInt num_free_elements = model_->columns()[subset].size();
295 for (const ElementIndex element : model_->columns()[subset]) {
296 if (coverage_[element] != 0) {
298 }
299 }
300 return num_free_elements;
301}
302
303void SetCoverInvariant::Flip(SubsetIndex subset,
304 ConsistencyLevel target_consistency) {
305 if (!is_selected_[subset]) {
306 Select(subset, target_consistency);
307 } else {
308 Deselect(subset, target_consistency);
309 }
310}
311
312void SetCoverInvariant::Select(SubsetIndex subset,
313 ConsistencyLevel target_consistency) {
314 const bool update_redundancy_info = target_consistency >= CL::kRedundancy;
315 if (update_redundancy_info) {
317 }
318 consistency_level_ = std::min(consistency_level_, target_consistency);
319 DVLOG(1) << "Selecting subset " << subset;
320 DCHECK(!is_selected_[subset]);
321 DCHECK(CheckConsistency(target_consistency));
322 trace_.push_back(SetCoverDecision(subset, true));
323 is_selected_[subset] = true;
324 const SubsetCostVector& subset_costs = model_->subset_costs();
325 cost_ += subset_costs[subset];
326 const SparseColumnView& columns = model_->columns();
327 const SparseRowView& rows = model_->rows();
328 // Fast path for kCostAndCoverage.
329 if (target_consistency == CL::kCostAndCoverage) {
330 for (const ElementIndex element : columns[subset]) {
331 ++coverage_[element];
332 }
333 return;
334 }
335 for (const ElementIndex element : columns[subset]) {
336 if (coverage_[element] == 0) {
337 // `element` will be newly covered.
338 --num_uncovered_elements_;
339 for (const SubsetIndex impacted_subset : rows[element]) {
340 --num_free_elements_[impacted_subset];
341 }
342 } else if (update_redundancy_info && coverage_[element] == 1) {
343 // `element` will be newly overcovered.
344 for (const SubsetIndex impacted_subset : rows[element]) {
345 --num_non_overcovered_elements_[impacted_subset];
346 if (num_non_overcovered_elements_[impacted_subset] == 0) {
347 // All the elements in impacted_subset are now overcovered, so it
348 // is removable. Note that this happens only when the last element
349 // of impacted_subset becomes overcovered.
350 DCHECK(!is_redundant_[impacted_subset]);
351 if (is_selected_[impacted_subset]) {
352 newly_removable_subsets_.push_back(impacted_subset);
353 }
354 is_redundant_[impacted_subset] = true;
355 }
356 }
357 }
358 // Update coverage. Notice the asymmetry with Deselect where coverage is
359 // **decremented** before being tested. This allows to have more
360 // symmetrical code for conditions.
361 ++coverage_[element];
362 }
363 if (update_redundancy_info) {
364 if (is_redundant_[subset]) {
365 newly_removable_subsets_.push_back(subset);
366 } else {
367 newly_non_removable_subsets_.push_back(subset);
368 }
369 }
370 DCHECK(CheckConsistency(target_consistency));
371}
372
373void SetCoverInvariant::Deselect(SubsetIndex subset,
374 ConsistencyLevel target_consistency) {
375 DCHECK(CheckConsistency(target_consistency));
376 const bool update_redundancy_info = target_consistency >= CL::kRedundancy;
377 if (update_redundancy_info) {
379 }
380 consistency_level_ = std::min(consistency_level_, target_consistency);
381 DVLOG(1) << "Deselecting subset " << subset;
382 // If already selected, then num_free_elements == 0.
383 DCHECK(is_selected_[subset]);
384 DCHECK_EQ(num_free_elements_[subset], 0);
385 trace_.push_back(SetCoverDecision(subset, false));
386 is_selected_[subset] = false;
387 const SubsetCostVector& subset_costs = model_->subset_costs();
388 cost_ -= subset_costs[subset];
389 const SparseColumnView& columns = model_->columns();
390 const SparseRowView& rows = model_->rows();
391 // Fast path for kCostAndCoverage.
392 if (target_consistency == CL::kCostAndCoverage) {
393 for (const ElementIndex element : columns[subset]) {
394 --coverage_[element];
395 }
396 return;
397 }
398 for (const ElementIndex element : columns[subset]) {
399 // Update coverage. Notice the asymmetry with Select where coverage is
400 // incremented after being tested.
401 --coverage_[element];
402 if (coverage_[element] == 0) {
403 // `element` is no longer covered.
404 ++num_uncovered_elements_;
405 for (const SubsetIndex impacted_subset : rows[element]) {
406 ++num_free_elements_[impacted_subset];
407 }
408 } else if (update_redundancy_info && coverage_[element] == 1) {
409 // `element` will be no longer overcovered.
410 for (const SubsetIndex impacted_subset : rows[element]) {
411 if (num_non_overcovered_elements_[impacted_subset] == 0) {
412 // There is one element of impacted_subset which is not overcovered.
413 // impacted_subset has just become non-removable.
414 DCHECK(is_redundant_[impacted_subset]);
415 if (is_selected_[impacted_subset]) {
416 newly_non_removable_subsets_.push_back(impacted_subset);
417 }
418 is_redundant_[impacted_subset] = false;
419 }
420 ++num_non_overcovered_elements_[impacted_subset];
421 }
422 }
423 }
424 // Since subset is now deselected, there is no need
425 // nor meaning in adding it a list of removable or non-removable subsets.
426 // This is a dissymmetry with Select.
427 DCHECK(CheckConsistency(target_consistency));
428}
429
430SetCoverSolutionResponse SetCoverInvariant::ExportSolutionAsProto() const {
431 SetCoverSolutionResponse message;
432 message.set_num_subsets(is_selected_.size());
433 Cost lower_bound = std::numeric_limits<Cost>::max();
434 for (const SubsetIndex subset : model_->SubsetRange()) {
435 if (is_selected_[subset]) {
436 message.add_subset(subset.value());
437 }
438 lower_bound = std::min(model_->subset_costs()[subset], lower_bound);
439 }
440 message.set_cost(cost_);
441 message.set_cost_lower_bound(lower_bound);
442 return message;
443}
444
446 const SetCoverSolutionResponse& message) {
447 is_selected_.resize(SubsetIndex(message.num_subsets()), false);
448 for (auto s : message.subset()) {
449 is_selected_[SubsetIndex(s)] = true;
450 }
451 Cost cost = message.cost();
452 CHECK(MathUtil::AlmostEquals(cost, cost_));
453}
454} // namespace operations_research
static bool AlmostEquals(const T x, const T y)
Definition mathutil.h:326
A helper class used to store the decisions made during a search.
bool ComputeIsRedundant(SubsetIndex subset) const
void Recompute(ConsistencyLevel target_consistency)
Recomputes the invariant to the given consistency level.
const SubsetToIntVector & num_free_elements() const
SetCoverSolutionResponse ExportSolutionAsProto() const
Returns the current solution as a proto.
bool CheckConsistency(ConsistencyLevel consistency) const
Checks the consistency of the invariant at the given consistency level.
void Select(SubsetIndex subset, ConsistencyLevel consistency)
void LoadSolution(const SubsetBoolVector &solution)
Loads the solution and recomputes the data in the invariant.
void ImportSolutionFromProto(const SetCoverSolutionResponse &message)
Imports the solution from a proto.
void Flip(SubsetIndex subset, ConsistencyLevel consistency)
const ElementToIntVector & coverage() const
Returns vector containing number of subsets covering each element.
BaseInt ComputeNumFreeElements(SubsetIndex subset) const
Computes the number of free (uncovered) elements in the given subset.
Cost cost() const
Returns the cost of current solution.
void Deselect(SubsetIndex subset, ConsistencyLevel consistency)
ElementToIntVector ComputeCoverageInFocus(absl::Span< const SubsetIndex > focus) const
void Clear()
Clears the invariant. Also called by Initialize.
void ClearRemovabilityInformation()
Clear the removability information.
util_intops::StrongIntRange< ElementIndex > ElementRange() const
const SparseRowView & rows() const
Row view of the set covering problem.
util_intops::StrongIntRange< SubsetIndex > SubsetRange() const
Access to the ranges of subsets and elements.
const SubsetCostVector & subset_costs() const
Vector of costs for each subset.
const SparseColumnView & columns() const
Column view of the set covering problem.
In SWIG mode, we don't want anything besides these top-level includes.
util_intops::StrongVector< SubsetIndex, Cost > SubsetCostVector
Select next search node to expand Select next item_i to add this new search node to the search Generate a new search node where item_i is not in the knapsack Check validity of this new partial solution(using propagators) - If valid
SetCoverInvariant::ConsistencyLevel CL
util_intops::StrongVector< SubsetIndex, BaseInt > SubsetToIntVector
util_intops::StrongVector< ElementIndex, BaseInt > ElementToIntVector
util_intops::StrongVector< SubsetIndex, bool > SubsetBoolVector
util_intops::StrongVector< SubsetIndex, SparseColumn > SparseColumnView
double Cost
Basic non-strict type for cost. The speed penalty for using double is ~2%.
util_intops::StrongVector< ElementIndex, SparseRow > SparseRowView
trees with all degrees equal w the current value of w