From e8540c1363ffbc69dc4b2242a9a3bde1ebdb9f2f Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Mon, 16 Oct 2017 11:20:54 +0200 Subject: [PATCH] add vrp constraint to the sat solver; remove old no_cycle code --- makefiles/Makefile.gen.mk | 36 +-- ortools/linear_solver/linear_expr.cc | 9 + ortools/linear_solver/linear_expr.h | 6 +- ortools/linear_solver/linear_solver.cc | 12 +- ortools/linear_solver/linear_solver.h | 4 - ortools/sat/circuit.cc | 283 +++++++++++++++++++++ ortools/sat/circuit.h | 149 +++++++++++ ortools/sat/clause.cc | 8 +- ortools/sat/cp_constraints.cc | 249 ------------------- ortools/sat/cp_constraints.h | 83 ------- ortools/sat/cp_model.proto | 22 ++ ortools/sat/cp_model_checker.cc | 61 +++++ ortools/sat/cp_model_solver.cc | 24 ++ ortools/sat/cp_model_utils.cc | 12 + ortools/sat/no_cycle.cc | 329 ------------------------- ortools/sat/no_cycle.h | 183 -------------- ortools/sat/python/cp_model.py | 1 + ortools/sat/sat_solver.h | 2 +- 18 files changed, 593 insertions(+), 880 deletions(-) create mode 100644 ortools/sat/circuit.cc create mode 100644 ortools/sat/circuit.h delete mode 100644 ortools/sat/no_cycle.cc delete mode 100644 ortools/sat/no_cycle.h diff --git a/makefiles/Makefile.gen.mk b/makefiles/Makefile.gen.mk index 67d750c230..bedbceb927 100644 --- a/makefiles/Makefile.gen.mk +++ b/makefiles/Makefile.gen.mk @@ -1526,6 +1526,7 @@ SAT_DEPS = \ SAT_LIB_OBJS = \ $(OBJ_DIR)/sat/all_different.$O \ $(OBJ_DIR)/sat/boolean_problem.$O \ + $(OBJ_DIR)/sat/circuit.$O \ $(OBJ_DIR)/sat/clause.$O \ $(OBJ_DIR)/sat/cp_constraints.$O \ $(OBJ_DIR)/sat/cp_model_checker.$O \ @@ -1543,7 +1544,6 @@ SAT_LIB_OBJS = \ $(OBJ_DIR)/sat/intervals.$O \ $(OBJ_DIR)/sat/linear_programming_constraint.$O \ $(OBJ_DIR)/sat/lp_utils.$O \ - $(OBJ_DIR)/sat/no_cycle.$O \ $(OBJ_DIR)/sat/optimization.$O \ $(OBJ_DIR)/sat/overload_checker.$O \ $(OBJ_DIR)/sat/pb_constraint.$O \ @@ -1579,6 +1579,16 @@ $(SRC_DIR)/ortools/sat/boolean_problem.h: \ $(SRC_DIR)/ortools/base/status.h \ $(SRC_DIR)/ortools/algorithms/sparse_permutation.h +$(SRC_DIR)/ortools/sat/circuit.h: \ + $(SRC_DIR)/ortools/sat/integer.h \ + $(SRC_DIR)/ortools/sat/model.h \ + $(SRC_DIR)/ortools/sat/sat_base.h \ + $(SRC_DIR)/ortools/base/integral_types.h \ + $(SRC_DIR)/ortools/base/int_type.h \ + $(SRC_DIR)/ortools/base/logging.h \ + $(SRC_DIR)/ortools/base/macros.h \ + $(SRC_DIR)/ortools/util/rev.h + $(SRC_DIR)/ortools/sat/clause.h: \ $(SRC_DIR)/ortools/sat/sat_base.h \ $(GEN_DIR)/ortools/sat/sat_parameters.pb.h \ @@ -1733,14 +1743,6 @@ $(SRC_DIR)/ortools/sat/model.h: \ $(SRC_DIR)/ortools/base/map_util.h \ $(SRC_DIR)/ortools/base/typeid.h -$(SRC_DIR)/ortools/sat/no_cycle.h: \ - $(SRC_DIR)/ortools/sat/sat_base.h \ - $(SRC_DIR)/ortools/base/integral_types.h \ - $(SRC_DIR)/ortools/base/int_type.h \ - $(SRC_DIR)/ortools/base/int_type_indexed_vector.h \ - $(SRC_DIR)/ortools/base/macros.h \ - $(SRC_DIR)/ortools/base/span.h - $(SRC_DIR)/ortools/sat/optimization.h: \ $(GEN_DIR)/ortools/sat/boolean_problem.pb.h \ $(SRC_DIR)/ortools/sat/integer.h \ @@ -1888,11 +1890,17 @@ $(OBJ_DIR)/sat/boolean_problem.$O: \ $(SRC_DIR)/ortools/base/map_util.h \ $(SRC_DIR)/ortools/base/stringprintf.h \ $(SRC_DIR)/ortools/algorithms/find_graph_symmetries.h \ - $(SRC_DIR)/ortools/graph/graph.h \ $(SRC_DIR)/ortools/graph/io.h \ $(SRC_DIR)/ortools/graph/util.h $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sortools$Ssat$Sboolean_problem.cc $(OBJ_OUT)$(OBJ_DIR)$Ssat$Sboolean_problem.$O +$(OBJ_DIR)/sat/circuit.$O: \ + $(SRC_DIR)/ortools/sat/circuit.cc \ + $(SRC_DIR)/ortools/sat/circuit.h \ + $(SRC_DIR)/ortools/sat/sat_solver.h \ + $(SRC_DIR)/ortools/base/map_util.h + $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sortools$Ssat$Scircuit.cc $(OBJ_OUT)$(OBJ_DIR)$Ssat$Scircuit.$O + $(OBJ_DIR)/sat/clause.$O: \ $(SRC_DIR)/ortools/sat/clause.cc \ $(SRC_DIR)/ortools/sat/clause.h \ @@ -1947,6 +1955,7 @@ $(OBJ_DIR)/sat/cp_model_search.$O: \ $(OBJ_DIR)/sat/cp_model_solver.$O: \ $(SRC_DIR)/ortools/sat/cp_model_solver.cc \ $(SRC_DIR)/ortools/sat/all_different.h \ + $(SRC_DIR)/ortools/sat/circuit.h \ $(SRC_DIR)/ortools/sat/cp_constraints.h \ $(SRC_DIR)/ortools/sat/cp_model_checker.h \ $(SRC_DIR)/ortools/sat/cp_model_presolve.h \ @@ -2083,12 +2092,6 @@ $(OBJ_DIR)/sat/lp_utils.$O: \ $(GEN_DIR)/ortools/glop/parameters.pb.h $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sortools$Ssat$Slp_utils.cc $(OBJ_OUT)$(OBJ_DIR)$Ssat$Slp_utils.$O -$(OBJ_DIR)/sat/no_cycle.$O: \ - $(SRC_DIR)/ortools/sat/no_cycle.cc \ - $(SRC_DIR)/ortools/sat/no_cycle.h \ - $(SRC_DIR)/ortools/base/logging.h - $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sortools$Ssat$Sno_cycle.cc $(OBJ_OUT)$(OBJ_DIR)$Ssat$Sno_cycle.$O - $(OBJ_DIR)/sat/optimization.$O: \ $(SRC_DIR)/ortools/sat/optimization.cc \ $(SRC_DIR)/ortools/sat/boolean_problem.h \ @@ -2725,6 +2728,7 @@ $(OBJ_DIR)/linear_solver/gurobi_interface.$O: \ $(OBJ_DIR)/linear_solver/linear_expr.$O: \ $(SRC_DIR)/ortools/linear_solver/linear_expr.cc \ $(SRC_DIR)/ortools/linear_solver/linear_expr.h \ + $(SRC_DIR)/ortools/linear_solver/linear_solver.h \ $(SRC_DIR)/ortools/base/logging.h $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sortools$Slinear_solver$Slinear_expr.cc $(OBJ_OUT)$(OBJ_DIR)$Slinear_solver$Slinear_expr.$O diff --git a/ortools/linear_solver/linear_expr.cc b/ortools/linear_solver/linear_expr.cc index b4258af9b7..1f8271124e 100644 --- a/ortools/linear_solver/linear_expr.cc +++ b/ortools/linear_solver/linear_expr.cc @@ -16,6 +16,7 @@ #include #include "ortools/base/logging.h" +#include "ortools/linear_solver/linear_solver.h" namespace operations_research { @@ -70,6 +71,14 @@ LinearExpr LinearExpr::NotVar(LinearExpr var) { return var; } +double LinearExpr::SolutionValue() const { + double solution = offset_; + for (const auto& pair : terms_) { + solution += pair.first->solution_value() * pair.second; + } + return solution; +} + LinearExpr operator+(LinearExpr lhs, const LinearExpr& rhs) { lhs += rhs; return lhs; diff --git a/ortools/linear_solver/linear_expr.h b/ortools/linear_solver/linear_expr.h index 610296a67a..f74d806df6 100644 --- a/ortools/linear_solver/linear_expr.h +++ b/ortools/linear_solver/linear_expr.h @@ -92,7 +92,7 @@ class MPVariable; // * Get the value of the quantity after solving, e.g. // // solver.Solve(); -// solver.SolutionValue(linear_expr); +// linear_expr.SolutionValue(); // // LinearExpr is allowed to delete variables with coefficient zero from the map, // but is not obligated to do so. @@ -122,6 +122,10 @@ class LinearExpr { return terms_; } + // Call only after calling MPSolver::Solve. Evaluates the value of this + // expression at the solution found. + double SolutionValue() const; + private: double offset_; std::unordered_map terms_; diff --git a/ortools/linear_solver/linear_solver.cc b/ortools/linear_solver/linear_solver.cc index 12a950ed94..5ffe3f938e 100644 --- a/ortools/linear_solver/linear_solver.cc +++ b/ortools/linear_solver/linear_solver.cc @@ -989,14 +989,6 @@ MPSolver::ResultStatus MPSolver::Solve(const MPSolverParameters& param) { return status; } -double MPSolver::SolutionValue(const LinearExpr& linear_expr) const { - double ans = linear_expr.offset(); - for (const auto& kv : linear_expr.terms()) { - ans += (kv.second * kv.first->solution_value()); - } - return ans; -} - void MPSolver::Write(const std::string& file_name) { interface_->Write(file_name); } namespace { @@ -1323,7 +1315,7 @@ bool MPSolverInterface::CheckSolutionIsSynchronized() const { } // Default version that can be overwritten by a solver-specific -// version to accomodate for the quirks of each solver. +// version to accommodate for the quirks of each solver. bool MPSolverInterface::CheckSolutionExists() const { if (result_status_ != MPSolver::OPTIMAL && result_status_ != MPSolver::FEASIBLE) { @@ -1335,7 +1327,7 @@ bool MPSolverInterface::CheckSolutionExists() const { } // Default version that can be overwritten by a solver-specific -// version to accomodate for the quirks of each solver. +// version to accommodate for the quirks of each solver. bool MPSolverInterface::CheckBestObjectiveBoundExists() const { if (result_status_ != MPSolver::OPTIMAL && result_status_ != MPSolver::FEASIBLE) { diff --git a/ortools/linear_solver/linear_solver.h b/ortools/linear_solver/linear_solver.h index e1104726e7..29155daeff 100644 --- a/ortools/linear_solver/linear_solver.h +++ b/ortools/linear_solver/linear_solver.h @@ -332,10 +332,6 @@ class MPSolver { // Solves the problem using the specified parameter values. ResultStatus Solve(const MPSolverParameters& param); - // Call only after calling MPSolver::Solve. Evaluates "linear_expr" for the - // variable values at the solution found by solving. - double SolutionValue(const LinearExpr& linear_expr) const; - // Writes the model using the solver internal write function. Currently only // available for Gurobi. void Write(const std::string& file_name); diff --git a/ortools/sat/circuit.cc b/ortools/sat/circuit.cc new file mode 100644 index 0000000000..e1be718473 --- /dev/null +++ b/ortools/sat/circuit.cc @@ -0,0 +1,283 @@ +// Copyright 2010-2014 Google +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#include "ortools/sat/circuit.h" + +#include +#include + +#include "ortools/base/map_util.h" +#include "ortools/sat/sat_solver.h" + +namespace operations_research { +namespace sat { + +CircuitPropagator::CircuitPropagator( + const std::vector>& graph, Options options, + Trail* trail) + : num_nodes_(graph.size()), + options_(options), + trail_(trail), + propagation_trail_index_(0) { + // TODO(user): add a way to properly handle trivially UNSAT cases. + // For now we just check that they don't occur at construction. + CHECK_GT(num_nodes_, 1) + << "Trivial or UNSAT constraint, shouldn't be constructed!"; + next_.resize(num_nodes_, -1); + prev_.resize(num_nodes_, -1); + next_literal_.resize(num_nodes_); + self_arcs_.resize(num_nodes_); + std::unordered_map literal_to_watch_index; + const VariablesAssignment& assignment = trail->Assignment(); + for (int tail = 0; tail < num_nodes_; ++tail) { + self_arcs_[tail] = graph[tail][tail]; + for (int head = 0; head < num_nodes_; ++head) { + const LiteralIndex index = graph[tail][head]; + // Note that we need to test for both "special" cases before we can + // call assignment.LiteralIsTrue() or LiteralIsFalse(). + if (index == kFalseLiteralIndex) continue; + if (index == kTrueLiteralIndex || + assignment.LiteralIsTrue(Literal(index))) { + CHECK_EQ(next_[tail], -1) + << "Trivially UNSAT or duplicate arcs while adding " << tail + << " -> " << head; + CHECK_EQ(prev_[head], -1) + << "Trivially UNSAT or duplicate arcs while adding " << tail + << " -> " << head; + AddArc(tail, head, kNoLiteralIndex); + continue; + } + if (assignment.LiteralIsFalse(Literal(index))) continue; + + int watch_index_ = FindWithDefault(literal_to_watch_index, index, -1); + if (watch_index_ == -1) { + watch_index_ = watch_index_to_literal_.size(); + literal_to_watch_index[index] = watch_index_; + watch_index_to_literal_.push_back(Literal(index)); + watch_index_to_arcs_.push_back(std::vector()); + } + watch_index_to_arcs_[watch_index_].push_back({tail, head}); + } + } +} + +void CircuitPropagator::SetLevel(int level) { + if (level == level_ends_.size()) return; + if (level > level_ends_.size()) { + while (level > level_ends_.size()) { + level_ends_.push_back(added_arcs_.size()); + } + return; + } + + // Backtrack. + for (int i = level_ends_[level]; i < added_arcs_.size(); ++i) { + const Arc arc = added_arcs_[i]; + next_[arc.tail] = -1; + prev_[arc.head] = -1; + } + added_arcs_.resize(level_ends_[level]); + level_ends_.resize(level); +} + +void CircuitPropagator::FillConflictFromCircuitAt(int start) { + std::vector* conflict = trail_->MutableConflict(); + conflict->clear(); + int node = start; + do { + CHECK_NE(node, -1); + if (next_literal_[node] != kNoLiteralIndex) { + conflict->push_back(Literal(next_literal_[node]).Negated()); + } + node = next_[node]; + } while (node != start); +} + +bool CircuitPropagator::Propagate() { return true; } + +// If multiple_subcircuit_through_zero is true, we never fill next_[0] and +// prev_[0]. +void CircuitPropagator::AddArc(int tail, int head, LiteralIndex literal_index) { + if (tail != 0 || !options_.multiple_subcircuit_through_zero) { + next_[tail] = head; + next_literal_[tail] = literal_index; + } + if (head != 0 || !options_.multiple_subcircuit_through_zero) { + prev_[head] = tail; + } +} + +bool CircuitPropagator::IncrementalPropagate( + const std::vector& watch_indices) { + for (const int w : watch_indices) { + const Literal literal = watch_index_to_literal_[w]; + for (const Arc arc : watch_index_to_arcs_[w]) { + // Get rid of the trivial conflicts: + // - At most one incoming and one ougtoing arc for each nodes. + if (next_[arc.tail] != -1) { + std::vector* conflict = trail_->MutableConflict(); + if (next_literal_[arc.tail] != kNoLiteralIndex) { + *conflict = {Literal(next_literal_[arc.tail]).Negated(), + literal.Negated()}; + } else { + *conflict = {literal.Negated()}; + } + return false; + } + if (prev_[arc.head] != -1) { + std::vector* conflict = trail_->MutableConflict(); + if (next_literal_[prev_[arc.head]] != kNoLiteralIndex) { + *conflict = {Literal(next_literal_[prev_[arc.head]]).Negated(), + literal.Negated()}; + } else { + *conflict = {literal.Negated()}; + } + return false; + } + + // Add the arc. + AddArc(arc.tail, arc.head, literal.Index()); + added_arcs_.push_back(arc); + + // Circuit? + in_circuit_.assign(num_nodes_, false); + in_circuit_[arc.tail] = true; + int size = 1; + int node = arc.head; + while (node != arc.tail && node != -1) { + in_circuit_[node] = true; + node = next_[node]; + size++; + } + + if (options_.multiple_subcircuit_through_zero) { + // If we reached zero, this is a valid path provided that we can + // reach the beginning of the path from zero. Note that we only check + // the basic case that the beginning of the path must have "open" arcs + // thanks to ExactlyOnePerRowAndPerColumn(). + if (node == 0 || node != arc.tail) continue; + + // We have a cycle not touching zero, this is a conflict. + FillConflictFromCircuitAt(arc.tail); + return false; + } + + if (node != arc.tail) continue; + + // We have one circuit. + if (size == num_nodes_) return true; + if (size == 1) continue; // self-arc. + + // HACK: we can reuse the conflict vector even though we don't have a + // conflict. + FillConflictFromCircuitAt(arc.tail); + BooleanVariable variable_with_same_reason = kNoBooleanVariable; + + // We can propagate all the other nodes to point to themselves. + // If this is not already the case, we have a conflict. + for (int node = 0; node < num_nodes_; ++node) { + if (in_circuit_[node] || next_[node] == node) continue; + if (next_[node] != -1) { + std::vector* conflict = trail_->MutableConflict(); + if (next_literal_[node] != kNoLiteralIndex) { + conflict->push_back(Literal(next_literal_[node]).Negated()); + } + return false; + } else if (self_arcs_[node] == kFalseLiteralIndex) { + return false; + } else { + DCHECK_NE(self_arcs_[node], kTrueLiteralIndex); + const Literal literal(self_arcs_[node]); + + // We may not have processed this literal yet. + if (trail_->Assignment().LiteralIsTrue(literal)) continue; + if (trail_->Assignment().LiteralIsFalse(literal)) { + std::vector* conflict = trail_->MutableConflict(); + conflict->push_back(literal); + return false; + } + + // Propagate. + if (variable_with_same_reason == kNoBooleanVariable) { + variable_with_same_reason = literal.Variable(); + const int index = trail_->Index(); + trail_->Enqueue(literal, AssignmentType::kCachedReason); + *trail_->GetVectorToStoreReason(index) = *trail_->MutableConflict(); + trail_->NotifyThatReasonIsCached(literal.Variable()); + } else { + trail_->EnqueueWithSameReasonAs(literal, variable_with_same_reason); + } + } + } + } + } + return true; +} + +void CircuitPropagator::RegisterWith(GenericLiteralWatcher* watcher) { + const int id = watcher->Register(this); + for (int w = 0; w < watch_index_to_literal_.size(); ++w) { + watcher->WatchLiteral(watch_index_to_literal_[w], id, w); + } + watcher->RegisterReversibleClass(id, this); + watcher->RegisterReversibleInt(id, &propagation_trail_index_); +} + +std::function ExactlyOnePerRowAndPerColumn( + const std::vector>& square_matrix, + bool ignore_row_and_col_zero) { + return [=](Model* model) { + int n = square_matrix.size(); + int num_trivially_false = 0; + Trail* trail = model->GetOrCreate(); + std::vector exactly_one_constraint; + for (const bool transpose : {false, true}) { + for (int i = ignore_row_and_col_zero ? 1 : 0; i < n; ++i) { + int num_true = 0; + exactly_one_constraint.clear(); + for (int j = 0; j < n; ++j) { + CHECK_EQ(n, square_matrix[i].size()); + const LiteralIndex index = + transpose ? square_matrix[j][i] : square_matrix[i][j]; + if (index == kFalseLiteralIndex) continue; + if (index == kTrueLiteralIndex) { + ++num_true; + continue; + } + exactly_one_constraint.push_back(Literal(index)); + } + if (num_true > 1) { + LOG(WARNING) << "UNSAT in ExactlyOnePerRowAndPerColumn()."; + return model->GetOrCreate()->NotifyThatModelIsUnsat(); + } + CHECK_LE(num_true, 1); + if (num_true == 1) { + for (const Literal l : exactly_one_constraint) { + if (!trail->Assignment().VariableIsAssigned(l.Variable())) { + ++num_trivially_false; + trail->EnqueueWithUnitReason(l.Negated()); + } + } + } else { + model->Add(ExactlyOneConstraint(exactly_one_constraint)); + } + } + } + if (num_trivially_false > 0) { + LOG(INFO) << "Num extra fixed literal: " << num_trivially_false; + } + }; +} + +} // namespace sat +} // namespace operations_research diff --git a/ortools/sat/circuit.h b/ortools/sat/circuit.h new file mode 100644 index 0000000000..2c35320ea8 --- /dev/null +++ b/ortools/sat/circuit.h @@ -0,0 +1,149 @@ +// Copyright 2010-2014 Google +// Licensed under the Apache License, Version 2.0 (the "License"); +// you may not use this file except in compliance with the License. +// You may obtain a copy of the License at +// +// http://www.apache.org/licenses/LICENSE-2.0 +// +// Unless required by applicable law or agreed to in writing, software +// distributed under the License is distributed on an "AS IS" BASIS, +// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +// See the License for the specific language governing permissions and +// limitations under the License. + +#ifndef OR_TOOLS_SAT_CIRCUIT_H_ +#define OR_TOOLS_SAT_CIRCUIT_H_ + +#include +#include +#include + +#include "ortools/base/integral_types.h" +#include "ortools/base/logging.h" +#include "ortools/base/macros.h" +#include "ortools/base/int_type.h" +#include "ortools/sat/integer.h" +#include "ortools/sat/model.h" +#include "ortools/sat/sat_base.h" +#include "ortools/util/rev.h" + +namespace operations_research { +namespace sat { + +// Initial version of the circuit/sub-circuit constraint. +// +// Nodes that are not in the unique allowed sub-circuit must point to themseves. +// A nodes that has no self-arc must thus be inside the sub-circuit. If there is +// no self-arc at all, then this constaint forces the circuit to go through all +// the nodes. +// +// Important: for correctness, this constraint requires to call +// ExactlyOnePerRowAndPerColumn() on the same graph. +class CircuitPropagator : PropagatorInterface, ReversibleInterface { + public: + struct Options { + // Hack for the VRP to allow for more than one sub-circuit and forces all + // the subcircuits to go through the node zero. + bool multiple_subcircuit_through_zero = false; + }; + + // The constraints take a dense representation of a graph on [0, n). Each arc + // being present when the given literal is true. The special values + // kTrueLiteralIndex and kFalseLiteralIndex can be used for arcs that are + // either always there or never there. + CircuitPropagator(const std::vector>& graph, + Options options, Trail* trail); + + void SetLevel(int level) final; + bool Propagate() final; + bool IncrementalPropagate(const std::vector& watch_indices) final; + void RegisterWith(GenericLiteralWatcher* watcher); + + private: + // Updates the structures when the given arc is added to the paths. + void AddArc(int tail, int head, LiteralIndex literal_index); + + // Clears and fills trail_->MutableConflict() with the literals of the arcs + // that form a cycle containing the given node. + void FillConflictFromCircuitAt(int start); + + const int num_nodes_; + const Options options_; + Trail* trail_; + + // Internal representation of the graph given at construction. Const. + struct Arc { + int tail; + int head; + }; + std::vector self_arcs_; + + std::vector watch_index_to_literal_; + std::vector> watch_index_to_arcs_; + + // Index in trail_ up to which we propagated all the assigned Literals. + int propagation_trail_index_; + + // Current partial chains of arc that are present. + std::vector next_; // -1 if not assigned yet. + std::vector prev_; // -1 if not assigned yet. + std::vector next_literal_; + + // Backtrack support for the partial chains of arcs, level_ends_[level] is an + // index in added_arcs_; + std::vector level_ends_; + std::vector added_arcs_; + + // Temporary vector. + std::vector in_circuit_; + + DISALLOW_COPY_AND_ASSIGN(CircuitPropagator); +}; + +// ============================================================================ +// Model based functions. +// ============================================================================ + +// Enforces that exactly one literal per rows and per columns is true. +// This only work for a square matrix (but could easily be generalized). +// +// If ignore_row_and_column_zero is true, this adds two less constraints by +// skipping the one for the row zero and for the column zero. Note however that +// the other constraints are not changed, i.e. matrix[0][5] is still counted +// in column 5. +std::function ExactlyOnePerRowAndPerColumn( + const std::vector>& square_matrix, + bool ignore_row_and_column_zero = false); + +inline std::function SubcircuitConstraint( + const std::vector>& graph) { + return [=](Model* model) { + if (graph.empty()) return; + model->Add(ExactlyOnePerRowAndPerColumn(graph)); + CircuitPropagator::Options options; + CircuitPropagator* constraint = + new CircuitPropagator(graph, options, model->GetOrCreate()); + constraint->RegisterWith(model->GetOrCreate()); + model->TakeOwnership(constraint); + }; +} + +inline std::function MultipleSubcircuitThroughZeroConstraint( + const std::vector>& graph) { + return [=](Model* model) { + if (graph.empty()) return; + model->Add(ExactlyOnePerRowAndPerColumn( + graph, /*ignore_row_and_column_zero=*/true)); + CircuitPropagator::Options options; + options.multiple_subcircuit_through_zero = true; + CircuitPropagator* constraint = + new CircuitPropagator(graph, options, model->GetOrCreate()); + constraint->RegisterWith(model->GetOrCreate()); + model->TakeOwnership(constraint); + }; +} + +} // namespace sat +} // namespace operations_research + +#endif // OR_TOOLS_SAT_CIRCUIT_H_ diff --git a/ortools/sat/clause.cc b/ortools/sat/clause.cc index 293a509cec..3aa5a754f4 100644 --- a/ortools/sat/clause.cc +++ b/ortools/sat/clause.cc @@ -536,7 +536,7 @@ void BinaryImplicationGraph::MinimizeConflictExperimental( void BinaryImplicationGraph::RemoveFixedVariables( int first_unprocessed_trail_index, const Trail& trail) { - const VariablesAssignment& assigment = trail.Assignment(); + const VariablesAssignment& assignment = trail.Assignment(); SCOPED_TIME_STAT(&stats_); is_marked_.ClearAndResize(LiteralIndex(implications_.size())); for (int i = first_unprocessed_trail_index; i < trail.Index(); ++i) { @@ -553,9 +553,9 @@ void BinaryImplicationGraph::RemoveFixedVariables( STLClearObject(&(implications_[true_literal.NegatedIndex()])); } for (const LiteralIndex i : is_marked_.PositionsSetAtLeastOnce()) { - RemoveIf(&implications_[i], - std::bind1st(std::mem_fun(&VariablesAssignment::LiteralIsTrue), - &assigment)); + RemoveIf(&implications_[i], [&assignment](const Literal& lit) { + return assignment.LiteralIsTrue(lit); + }); } } diff --git a/ortools/sat/cp_constraints.cc b/ortools/sat/cp_constraints.cc index 3ee279728e..63ec130b49 100644 --- a/ortools/sat/cp_constraints.cc +++ b/ortools/sat/cp_constraints.cc @@ -76,210 +76,6 @@ void BooleanXorPropagator::RegisterWith(GenericLiteralWatcher* watcher) { } } -CircuitPropagator::CircuitPropagator( - const std::vector>& graph, bool allow_subcircuit, - Trail* trail) - : num_nodes_(graph.size()), - allow_subcircuit_(allow_subcircuit), - trail_(trail), - propagation_trail_index_(0) { - // TODO(user): add a way to properly handle trivially UNSAT cases. - // For now we just check that they don't occur at construction. - CHECK_GT(num_nodes_, 1) - << "Trivial or UNSAT constraint, shouldn't be constructed!"; - next_.resize(num_nodes_, -1); - prev_.resize(num_nodes_, -1); - next_literal_.resize(num_nodes_); - self_arcs_.resize(num_nodes_); - std::unordered_map literal_to_watch_index; - const VariablesAssignment& assignment = trail->Assignment(); - for (int tail = 0; tail < num_nodes_; ++tail) { - self_arcs_[tail] = graph[tail][tail]; - for (int head = 0; head < num_nodes_; ++head) { - const LiteralIndex index = graph[tail][head]; - // Note that we need to test for both "special" cases before we can - // call assignment.LiteralIsTrue() or LiteralIsFalse(). - if (index == kFalseLiteralIndex) continue; - if (index == kTrueLiteralIndex || - assignment.LiteralIsTrue(Literal(index))) { - if (!allow_subcircuit) { - CHECK_NE(tail, head) << "Trivially UNSAT."; - } - CHECK_EQ(next_[tail], -1) - << "Trivially UNSAT or duplicate arcs while adding " << tail - << " -> " << head; - CHECK_EQ(prev_[head], -1) - << "Trivially UNSAT or duplicate arcs while adding " << tail - << " -> " << head; - next_[tail] = head; - prev_[head] = tail; - next_literal_[tail] = kNoLiteralIndex; - continue; - } - if (assignment.LiteralIsFalse(Literal(index))) continue; - - int watch_index_ = FindWithDefault(literal_to_watch_index, index, -1); - if (watch_index_ == -1) { - watch_index_ = watch_index_to_literal_.size(); - literal_to_watch_index[index] = watch_index_; - watch_index_to_literal_.push_back(Literal(index)); - watch_index_to_arcs_.push_back(std::vector()); - } - watch_index_to_arcs_[watch_index_].push_back({tail, head}); - } - } -} - -void CircuitPropagator::SetLevel(int level) { - if (level == level_ends_.size()) return; - if (level > level_ends_.size()) { - while (level > level_ends_.size()) { - level_ends_.push_back(added_arcs_.size()); - } - return; - } - - // Backtrack. - for (int i = level_ends_[level]; i < added_arcs_.size(); ++i) { - const Arc arc = added_arcs_[i]; - next_[arc.tail] = -1; - prev_[arc.head] = -1; - } - added_arcs_.resize(level_ends_[level]); - level_ends_.resize(level); -} - -void CircuitPropagator::FillConflictFromCircuitAt(int start) { - std::vector* conflict = trail_->MutableConflict(); - conflict->clear(); - int node = start; - do { - CHECK_NE(node, -1); - if (next_literal_[node] != kNoLiteralIndex) { - conflict->push_back(Literal(next_literal_[node]).Negated()); - } - node = next_[node]; - } while (node != start); -} - -bool CircuitPropagator::Propagate() { return true; } - -bool CircuitPropagator::IncrementalPropagate( - const std::vector& watch_indices) { - for (const int w : watch_indices) { - const Literal literal = watch_index_to_literal_[w]; - for (const Arc arc : watch_index_to_arcs_[w]) { - // Get rid of the trivial conflicts: - // - At most one incoming and one ougtoing arc for each nodes. - // - No self arc if allow_subcircuit_ is false - std::vector* conflict = trail_->MutableConflict(); - if (arc.tail == arc.head && !allow_subcircuit_) { - *conflict = {literal.Negated()}; - return false; - } - if (next_[arc.tail] != -1) { - if (next_literal_[arc.tail] != kNoLiteralIndex) { - *conflict = {Literal(next_literal_[arc.tail]).Negated(), - literal.Negated()}; - } else { - *conflict = {literal.Negated()}; - } - return false; - } - if (prev_[arc.head] != -1) { - if (next_literal_[prev_[arc.head]] != kNoLiteralIndex) { - *conflict = {Literal(next_literal_[prev_[arc.head]]).Negated(), - literal.Negated()}; - } else { - *conflict = {literal.Negated()}; - } - return false; - } - - // Add the arc. - added_arcs_.push_back(arc); - next_[arc.tail] = arc.head; - prev_[arc.head] = arc.tail; - next_literal_[arc.tail] = literal.Index(); - - // Circuit? - if (allow_subcircuit_) { - in_circuit_.assign(num_nodes_, false); - in_circuit_[arc.tail] = true; - } - int node = arc.head; - int size = 1; - while (node != arc.tail && node != -1) { - if (allow_subcircuit_) in_circuit_[node] = true; - node = next_[node]; - size++; - } - if (node != arc.tail) continue; - - // We have one circuit. - if (size == num_nodes_) return true; - if (!allow_subcircuit_) { - FillConflictFromCircuitAt(arc.tail); - return false; - } - - if (size == 1) continue; // self-arc. - - // HACK: we can reuse the conflict vector even though we don't have a - // conflict. - FillConflictFromCircuitAt(arc.tail); - BooleanVariable variable_with_same_reason = kNoBooleanVariable; - - // We can propagate all the other nodes to point to themselves. - // If this is not already the case, we have a conflict. - for (int node = 0; node < num_nodes_; ++node) { - if (in_circuit_[node] || next_[node] == node) continue; - if (next_[node] != -1) { - std::vector* conflict = trail_->MutableConflict(); - if (next_literal_[node] != kNoLiteralIndex) { - conflict->push_back(Literal(next_literal_[node]).Negated()); - } - return false; - } else if (self_arcs_[node] == kFalseLiteralIndex) { - return false; - } else { - DCHECK_NE(self_arcs_[node], kTrueLiteralIndex); - const Literal literal(self_arcs_[node]); - - // We may not have processed this literal yet. - if (trail_->Assignment().LiteralIsTrue(literal)) continue; - if (trail_->Assignment().LiteralIsFalse(literal)) { - std::vector* conflict = trail_->MutableConflict(); - conflict->push_back(literal); - return false; - } - - // Propagate. - if (variable_with_same_reason == kNoBooleanVariable) { - variable_with_same_reason = literal.Variable(); - const int index = trail_->Index(); - trail_->Enqueue(literal, AssignmentType::kCachedReason); - *trail_->GetVectorToStoreReason(index) = *trail_->MutableConflict(); - trail_->NotifyThatReasonIsCached(literal.Variable()); - } else { - trail_->EnqueueWithSameReasonAs(literal, variable_with_same_reason); - } - } - } - } - } - return true; -} - -void CircuitPropagator::RegisterWith(GenericLiteralWatcher* watcher) { - const int id = watcher->Register(this); - for (int w = 0; w < watch_index_to_literal_.size(); ++w) { - watcher->WatchLiteral(watch_index_to_literal_[w], id, w); - } - watcher->RegisterReversibleClass(id, this); - watcher->RegisterReversibleInt(id, &propagation_trail_index_); -} - // ----- CpPropagator ----- CpPropagator::CpPropagator(IntegerTrail* integer_trail) @@ -685,50 +481,5 @@ void OneOfVarMinPropagator::RegisterWith(GenericLiteralWatcher* watcher) { for (const IntegerVariable v : vars_) watcher->WatchLowerBound(v, id); } -std::function ExactlyOnePerRowAndPerColumn( - const std::vector>& square_matrix) { - return [=](Model* model) { - int n = square_matrix.size(); - int num_trivially_false = 0; - Trail* trail = model->GetOrCreate(); - std::vector exactly_one_constraint; - for (const bool transpose : {false, true}) { - for (int i = 0; i < n; ++i) { - int num_true = 0; - exactly_one_constraint.clear(); - for (int j = 0; j < n; ++j) { - CHECK_EQ(n, square_matrix[i].size()); - const LiteralIndex index = - transpose ? square_matrix[j][i] : square_matrix[i][j]; - if (index == kFalseLiteralIndex) continue; - if (index == kTrueLiteralIndex) { - ++num_true; - continue; - } - exactly_one_constraint.push_back(Literal(index)); - } - if (num_true > 1) { - LOG(WARNING) << "UNSAT in ExactlyOnePerRowAndPerColumn()."; - return model->GetOrCreate()->NotifyThatModelIsUnsat(); - } - CHECK_LE(num_true, 1); - if (num_true == 1) { - for (const Literal l : exactly_one_constraint) { - if (!trail->Assignment().VariableIsAssigned(l.Variable())) { - ++num_trivially_false; - trail->EnqueueWithUnitReason(l.Negated()); - } - } - } else { - model->Add(ExactlyOneConstraint(exactly_one_constraint)); - } - } - } - if (num_trivially_false > 0) { - LOG(INFO) << "Num extra fixed literal: " << num_trivially_false; - } - }; -} - } // namespace sat } // namespace operations_research diff --git a/ortools/sat/cp_constraints.h b/ortools/sat/cp_constraints.h index 7e07db3ade..232e36224d 100644 --- a/ortools/sat/cp_constraints.h +++ b/ortools/sat/cp_constraints.h @@ -57,60 +57,6 @@ class BooleanXorPropagator : public PropagatorInterface { DISALLOW_COPY_AND_ASSIGN(BooleanXorPropagator); }; -// Initial version of the circuit/sub-circuit constraint. -// It mainly report conflicts and do not propagate much. -class CircuitPropagator : PropagatorInterface, ReversibleInterface { - public: - // The constraints take a dense representation of a graph on [0, n). Each arc - // being present when the given literal is true. The special values - // kTrueLiteralIndex and kFalseLiteralIndex can be used for arcs that are - // either always there or never there. - CircuitPropagator(const std::vector>& graph, - bool allow_subcircuit, Trail* trail); - - void SetLevel(int level) final; - bool Propagate() final; - bool IncrementalPropagate(const std::vector& watch_indices) final; - void RegisterWith(GenericLiteralWatcher* watcher); - - private: - // Clears and fills trail_->MutableConflict() with the literals of the arcs - // that form a cycle containing the given node. - void FillConflictFromCircuitAt(int start); - - const int num_nodes_; - const bool allow_subcircuit_; - Trail* trail_; - - // Internal representation of the graph given at construction. Const. - struct Arc { - int tail; - int head; - }; - std::vector self_arcs_; - - std::vector watch_index_to_literal_; - std::vector> watch_index_to_arcs_; - - // Index in trail_ up to which we propagated all the assigned Literals. - int propagation_trail_index_; - - // Current partial chains of arc that are present. - std::vector next_; // -1 if not assigned yet. - std::vector prev_; // -1 if not assigned yet. - std::vector next_literal_; - - // Backtrack support for the partial chains of arcs, level_ends_[level] is an - // index in added_arcs_; - std::vector level_ends_; - std::vector added_arcs_; - - // Temporary vector. - std::vector in_circuit_; - - DISALLOW_COPY_AND_ASSIGN(CircuitPropagator); -}; - // Base class to help writing CP inspired constraints. class CpPropagator : public PropagatorInterface { public: @@ -338,35 +284,6 @@ inline std::function StrictNonOverlappingFixedSizeRectangles( }; } -// Enforces that exactly one literal per rows and per columns is true. -// This only work for a square matrix (but could easily be generalized). -std::function ExactlyOnePerRowAndPerColumn( - const std::vector>& square_matrix); - -inline std::function CircuitConstraint( - const std::vector>& graph) { - return [=](Model* model) { - if (graph.empty()) return; - model->Add(ExactlyOnePerRowAndPerColumn(graph)); - CircuitPropagator* constraint = new CircuitPropagator( - graph, /*allow_subcircuit=*/false, model->GetOrCreate()); - constraint->RegisterWith(model->GetOrCreate()); - model->TakeOwnership(constraint); - }; -} - -inline std::function SubcircuitConstraint( - const std::vector>& graph) { - return [=](Model* model) { - if (graph.empty()) return; - model->Add(ExactlyOnePerRowAndPerColumn(graph)); - CircuitPropagator* constraint = new CircuitPropagator( - graph, /*allow_subcircuit=*/true, model->GetOrCreate()); - constraint->RegisterWith(model->GetOrCreate()); - model->TakeOwnership(constraint); - }; -} - // The target variable is equal to exactly one of the candidate variable. The // equality is controlled by the given "selector" literals. // diff --git a/ortools/sat/cp_model.proto b/ortools/sat/cp_model.proto index 1e44912e91..5bb2458c99 100644 --- a/ortools/sat/cp_model.proto +++ b/ortools/sat/cp_model.proto @@ -152,6 +152,27 @@ message CircuitConstraintProto { repeated int32 nexts = 2; } +// The "VRP" constraint. +// +// The direct graph where arc #i (from tails[i] to head[i]) is present iff +// literals[i] is true must satisfy this set of properties: +// - #incoming arcs == 1 except for node 0. +// - #outgoing arcs == 1 except for node 0. +// - for node zero, #incoming arcs == #outgoing arcs. +// - There is no duplicate arcs. +// - Self-arc are allowed except for node 0. +// - There is no cycle in this graph, except through node 0. +// +// TODO(user): It is probably possible to generalize this constraint to a +// no-cycle in a general graph, or a no-cycle with sum incoming <= 1 and sum +// outgoing <= 1 (more efficient implementation). On the other hand, having this +// specific constraint allow us to add specific "cuts" to a VRP problem. +message RoutesConstraintProto { + repeated int32 tails = 1; + repeated int32 heads = 2; + repeated int32 literals = 3; +} + // The values of the n-tuple formed by the given variables can only be one of // the listed n-tuples in values. The n-tuples are encoded in a flattened way: // [tuple0_v0, tuple0_v1, ..., tuple0_v{n-1}, tuple1_v0, ...]. @@ -219,6 +240,7 @@ message ConstraintProto { AllDifferentConstraintProto all_diff = 13; ElementConstraintProto element = 14; CircuitConstraintProto circuit = 15; + RoutesConstraintProto routes = 23; TableConstraintProto table = 16; AutomataConstraintProto automata = 17; InverseConstraintProto inverse = 18; diff --git a/ortools/sat/cp_model_checker.cc b/ortools/sat/cp_model_checker.cc index 9b88ddcb56..3e5b305ef2 100644 --- a/ortools/sat/cp_model_checker.cc +++ b/ortools/sat/cp_model_checker.cc @@ -522,6 +522,64 @@ class ConstraintChecker { return num_visited + num_inactive == num_nodes; } + bool RoutesConstraintIsFeasible(const CpModelProto& model, + const ConstraintProto& ct) { + const int num_arcs = ct.routes().tails_size(); + int num_used_arcs = 0; + int num_self_arcs = 0; + int num_nodes = 0; + std::vector tail_to_head; + std::vector depot_nexts; + for (int i = 0; i < num_arcs; ++i) { + const int tail = ct.routes().tails(i); + const int head = ct.routes().heads(i); + num_nodes = std::max(num_nodes, 1 + tail); + num_nodes = std::max(num_nodes, 1 + head); + tail_to_head.resize(num_nodes, -1); + if (LiteralIsTrue(ct.routes().literals(i))) { + if (tail == head) { + if (tail == 0) return false; + ++num_self_arcs; + continue; + } + ++num_used_arcs; + if (tail == 0) { + depot_nexts.push_back(head); + } else { + if (tail_to_head[tail] != -1) return false; + tail_to_head[tail] = head; + } + } + } + + // Make sure each routes from the depot go back to it, and count such arcs. + int count = 0; + for (int start : depot_nexts) { + ++count; + while (start != 0) { + if (tail_to_head[start] == -1) return false; + start = tail_to_head[start]; + ++count; + } + } + + if (count != num_used_arcs) { + VLOG(1) << "count: " << count << " != num_used_arcs:" << num_used_arcs; + return false; + } + + // Each routes cover as many node as there is arcs, but this way we count + // multiple times the depot. So the number of nodes covered are: + // count - depot_nexts.size() + 1. + // And this number + the self arcs should be num_nodes. + if (count - depot_nexts.size() + 1 + num_self_arcs != num_nodes) { + VLOG(1) << "Not all nodes are covered!"; + return false; + } + + return true; + } + bool InverseConstraintIsFeasible(const CpModelProto& model, const ConstraintProto& ct) { const int num_variables = ct.inverse().f_direct_size(); @@ -621,6 +679,9 @@ bool SolutionIsFeasible(const CpModelProto& model, case ConstraintProto::ConstraintCase::kCircuit: is_feasible = checker.CircuitConstraintIsFeasible(model, ct); break; + case ConstraintProto::ConstraintCase::kRoutes: + is_feasible = checker.RoutesConstraintIsFeasible(model, ct); + break; case ConstraintProto::ConstraintCase::kInverse: is_feasible = checker.InverseConstraintIsFeasible(model, ct); break; diff --git a/ortools/sat/cp_model_solver.cc b/ortools/sat/cp_model_solver.cc index 552464417d..b7b8bf8020 100644 --- a/ortools/sat/cp_model_solver.cc +++ b/ortools/sat/cp_model_solver.cc @@ -37,6 +37,7 @@ #include "ortools/base/stl_util.h" #include "ortools/graph/connectivity.h" #include "ortools/sat/all_different.h" +#include "ortools/sat/circuit.h" #include "ortools/sat/cp_constraints.h" #include "ortools/sat/cp_model_checker.h" #include "ortools/sat/cp_model_presolve.h" @@ -1239,6 +1240,26 @@ void LoadCircuitConstraint(const ConstraintProto& ct, ModelWithMapping* m) { m->Add(SubcircuitConstraint(graph)); } +// TODO(user): provide a sparse API. +void LoadRoutesConstraint(const ConstraintProto& ct, ModelWithMapping* m) { + int num_nodes = 0; + const RoutesConstraintProto& arg = ct.routes(); + for (const int32 tail : arg.tails()) { + num_nodes = std::max(num_nodes, tail + 1); + } + for (const int32 head : arg.heads()) { + num_nodes = std::max(num_nodes, head + 1); + } + std::vector> graph( + num_nodes, std::vector(num_nodes, kFalseLiteralIndex)); + + const int num_arcs = arg.tails_size(); + for (int i = 0; i < num_arcs; ++i) { + graph[arg.tails(i)][arg.heads(i)] = m->Literal(arg.literals(i)).Index(); + } + m->Add(MultipleSubcircuitThroughZeroConstraint(graph)); +} + void LoadInverseConstraint(const ConstraintProto& ct, ModelWithMapping* m) { // Fully encode both arrays of variables, encode the constraint using Boolean // equalities: f_direct[i] == j <=> f_inverse[j] == i. @@ -1504,6 +1525,9 @@ bool LoadConstraint(const ConstraintProto& ct, ModelWithMapping* m) { case ConstraintProto::ConstraintProto::kCircuit: LoadCircuitConstraint(ct, m); return true; + case ConstraintProto::ConstraintProto::kRoutes: + LoadRoutesConstraint(ct, m); + return true; case ConstraintProto::ConstraintProto::kInverse: LoadInverseConstraint(ct, m); return true; diff --git a/ortools/sat/cp_model_utils.cc b/ortools/sat/cp_model_utils.cc index bb1f13f115..eacfbf9310 100644 --- a/ortools/sat/cp_model_utils.cc +++ b/ortools/sat/cp_model_utils.cc @@ -71,6 +71,9 @@ void AddReferencesUsedByConstraint(const ConstraintProto& ct, case ConstraintProto::ConstraintCase::kCircuit: AddIndices(ct.circuit().nexts(), &output->variables); break; + case ConstraintProto::ConstraintCase::kRoutes: + AddIndices(ct.routes().literals(), &output->literals); + break; case ConstraintProto::ConstraintCase::kInverse: AddIndices(ct.inverse().f_direct(), &output->variables); AddIndices(ct.inverse().f_inverse(), &output->variables); @@ -147,6 +150,9 @@ void ApplyToAllLiteralIndices(const std::function& f, break; case ConstraintProto::ConstraintCase::kCircuit: break; + case ConstraintProto::ConstraintCase::kRoutes: + APPLY_TO_REPEATED_FIELD(routes, literals); + break; case ConstraintProto::ConstraintCase::kInverse: break; case ConstraintProto::ConstraintCase::kTable: @@ -209,6 +215,8 @@ void ApplyToAllVariableIndices(const std::function& f, case ConstraintProto::ConstraintCase::kCircuit: APPLY_TO_REPEATED_FIELD(circuit, nexts); break; + case ConstraintProto::ConstraintCase::kRoutes: + break; case ConstraintProto::ConstraintCase::kInverse: APPLY_TO_REPEATED_FIELD(inverse, f_direct); APPLY_TO_REPEATED_FIELD(inverse, f_inverse); @@ -264,6 +272,8 @@ void ApplyToAllIntervalIndices(const std::function& f, break; case ConstraintProto::ConstraintCase::kCircuit: break; + case ConstraintProto::ConstraintCase::kRoutes: + break; case ConstraintProto::ConstraintCase::kInverse: break; case ConstraintProto::ConstraintCase::kTable: @@ -316,6 +326,8 @@ std::string ConstraintCaseName(ConstraintProto::ConstraintCase constraint_case) return "kElement"; case ConstraintProto::ConstraintCase::kCircuit: return "kCircuit"; + case ConstraintProto::ConstraintCase::kRoutes: + return "kRoutes"; case ConstraintProto::ConstraintCase::kInverse: return "kInverse"; case ConstraintProto::ConstraintCase::kTable: diff --git a/ortools/sat/no_cycle.cc b/ortools/sat/no_cycle.cc deleted file mode 100644 index 1ad1adab38..0000000000 --- a/ortools/sat/no_cycle.cc +++ /dev/null @@ -1,329 +0,0 @@ -// Copyright 2010-2014 Google -// Licensed under the Apache License, Version 2.0 (the "License"); -// you may not use this file except in compliance with the License. -// You may obtain a copy of the License at -// -// http://www.apache.org/licenses/LICENSE-2.0 -// -// Unless required by applicable law or agreed to in writing, software -// distributed under the License is distributed on an "AS IS" BASIS, -// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. -// See the License for the specific language governing permissions and -// limitations under the License. - -#include "ortools/sat/no_cycle.h" - -#include -#include - -#include "ortools/base/logging.h" - -namespace operations_research { -namespace sat { - -bool NoCyclePropagator::Propagate(Trail* trail) { - if (problem_is_unsat_) { - trail->MutableConflict()->clear(); - return false; - } - - if (!initialization_is_done_) { - CHECK_EQ(0, trail->CurrentDecisionLevel()); - initialization_is_done_ = true; - - // Propagate all that can be propagated using the fixed arcs. - for (int node = 0; node < graph_.size(); ++node) { - FillNodeIsReachedWithAntecedentOf(node); - for (const Arc arc : potential_graph_[node]) { - if (!node_is_reached_[arc.head]) continue; - - // We know that l must be false. - const Literal l(arc.literal_index); - if (trail->Assignment().VariableIsAssigned(l.Variable())) { - if (trail->Assignment().LiteralIsTrue(l)) { - // The problem is UNSAT. - trail->MutableConflict()->clear(); - return false; - } - } else { - trail->EnqueueWithUnitReason(l.Negated()); - } - } - } - } - - while (propagation_trail_index_ < trail->Index()) { - const Literal literal = (*trail)[propagation_trail_index_++]; - if (literal.Index() >= potential_arcs_.size()) continue; - for (const std::pair& p : potential_arcs_[literal.Index()]) { - { - // Remove this arc from the potential_graph_. - // - // TODO(user): this is not super efficient, but it helps speeding up - // the propagation here and in the makespan constraint. - int new_size = 0; - std::vector& ref = potential_graph_[p.first]; - const int size = ref.size(); - while (new_size < size && - ref[new_size].literal_index != literal.Index()) { - ++new_size; - } - for (int i = new_size; i < size; ++i) { - if (ref[i].literal_index != literal.Index()) { - ref[new_size++] = ref[i]; - } - } - ref.resize(new_size); - } - - if (!include_propagated_arcs_in_graph_ && - trail->AssignmentType(literal.Variable()) == propagator_id_) { - continue; - } - - // Warning: The order of the following 3 lines matter! - if (IsReachable(p.first, p.second)) continue; - newly_reachable_ = NewlyReachable(p.second, &node_is_reached_); - FillNodeIsReachedWithAntecedentOf(p.first); - - // Do nothing if we reached the threshold on the number of arcs. - if (num_arcs_ == num_arcs_threshold_) continue; - - if (node_is_reached_[p.second]) { // Conflict. - // Note that this modifies node_is_reached_ and reached_nodes_, but - // since we abort aftewards, it is fine. - FindReasonForPath(*trail, p.second, p.first, propagation_trail_index_, - trail->MutableConflict()); - trail->MutableConflict()->push_back(literal.Negated()); - return false; - } - - ++num_arcs_; - graph_[p.first].push_back({p.second, literal.Index()}); - reverse_graph_[p.second].push_back({p.first, literal.Index()}); - - for (const int node : newly_reachable_) { - for (const Arc arc : potential_graph_[node]) { - CHECK_NE(arc.literal_index, kNoLiteralIndex); - if (node_is_reached_[arc.head]) { - const Literal l(arc.literal_index); - if (trail->Assignment().VariableIsAssigned(l.Variable())) { - // TODO(user): we could detect a conflict earlier if the literal - // l is already assigned to true. - continue; - } - - // Save the information needed for the lazy-explanation and enqueue - // the fact that "arc" cannot be in the graph. - const int trail_index = trail->Index(); - if (trail_index >= reason_arc_.size()) { - reason_arc_.resize(trail_index + 1); - reason_trail_limit_.resize(trail_index + 1); - } - reason_arc_[trail_index] = {node, arc.head}; - reason_trail_limit_[trail_index] = propagation_trail_index_; - trail->Enqueue(l.Negated(), propagator_id_); - } - } - } - } - } - return true; -} - -void NoCyclePropagator::Untrail(const Trail& trail, int trail_index) { - while (propagation_trail_index_ > trail_index) { - const Literal literal = trail[--propagation_trail_index_]; - if (literal.Index() >= potential_arcs_.size()) continue; - for (const std::pair& p : potential_arcs_[literal.Index()]) { - DCHECK_LT(p.first, graph_.size()); - DCHECK_GE(p.first, 0); - - potential_graph_[p.first].push_back({p.second, literal.Index()}); - - // We only remove this arc if it was added. That is if it is the last arc - // in graph_[p.first]. - if (graph_[p.first].empty()) continue; - if (graph_[p.first].back().literal_index != literal.Index()) continue; - - --num_arcs_; - graph_[p.first].pop_back(); - reverse_graph_[p.second].pop_back(); - } - } -} - -// TODO(user): If one literal propagate many arcs, and more than one is needed -// to form a cycle, this will not work properly. -gtl::Span NoCyclePropagator::Reason(const Trail& trail, - int trail_index) const { - const int source = reason_arc_[trail_index].second; - const int target = reason_arc_[trail_index].first; - const int trail_limit = reason_trail_limit_[trail_index]; - std::vector* const reason = - trail.GetVectorToStoreReason(trail_index); - - // Note that this modify node_is_reached_ and reached_nodes_. - FindReasonForPath(trail, source, target, trail_limit, reason); - return *reason; -} - -namespace { - -// This sets the given vector of Booleans to all false using a vector of its -// positions at true in order to exploit sparsity. -// -// TODO(user): Add a test depending on position.size() or use directly -// util::operations_research::SparseBitset. -void ResetBitsetWithPosition(int new_size, std::vector* bitset, - std::vector* true_positions) { - bitset->resize(new_size, false); - for (const int i : *true_positions) { - DCHECK((*bitset)[i]); - (*bitset)[i] = false; - } - true_positions->clear(); - DCHECK( - std::none_of(bitset->begin(), bitset->end(), [](bool v) { return v; })); -} - -} // namespace - -// We use a BFS to try to minimize the reason. -void NoCyclePropagator::FindReasonForPath(const Trail& trail, int source, - int target, int trail_limit, - std::vector* reason) const { - CHECK_NE(source, target); - ResetBitsetWithPosition(graph_.size(), &node_is_reached_, &reached_nodes_); - - // This is the same code as IsReachable() below, except that we need to - // remember the path taken to the target and we work on a subgraph. - reached_nodes_.push_back(source); - parent_index_with_literal_.clear(); - parent_index_with_literal_.push_back({0, LiteralIndex(-1)}); - node_is_reached_[source] = true; - - int i = 0; - for (; i < reached_nodes_.size(); ++i) { - const int node = reached_nodes_[i]; - if (node == target) break; - - // Only consider arc whose literal was assigned before trail_limit. The arcs - // in graph_[node] are ordered by increasing trail index, so it is okay to - // abort as soon as an arc was added after trail_limit. - for (const Arc& arc : graph_[node]) { - if (arc.literal_index != kNoLiteralIndex) { - const BooleanVariable var = Literal(arc.literal_index).Variable(); - if (trail.Info(var).trail_index >= trail_limit) break; - } - if (node_is_reached_[arc.head]) continue; - node_is_reached_[arc.head] = true; - reached_nodes_.push_back(arc.head); - parent_index_with_literal_.push_back({i, arc.literal_index}); - } - } - - // Follow the path backward and fill the reason. - CHECK_LT(i, reached_nodes_.size()) << "The target is not reachable!"; - reason->clear(); - while (i != 0) { - const LiteralIndex literal_index = parent_index_with_literal_[i].second; - if (literal_index != kNoLiteralIndex) { - reason->push_back(Literal(literal_index).Negated()); - } - i = parent_index_with_literal_[i].first; - } -} - -bool NoCyclePropagator::IsReachable(int source, int destination) const { - if (source == destination) return true; - ResetBitsetWithPosition(graph_.size(), &node_is_reached_, &reached_nodes_); - reached_nodes_.push_back(source); - node_is_reached_[source] = true; - for (int i = 0; i < reached_nodes_.size(); ++i) { - for (const Arc& arc : graph_[reached_nodes_[i]]) { - if (arc.head == destination) return true; - if (node_is_reached_[arc.head]) continue; - node_is_reached_[arc.head] = true; - reached_nodes_.push_back(arc.head); - } - } - return false; -} - -void NoCyclePropagator::FillNodeIsReachedWithAntecedentOf(int source) const { - ResetBitsetWithPosition(graph_.size(), &node_is_reached_, &reached_nodes_); - reached_nodes_.push_back(source); - node_is_reached_[source] = true; - for (int i = 0; i < reached_nodes_.size(); ++i) { - for (const Arc& arc : reverse_graph_[reached_nodes_[i]]) { - if (node_is_reached_[arc.head]) continue; - node_is_reached_[arc.head] = true; - reached_nodes_.push_back(arc.head); - } - } -} - -std::vector NoCyclePropagator::NewlyReachable( - int source, std::vector* already_reached) const { - if ((*already_reached)[source]) return std::vector(); - std::vector result; - result.push_back(source); - (*already_reached)[source] = true; - for (int i = 0; i < result.size(); ++i) { - for (const Arc& arc : graph_[result[i]]) { - if ((*already_reached)[arc.head]) continue; - (*already_reached)[arc.head] = true; - result.push_back(arc.head); - } - } - - // Restore already_reached to its original value. - for (const int node : result) (*already_reached)[node] = false; - return result; -} - -void NoCyclePropagator::AdjustSizes(int tail, int head, int literal_index) { - CHECK_NE(tail, head); - CHECK(!initialization_is_done_); - CHECK_EQ(0, propagation_trail_index_); - - const int m = std::max(tail, head); - DCHECK_GE(m, 0); - if (m >= graph_.size()) { - graph_.resize(m + 1); - potential_graph_.resize(m + 1); - reverse_graph_.resize(m + 1); - } - DCHECK_GE(literal_index, 0); - if (literal_index >= potential_arcs_.size()) { - potential_arcs_.resize(literal_index + 1); - } -} - -void NoCyclePropagator::AddArc(int tail, int head) { - AdjustSizes(tail, head, 0); - - // Deal with the corner case of a cycle with the fixed arcs. - if (problem_is_unsat_ || IsReachable(head, tail)) { - problem_is_unsat_ = true; - return; - } - - // TODO(user): test IsReachable(tail, head) and do not add the arc if true? - ++num_arcs_; - graph_[tail].push_back({head, LiteralIndex(-1)}); - reverse_graph_[head].push_back({tail, LiteralIndex(-1)}); -} - -void NoCyclePropagator::AddPotentialArc(int tail, int head, Literal literal) { - AdjustSizes(tail, head, literal.Index().value()); - potential_arcs_[literal.Index()].push_back({tail, head}); - potential_graph_[tail].push_back({head, literal.Index()}); - CHECK_EQ(1, potential_arcs_[literal.Index()].size()) - << "We don't support multiple arcs associated to the same literal. " - << "However, it should be fairly easy to support this case."; -} - -} // namespace sat -} // namespace operations_research diff --git a/ortools/sat/no_cycle.h b/ortools/sat/no_cycle.h deleted file mode 100644 index dc0d3e4782..0000000000 --- a/ortools/sat/no_cycle.h +++ /dev/null @@ -1,183 +0,0 @@ -// Copyright 2010-2014 Google -// Licensed under the Apache License, Version 2.0 (the "License"); -// you may not use this file except in compliance with the License. -// You may obtain a copy of the License at -// -// http://www.apache.org/licenses/LICENSE-2.0 -// -// Unless required by applicable law or agreed to in writing, software -// distributed under the License is distributed on an "AS IS" BASIS, -// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. -// See the License for the specific language governing permissions and -// limitations under the License. - -#ifndef OR_TOOLS_SAT_NO_CYCLE_H_ -#define OR_TOOLS_SAT_NO_CYCLE_H_ - -#include -#include -#include - -#include "ortools/base/integral_types.h" -#include "ortools/base/macros.h" -#include "ortools/base/span.h" -#include "ortools/base/int_type.h" -#include "ortools/base/int_type_indexed_vector.h" -#include "ortools/sat/sat_base.h" - -namespace operations_research { -namespace sat { - -// The "no-cycle" constraint. -// -// Each arc will be associated to a literal and this propagator will make sure -// that there is no cycle in the graph with only the arcs whose associated -// literal is set to true. -class NoCyclePropagator : public SatPropagator { - public: - NoCyclePropagator() - : SatPropagator("NoCyclePropagator"), - num_arcs_(0), - problem_is_unsat_(false), - initialization_is_done_(false), - num_arcs_threshold_(std::numeric_limits::max()), - include_propagated_arcs_in_graph_(true) {} - ~NoCyclePropagator() final {} - - bool Propagate(Trail* trail) final; - void Untrail(const Trail& trail, int trail_index) final; - gtl::Span Reason(const Trail& trail, - int trail_index) const final; - - // Stops doing anything when the number of arcs in the graph becomes greater - // that the given value. This allows to use this class to model a circuit - // constraint on n nodes: we don't want any cycle, but it is okay to have one - // when we add the n-th arc. Of course we also need to make sure that each - // node as an unique successor using at-most-one constraints. - void AllowCycleWhenNumArcsIsGreaterThan(int64 value) { - num_arcs_threshold_ = value; - } - - // If this is false, then we don't track inside our graphs the arcs that we - // propagated. This is meant to be turned on if an arc and its reverse are - // controlled by a literal and its negation. When this is the case, then we - // know that all the arcs propagated by this class don't change the - // reachability of the graph. - void IncludePropagatedArcsInGraph(bool value) { - include_propagated_arcs_in_graph_ = value; - } - - // Adds a "constant" arc to the graph. - // Self-arc are not allowed (it would create a trivial cycle). - void AddArc(int tail, int head); - - // Registers an arc that will be present in the graph iff 'literal' is true. - // Self-arc are not allowed (it would fix the given literal to false). - // - // TODO(user): support more than one arc associated to the same literal. - void AddPotentialArc(int tail, int head, Literal literal); - - // Getters for the current graph. This is only in sync with the trail iff - // SatPropagator::PropagationIsDone() is true. - // - // Note that these graphs will NOT contain all the arcs but will correctly - // encode the reachability of every node. More specifically, when an arc (tail - // -> head) is about to be added but a path from tail to head already exists - // in the graph, this arc will not be added. - struct Arc { - int head; - LiteralIndex literal_index; - }; - const std::vector>& Graph() const { return graph_; } - const std::vector>& ReverseGraph() const { - return reverse_graph_; - } - - // Getters for the "potential" arcs. That is the arcs that could be added to - // the graph or not depending on their associated literal value. Note that - // some already added arcs may not appear here for optimization purposes. - const std::vector>& PotentialGraph() const { - return potential_graph_; - } - const ITIVector>>& - PotentialArcs() const { - return potential_arcs_; - } - - private: - // Adjust the internal data structures when a new arc is added. - void AdjustSizes(int tail, int head, int literal_index); - - // Returns true if destination is reachable from source in graph_. - // Warning: this modifies node_is_reached_ and reached_nodes_. - bool IsReachable(int source, int destination) const; - - // Returns the set of node from which source can be reached (included). - // Warning: this modifies node_is_reached_ and reached_nodes_. - void FillNodeIsReachedWithAntecedentOf(int source) const; - - // Returns the vector of node that are reachable from source (included), but - // not in the given already_reached. The already_reached vector is not const - // because this function temporarily modifies it before restoring it to its - // original value for performance reason. - std::vector NewlyReachable(int source, - std::vector* already_reached) const; - - // Finds a path from source to target and output its reason. - // Only the arcs whose associated literal is assigned before the given - // trail_limit are considered. - // - // Warning: this modifies node_is_reached_ and reached_nodes_. - void FindReasonForPath(const Trail& trail, int source, int target, - int trail_limit, std::vector* reason) const; - - // The number of arcs in graph_ and reverse_graph_. - int64 num_arcs_; - - // Just used to detect the corner case of a cycle with fixed arcs. - bool problem_is_unsat_; - bool initialization_is_done_; - - // Control the options of this class. - int64 num_arcs_threshold_; - bool include_propagated_arcs_in_graph_; - - // The current graph wich is kept in sync with the literal trail. For each - // node, graph_[node] list the pair (head, literal_index) of the outgoing - // arcs. - // - // Important: this will always be kept acyclic. - std::vector> graph_; - std::vector> reverse_graph_; - - // The graph formed by all the potential arcs in the same format as graph_. - std::vector> potential_graph_; - - // The set of potential arc (tail, head) indexed by literal_index. - // - // TODO(user): Introduce a struct with .tail and .head to make the code more - // readable. - ITIVector>> potential_arcs_; - - // Temporary vectors used by the various BFS computations. We always have: - // node_is_reached_[node] is true iff reached_nodes_ contains node. - mutable std::vector reached_nodes_; - mutable std::vector node_is_reached_; - - // Temporary vector to hold the result of NewlyReachable(). - std::vector newly_reachable_; - - // Temporary vector used by FindReasonForPath(). - mutable std::vector> parent_index_with_literal_; - - // The arc that caused the literal at a given trail_index to be propagated. - std::vector> reason_arc_; - std::vector reason_trail_limit_; - - DISALLOW_COPY_AND_ASSIGN(NoCyclePropagator); -}; - -} // namespace sat -} // namespace operations_research - -#endif // OR_TOOLS_SAT_NO_CYCLE_H_ diff --git a/ortools/sat/python/cp_model.py b/ortools/sat/python/cp_model.py index 18be30faea..1373eed63a 100644 --- a/ortools/sat/python/cp_model.py +++ b/ortools/sat/python/cp_model.py @@ -34,6 +34,7 @@ from ortools.sat import pywrapsat INT_MIN = -9223372036854775808 # hardcoded to be platform independent. INT_MAX = 9223372036854775807 + # Cp Solver status (exported to avoid importing cp_model_cp2). UNKNOWN = cp_model_pb2.UNKNOWN MODEL_INVALID = cp_model_pb2.MODEL_INVALID diff --git a/ortools/sat/sat_solver.h b/ortools/sat/sat_solver.h index fa460e4d20..9eaecbbe94 100644 --- a/ortools/sat/sat_solver.h +++ b/ortools/sat/sat_solver.h @@ -636,7 +636,7 @@ class SatSolver { // as decaying all the variable activities, but it is a lot more efficient. void UpdateVariableActivityIncrement(); - // Activity managment for clauses. This work the same way at the ones for + // Activity management for clauses. This work the same way at the ones for // variables, but with different parameters. void BumpReasonActivities(const std::vector& literals); void BumpClauseActivity(SatClause* clause);