From 5cef83b2180009304ac6f28dd895e8af5dc97bfb Mon Sep 17 00:00:00 2001 From: "laurent.perron@gmail.com" Date: Fri, 2 May 2014 10:33:52 +0000 Subject: [PATCH] add previous flatzinc constraints; store flattening mapping --- makefiles/Makefile.cpp.mk | 4 + src/flatzinc2/constraints.cc | 83 +-- src/flatzinc2/flatzinc_constraints.cc | 733 ++++++++++++++++++++++++++ src/flatzinc2/flatzinc_constraints.h | 46 ++ src/flatzinc2/presolve.cc | 56 +- src/flatzinc2/presolve.h | 21 + 6 files changed, 888 insertions(+), 55 deletions(-) create mode 100644 src/flatzinc2/flatzinc_constraints.cc create mode 100644 src/flatzinc2/flatzinc_constraints.h diff --git a/makefiles/Makefile.cpp.mk b/makefiles/Makefile.cpp.mk index 9d2c831c02..2905e32c92 100644 --- a/makefiles/Makefile.cpp.mk +++ b/makefiles/Makefile.cpp.mk @@ -905,6 +905,7 @@ endif FLATZINC2_LIB_OBJS=\ $(OBJ_DIR)/flatzinc2/constraints.$O\ + $(OBJ_DIR)/flatzinc2/flatzinc_constraints.$O\ $(OBJ_DIR)/flatzinc2/model.$O\ $(OBJ_DIR)/flatzinc2/parallel_support.$O\ $(OBJ_DIR)/flatzinc2/parser.$O\ @@ -926,6 +927,9 @@ $(GEN_DIR)/flatzinc2/parser.tab.hh: $(GEN_DIR)/flatzinc2/parser.tab.cc $(OBJ_DIR)/flatzinc2/constraints.$O:$(SRC_DIR)/flatzinc2/constraints.cc $(SRC_DIR)/flatzinc2/model.h $(SRC_DIR)/flatzinc2/solver.h $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sflatzinc2$Sconstraints.cc $(OBJ_OUT)$(OBJ_DIR)$Sflatzinc2$Sconstraints.$O +$(OBJ_DIR)/flatzinc2/flatzinc_constraints.$O:$(SRC_DIR)/flatzinc2/flatzinc_constraints.cc $(SRC_DIR)/flatzinc2/model.h $(SRC_DIR)/flatzinc2/solver.h + $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sflatzinc2$Sflatzinc_constraints.cc $(OBJ_OUT)$(OBJ_DIR)$Sflatzinc2$Sflatzinc_constraints.$O + $(OBJ_DIR)/flatzinc2/model.$O:$(SRC_DIR)/flatzinc2/model.cc $(SRC_DIR)/flatzinc2/model.h $(CCC) $(CFLAGS) -c $(SRC_DIR)$Sflatzinc2$Smodel.cc $(OBJ_OUT)$(OBJ_DIR)$Sflatzinc2$Smodel.$O diff --git a/src/flatzinc2/constraints.cc b/src/flatzinc2/constraints.cc index 7c70e1cc16..63027064c6 100644 --- a/src/flatzinc2/constraints.cc +++ b/src/flatzinc2/constraints.cc @@ -15,6 +15,7 @@ #include "base/integral_types.h" #include "base/logging.h" #include "base/hash.h" +#include "flatzinc2/flatzinc_constraints.h" #include "flatzinc2/model.h" #include "flatzinc2/search.h" #include "flatzinc2/solver.h" @@ -388,17 +389,17 @@ void ExtractIntLinEq(FzSolver* fzsolver, FzConstraint* ct) { fzsolver->SetExtracted(ct->target_variable, target); } } else { - Constraint* ct = nullptr; + Constraint* constraint = nullptr; switch (size) { case 0: { - ct = rhs == 0 ? solver->MakeTrueConstraint() - : solver->MakeFalseConstraint(); + constraint = rhs == 0 ? solver->MakeTrueConstraint() + : solver->MakeFalseConstraint(); break; } case 1: { IntExpr* const e1 = fzsolver->Extract(fzvars[0]); const int64 c1 = coefficients[0]; - ct = solver->MakeEquality(solver->MakeProd(e1, c1), rhs); + constraint = solver->MakeEquality(solver->MakeProd(e1, c1), rhs); break; } case 2: { @@ -408,20 +409,20 @@ void ExtractIntLinEq(FzSolver* fzsolver, FzConstraint* ct) { const int64 c2 = coefficients[1]; if (c1 > 0) { if (c2 > 0) { - ct = solver->MakeEquality( + constraint = solver->MakeEquality( solver->MakeProd(e1, c1), solver->MakeDifference(rhs, solver->MakeProd(e2, c2))); } else { - ct = solver->MakeEquality( + constraint = solver->MakeEquality( solver->MakeProd(e1, c1), solver->MakeSum(solver->MakeProd(e2, -c2), rhs)); } } else if (c2 > 0) { - ct = solver->MakeEquality( + constraint = solver->MakeEquality( solver->MakeProd(e2, c2), solver->MakeSum(solver->MakeProd(e1, -c1), rhs)); } else { - ct = solver->MakeEquality( + constraint = solver->MakeEquality( solver->MakeProd(e1, -c1), solver->MakeDifference(-rhs, solver->MakeProd(e2, -c2))); } @@ -434,31 +435,40 @@ void ExtractIntLinEq(FzSolver* fzsolver, FzConstraint* ct) { const int64 c1 = coefficients[0]; const int64 c2 = coefficients[1]; const int64 c3 = coefficients[2]; - if (c1 < 0 && c2 > 0 && c3 > 0) { - ct = solver->MakeEquality( - solver->MakeSum(solver->MakeProd(e2, c2), - solver->MakeProd(e3, c3)), - solver->MakeSum(solver->MakeProd(e1, -c1), rhs)); - } else if (c1 > 0 && c2 < 0 && c3 > 0) { - ct = solver->MakeEquality( - solver->MakeSum(solver->MakeProd(e1, c1), - solver->MakeProd(e3, c3)), - solver->MakeSum(solver->MakeProd(e2, -c2), rhs)); - } else if (c1 > 0 && c2 > 0 && c3 < 0) { - ct = solver->MakeEquality( - solver->MakeSum(solver->MakeProd(e1, c1), - solver->MakeProd(e2, c2)), - solver->MakeSum(solver->MakeProd(e3, -c3), rhs)); - } else if (c1 < 0 && c2 < 0 && c3 > 0) { - ct = solver->MakeEquality( - solver->MakeSum(solver->MakeProd(e1, -c1), - solver->MakeProd(e2, -c2)), - solver->MakeSum(solver->MakeProd(e3, c3), -rhs)); + if (ct->strong_propagation) { + std::vector variables; + variables.push_back(e1->Var()); + variables.push_back(e2->Var()); + variables.push_back(e3->Var()); + constraint = + MakeStrongScalProdEquality(solver, variables, coefficients, rhs); } else { - ct = solver->MakeEquality( - solver->MakeSum(solver->MakeProd(e1, c1), - solver->MakeProd(e2, c2)), - solver->MakeDifference(rhs, solver->MakeProd(e3, c3))); + if (c1 < 0 && c2 > 0 && c3 > 0) { + constraint = solver->MakeEquality( + solver->MakeSum(solver->MakeProd(e2, c2), + solver->MakeProd(e3, c3)), + solver->MakeSum(solver->MakeProd(e1, -c1), rhs)); + } else if (c1 > 0 && c2 < 0 && c3 > 0) { + constraint = solver->MakeEquality( + solver->MakeSum(solver->MakeProd(e1, c1), + solver->MakeProd(e3, c3)), + solver->MakeSum(solver->MakeProd(e2, -c2), rhs)); + } else if (c1 > 0 && c2 > 0 && c3 < 0) { + constraint = solver->MakeEquality( + solver->MakeSum(solver->MakeProd(e1, c1), + solver->MakeProd(e2, c2)), + solver->MakeSum(solver->MakeProd(e3, -c3), rhs)); + } else if (c1 < 0 && c2 < 0 && c3 > 0) { + constraint = solver->MakeEquality( + solver->MakeSum(solver->MakeProd(e1, -c1), + solver->MakeProd(e2, -c2)), + solver->MakeSum(solver->MakeProd(e3, c3), -rhs)); + } else { + constraint = solver->MakeEquality( + solver->MakeSum(solver->MakeProd(e1, c1), + solver->MakeProd(e2, c2)), + solver->MakeDifference(rhs, solver->MakeProd(e3, c3))); + } } break; } @@ -481,12 +491,13 @@ void ExtractIntLinEq(FzSolver* fzsolver, FzConstraint* ct) { // PostBooleanSumInRange(model, spec, variables, rhs, rhs); // return; // } else { - ct = solver->MakeScalProdEquality(variables, new_coefficients, rhs); - // } + constraint = + solver->MakeScalProdEquality(variables, new_coefficients, rhs); + // } } } - FZVLOG << " - posted " << ct->DebugString() << FZENDL; - solver->AddConstraint(ct); + FZVLOG << " - posted " << constraint->DebugString() << FZENDL; + solver->AddConstraint(constraint); } } diff --git a/src/flatzinc2/flatzinc_constraints.cc b/src/flatzinc2/flatzinc_constraints.cc new file mode 100644 index 0000000000..3ad1dec448 --- /dev/null +++ b/src/flatzinc2/flatzinc_constraints.cc @@ -0,0 +1,733 @@ +// Copyright 2010-2013 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 "flatzinc2/flatzinc_constraints.h" + +#include "base/commandlineflags.h" +#include "constraint_solver/constraint_solveri.h" +#include "util/string_array.h" + +DECLARE_bool(cp_trace_search); +DECLARE_bool(cp_trace_propagation); +DECLARE_bool(use_sat); + +namespace operations_research { +namespace { +class BooleanSumOdd : public Constraint { + public: + BooleanSumOdd(Solver* const s, const std::vector& vars) + : Constraint(s), + vars_(vars), + num_possible_true_vars_(0), + num_always_true_vars_(0) {} + + virtual ~BooleanSumOdd() {} + + virtual void Post() { + for (int i = 0; i < vars_.size(); ++i) { + if (!vars_[i]->Bound()) { + Demon* const u = MakeConstraintDemon1( + solver(), this, &BooleanSumOdd::Update, "Update", i); + vars_[i]->WhenBound(u); + } + } + } + + virtual void InitialPropagate() { + int num_always_true = 0; + int num_possible_true = 0; + int possible_true_index = -1; + for (int i = 0; i < vars_.size(); ++i) { + const IntVar* const var = vars_[i]; + if (var->Min() == 1) { + num_always_true++; + num_possible_true++; + } else if (var->Max() == 1) { + num_possible_true++; + possible_true_index = i; + } + } + if (num_always_true == num_possible_true && num_possible_true % 2 == 0) { + solver()->Fail(); + } else if (num_possible_true == num_always_true + 1) { + DCHECK_NE(-1, possible_true_index); + if (num_possible_true % 2 == 1) { + vars_[possible_true_index]->SetMin(1); + } else { + vars_[possible_true_index]->SetMax(0); + } + } + num_possible_true_vars_.SetValue(solver(), num_possible_true); + num_always_true_vars_.SetValue(solver(), num_always_true); + } + + void Update(int index) { + DCHECK(vars_[index]->Bound()); + const int64 value = vars_[index]->Min(); // Faster than Value(). + if (value == 0) { + num_possible_true_vars_.Decr(solver()); + } else { + DCHECK_EQ(1, value); + num_always_true_vars_.Incr(solver()); + } + if (num_always_true_vars_.Value() == num_possible_true_vars_.Value() && + num_possible_true_vars_.Value() % 2 == 0) { + solver()->Fail(); + } else if (num_possible_true_vars_.Value() == + num_always_true_vars_.Value() + 1) { + int possible_true_index = -1; + for (int i = 0; i < vars_.size(); ++i) { + if (!vars_[i]->Bound()) { + possible_true_index = i; + break; + } + } + if (possible_true_index != -1) { + if (num_possible_true_vars_.Value() % 2 == 1) { + vars_[possible_true_index]->SetMin(1); + } else { + vars_[possible_true_index]->SetMax(0); + } + } + } + } + + virtual std::string DebugString() const { + return StringPrintf("BooleanSumOdd([%s])", + JoinDebugStringPtr(vars_, ", ").c_str()); + } + + virtual void Accept(ModelVisitor* const visitor) const { + visitor->BeginVisitConstraint(ModelVisitor::kSumEqual, this); + visitor->VisitIntegerVariableArrayArgument(ModelVisitor::kVarsArgument, + vars_); + visitor->EndVisitConstraint(ModelVisitor::kSumEqual, this); + } + + private: + const std::vector vars_; + NumericalRev num_possible_true_vars_; + NumericalRev num_always_true_vars_; +}; + +class VariableParity : public Constraint { + public: + VariableParity(Solver* const s, IntVar* const var, bool odd) + : Constraint(s), var_(var), odd_(odd) {} + + virtual ~VariableParity() {} + + virtual void Post() { + if (!var_->Bound()) { + Demon* const u = solver()->MakeConstraintInitialPropagateCallback(this); + var_->WhenRange(u); + } + } + + virtual void InitialPropagate() { + const int64 vmax = var_->Max(); + const int64 vmin = var_->Min(); + int64 new_vmax = vmax; + int64 new_vmin = vmin; + if (odd_) { + if (vmax % 2 == 0) { + new_vmax--; + } + if (vmin % 2 == 0) { + new_vmin++; + } + } else { + if (vmax % 2 == 1) { + new_vmax--; + } + if (vmin % 2 == 1) { + new_vmin++; + } + } + var_->SetRange(new_vmin, new_vmax); + } + + virtual std::string DebugString() const { + return StringPrintf("VarParity(%s, %d)", var_->DebugString().c_str(), odd_); + } + + virtual void Accept(ModelVisitor* const visitor) const { + visitor->BeginVisitConstraint("VarParity", this); + visitor->VisitIntegerExpressionArgument(ModelVisitor::kVariableArgument, + var_); + visitor->VisitIntegerArgument(ModelVisitor::kValuesArgument, odd_); + visitor->EndVisitConstraint("VarParity", this); + } + + private: + IntVar* const var_; + const bool odd_; +}; + +class IsBooleanSumInRange : public Constraint { + public: + IsBooleanSumInRange(Solver* const s, const std::vector& vars, + int64 range_min, int64 range_max, IntVar* const target) + : Constraint(s), + vars_(vars), + range_min_(range_min), + range_max_(range_max), + target_(target), + num_possible_true_vars_(0), + num_always_true_vars_(0) {} + + virtual ~IsBooleanSumInRange() {} + + virtual void Post() { + for (int i = 0; i < vars_.size(); ++i) { + if (!vars_[i]->Bound()) { + Demon* const u = MakeConstraintDemon1( + solver(), this, &IsBooleanSumInRange::Update, "Update", i); + vars_[i]->WhenBound(u); + } + } + if (!target_->Bound()) { + Demon* const u = MakeConstraintDemon0( + solver(), this, &IsBooleanSumInRange::UpdateTarget, "UpdateTarget"); + target_->WhenBound(u); + } + } + + virtual void InitialPropagate() { + int num_always_true = 0; + int num_possible_true = 0; + int possible_true_index = -1; + for (int i = 0; i < vars_.size(); ++i) { + const IntVar* const var = vars_[i]; + if (var->Min() == 1) { + num_always_true++; + num_possible_true++; + } else if (var->Max() == 1) { + num_possible_true++; + possible_true_index = i; + } + } + num_possible_true_vars_.SetValue(solver(), num_possible_true); + num_always_true_vars_.SetValue(solver(), num_always_true); + UpdateTarget(); + } + + void UpdateTarget() { + if (num_always_true_vars_.Value() > range_max_ || + num_possible_true_vars_.Value() < range_min_) { + inactive_.Switch(solver()); + target_->SetValue(0); + } else if (num_always_true_vars_.Value() >= range_min_ && + num_possible_true_vars_.Value() <= range_max_) { + inactive_.Switch(solver()); + target_->SetValue(1); + } else if (target_->Min() == 1) { + if (num_possible_true_vars_.Value() == range_min_) { + PushAllUnboundToOne(); + } else if (num_always_true_vars_.Value() == range_max_) { + PushAllUnboundToZero(); + } + } else if (target_->Max() == 0) { + if (num_possible_true_vars_.Value() == range_max_ + 1 && + num_always_true_vars_.Value() >= range_min_) { + PushAllUnboundToOne(); + } else if (num_always_true_vars_.Value() == range_min_ - 1 && + num_possible_true_vars_.Value() <= range_max_) { + PushAllUnboundToZero(); + } + } + } + + void Update(int index) { + if (!inactive_.Switched()) { + DCHECK(vars_[index]->Bound()); + const int64 value = vars_[index]->Min(); // Faster than Value(). + if (value == 0) { + num_possible_true_vars_.Decr(solver()); + } else { + DCHECK_EQ(1, value); + num_always_true_vars_.Incr(solver()); + } + UpdateTarget(); + } + } + + virtual std::string DebugString() const { + return StringPrintf("Sum([%s]) in [%" GG_LL_FORMAT "d..%" GG_LL_FORMAT + "d] == %s", + JoinDebugStringPtr(vars_, ", ").c_str(), range_min_, + range_max_, target_->DebugString().c_str()); + } + + virtual void Accept(ModelVisitor* const visitor) const { + visitor->BeginVisitConstraint(ModelVisitor::kSumEqual, this); + visitor->VisitIntegerVariableArrayArgument(ModelVisitor::kVarsArgument, + vars_); + visitor->EndVisitConstraint(ModelVisitor::kSumEqual, this); + } + + private: + void PushAllUnboundToZero() { + inactive_.Switch(solver()); + int true_vars = 0; + for (int i = 0; i < vars_.size(); ++i) { + if (vars_[i]->Min() == 0) { + vars_[i]->SetValue(0); + } else { + true_vars++; + } + } + target_->SetValue(true_vars >= range_min_ && true_vars <= range_max_); + } + + void PushAllUnboundToOne() { + inactive_.Switch(solver()); + int true_vars = 0; + for (int i = 0; i < vars_.size(); ++i) { + if (vars_[i]->Max() == 1) { + vars_[i]->SetValue(1); + true_vars++; + } + } + target_->SetValue(true_vars >= range_min_ && true_vars <= range_max_); + } + + const std::vector vars_; + const int64 range_min_; + const int64 range_max_; + IntVar* const target_; + NumericalRev num_possible_true_vars_; + NumericalRev num_always_true_vars_; + RevSwitch inactive_; +}; + +class BooleanSumInRange : public Constraint { + public: + BooleanSumInRange(Solver* const s, const std::vector& vars, + int64 range_min, int64 range_max) + : Constraint(s), + vars_(vars), + range_min_(range_min), + range_max_(range_max), + num_possible_true_vars_(0), + num_always_true_vars_(0) {} + + virtual ~BooleanSumInRange() {} + + virtual void Post() { + for (int i = 0; i < vars_.size(); ++i) { + if (!vars_[i]->Bound()) { + Demon* const u = MakeConstraintDemon1( + solver(), this, &BooleanSumInRange::Update, "Update", i); + vars_[i]->WhenBound(u); + } + } + } + + virtual void InitialPropagate() { + int num_always_true = 0; + int num_possible_true = 0; + int possible_true_index = -1; + for (int i = 0; i < vars_.size(); ++i) { + const IntVar* const var = vars_[i]; + if (var->Min() == 1) { + num_always_true++; + num_possible_true++; + } else if (var->Max() == 1) { + num_possible_true++; + possible_true_index = i; + } + } + num_possible_true_vars_.SetValue(solver(), num_possible_true); + num_always_true_vars_.SetValue(solver(), num_always_true); + Check(); + } + + void Check() { + if (num_always_true_vars_.Value() > range_max_ || + num_possible_true_vars_.Value() < range_min_) { + solver()->Fail(); + } else if (num_always_true_vars_.Value() >= range_min_ && + num_possible_true_vars_.Value() <= range_max_) { + // Inhibit. + } else { + if (num_possible_true_vars_.Value() == range_min_) { + PushAllUnboundToOne(); + } else if (num_always_true_vars_.Value() == range_max_) { + PushAllUnboundToZero(); + } + } + } + + void Update(int index) { + DCHECK(vars_[index]->Bound()); + const int64 value = vars_[index]->Min(); // Faster than Value(). + if (value == 0) { + num_possible_true_vars_.Decr(solver()); + } else { + DCHECK_EQ(1, value); + num_always_true_vars_.Incr(solver()); + } + Check(); + } + + virtual std::string DebugString() const { + return StringPrintf("Sum([%s]) in [%" GG_LL_FORMAT "d..%" GG_LL_FORMAT "d]", + JoinDebugStringPtr(vars_, ", ").c_str(), range_min_, + range_max_); + } + + virtual void Accept(ModelVisitor* const visitor) const { + visitor->BeginVisitConstraint(ModelVisitor::kSumEqual, this); + visitor->VisitIntegerVariableArrayArgument(ModelVisitor::kVarsArgument, + vars_); + visitor->EndVisitConstraint(ModelVisitor::kSumEqual, this); + } + + private: + void PushAllUnboundToZero() { + for (int i = 0; i < vars_.size(); ++i) { + if (vars_[i]->Min() == 0) { + vars_[i]->SetValue(0); + } else { + } + } + } + + void PushAllUnboundToOne() { + for (int i = 0; i < vars_.size(); ++i) { + if (vars_[i]->Max() == 1) { + vars_[i]->SetValue(1); + } + } + } + + const std::vector vars_; + const int64 range_min_; + const int64 range_max_; + NumericalRev num_possible_true_vars_; + NumericalRev num_always_true_vars_; +}; + +// ----- Inverse constraint ----- + +// Variable demand cumulative time table + +class VariableCumulativeTask : public BaseObject { + public: + VariableCumulativeTask(IntVar* const start, IntVar* const duration, + IntVar* const demand) + : start_(start), duration_(duration), demand_(demand) {} + // Copy constructor. Cannot be explicit, because we want to pass instances of + // VariableCumulativeTask by copy. + VariableCumulativeTask(const VariableCumulativeTask& task) + : start_(task.start_), duration_(task.duration_), demand_(task.demand_) {} + + IntVar* start() const { return start_; } + IntVar* duration() const { return duration_; } + IntVar* demand() const { return demand_; } + int64 StartMin() const { return start_->Min(); } + int64 StartMax() const { return start_->Max(); } + int64 EndMin() const { return start_->Min() + duration_->Min(); } + virtual std::string DebugString() const { + return StringPrintf("Task{ start: %s, duration: %s, demand: %s }", + start_->DebugString().c_str(), + duration_->DebugString().c_str(), + demand_->DebugString().c_str()); + } + + private: + IntVar* const start_; + IntVar* const duration_; + IntVar* const demand_; +}; + +struct ProfileDelta { + ProfileDelta(int64 _time, int64 _delta) : time(_time), delta(_delta) {} + int64 time; + int64 delta; +}; + +bool TimeLessThan(const ProfileDelta& delta1, const ProfileDelta& delta2) { + return delta1.time < delta2.time; +} + +bool TaskStartMinLessThan(const VariableCumulativeTask* const task1, + const VariableCumulativeTask* const task2) { + return (task1->StartMin() < task2->StartMin()); +} + +class VariableCumulativeTimeTable : public Constraint { + public: + VariableCumulativeTimeTable(Solver* const solver, + const std::vector& tasks, + IntVar* const capacity) + : Constraint(solver), by_start_min_(tasks), capacity_(capacity) { + // There may be up to 2 delta's per interval (one on each side), + // plus two sentinels + const int profile_max_size = 2 * NumTasks() + 2; + profile_non_unique_time_.reserve(profile_max_size); + profile_unique_time_.reserve(profile_max_size); + } + + virtual void InitialPropagate() { + BuildProfile(); + PushTasks(); + } + + virtual void Post() { + Demon* const demon = + solver()->MakeDelayedConstraintInitialPropagateCallback(this); + for (int i = 0; i < NumTasks(); ++i) { + by_start_min_[i]->start()->WhenRange(demon); + by_start_min_[i]->duration()->WhenRange(demon); + by_start_min_[i]->demand()->WhenRange(demon); + } + capacity_->WhenRange(demon); + } + + int NumTasks() const { return by_start_min_.size(); } + + void Accept(ModelVisitor* const visitor) const { + LOG(FATAL) << "Should not be visited"; + } + + virtual std::string DebugString() const { + return StringPrintf("VariableCumulativeTimeTable([%s], capacity = %s)", + JoinDebugStringPtr(by_start_min_, ", ").c_str(), + capacity_->DebugString().c_str()); + } + + private: + // Build the usage profile. Runs in O(n log n). + void BuildProfile() { + // Build profile with non unique time + profile_non_unique_time_.clear(); + for (int i = 0; i < NumTasks(); ++i) { + const VariableCumulativeTask* task = by_start_min_[i]; + const int64 start_max = task->StartMax(); + const int64 end_min = task->EndMin(); + const int64 demand_min = task->demand()->Min(); + if (start_max < end_min && demand_min > 0) { + profile_non_unique_time_.push_back( + ProfileDelta(start_max, +demand_min)); + profile_non_unique_time_.push_back(ProfileDelta(end_min, -demand_min)); + } + } + // Sort + std::sort(profile_non_unique_time_.begin(), profile_non_unique_time_.end(), + TimeLessThan); + // Build profile with unique times + profile_unique_time_.clear(); + profile_unique_time_.push_back(ProfileDelta(kint64min, 0)); + for (int i = 0; i < profile_non_unique_time_.size(); ++i) { + const ProfileDelta& profile_delta = profile_non_unique_time_[i]; + if (profile_delta.time == profile_unique_time_.back().time) { + profile_unique_time_.back().delta += profile_delta.delta; + } else { + profile_unique_time_.push_back(profile_delta); + } + } + // Re-scan profile to compute min usage, max usage, and check + // final usage to be 0. + int64 usage = 0; + int64 max_required_usage = 0; + const int64 max_capacity = capacity_->Max(); + for (int i = 0; i < profile_unique_time_.size(); ++i) { + const ProfileDelta& profile_delta = profile_unique_time_[i]; + usage += profile_delta.delta; + max_required_usage = std::max(max_required_usage, usage); + if (usage > max_capacity) { + solver()->Fail(); + } + } + DCHECK_EQ(0, usage); + profile_unique_time_.push_back(ProfileDelta(kint64max, 0)); + + // Propagate on capacity. + capacity_->SetMin(max_required_usage); + } + + // Update the start min for all tasks. Runs in O(n^2) and Omega(n). + void PushTasks() { + std::sort(by_start_min_.begin(), by_start_min_.end(), TaskStartMinLessThan); + int64 usage = 0; + int profile_index = 0; + for (int task_index = 0; task_index < NumTasks(); ++task_index) { + VariableCumulativeTask* const task = by_start_min_[task_index]; + if (task->duration()->Min() > 0) { + while (task->StartMin() > profile_unique_time_[profile_index].time) { + DCHECK(profile_index < profile_unique_time_.size()); + ++profile_index; + usage += profile_unique_time_[profile_index].delta; + } + PushTask(task, profile_index, usage); + } + } + } + + // Push the given task to new_start_min, defined as the smallest integer such + // that the profile usage for all tasks, excluding the current one, does not + // exceed capacity_ - task->demand()->Min() on the interval + // [new_start_min, new_start_min + task->interval()->DurationMin() ). + void PushTask(VariableCumulativeTask* const task, int profile_index, + int64 usage) { + // Init + const int64 demand_min = task->demand()->Min(); + const int64 adjusted_demand = demand_min == 0 ? 1 : demand_min; + const bool is_adjusted = demand_min == 0; + const int64 residual_capacity = capacity_->Max() - adjusted_demand; + const int64 duration_min = task->duration()->Min(); + const ProfileDelta& first_prof_delta = profile_unique_time_[profile_index]; + + int64 new_start_min = task->StartMin(); + + DCHECK_GE(first_prof_delta.time, task->StartMin()); + // The check above is with a '>='. Let's first treat the '>' case + if (first_prof_delta.time > task->StartMin()) { + // There was no profile delta at a time between interval->StartMin() + // (included) and the current one. + // As we don't delete delta's of 0 value, this means the current task + // does not contribute to the usage before: + DCHECK((task->StartMax() >= first_prof_delta.time) || + (task->StartMax() >= task->EndMin())); + // The 'usage' given in argument is valid at first_prof_delta.time. To + // compute the usage at the start min, we need to remove the last delta. + const int64 usage_at_start_min = usage - first_prof_delta.delta; + if (usage_at_start_min > residual_capacity) { + new_start_min = profile_unique_time_[profile_index].time; + } + } + + // Influence of current task + const int64 start_max = task->StartMax(); + const int64 end_min = task->EndMin(); + // TODO(user): remove delta_start and delta_end. + ProfileDelta delta_start(start_max, 0); + ProfileDelta delta_end(end_min, 0); + if (start_max < end_min) { + delta_start.delta = +demand_min; + delta_end.delta = -demand_min; + } + while (profile_unique_time_[profile_index].time < + duration_min + new_start_min) { + const ProfileDelta& profile_delta = profile_unique_time_[profile_index]; + DCHECK(profile_index < profile_unique_time_.size()); + // Compensate for current task + if (profile_delta.time == delta_start.time) { + usage -= delta_start.delta; + } + if (profile_delta.time == delta_end.time) { + usage -= delta_end.delta; + } + // Increment time + ++profile_index; + DCHECK(profile_index < profile_unique_time_.size()); + // Does it fit? + if (usage > residual_capacity) { + new_start_min = profile_unique_time_[profile_index].time; + } + usage += profile_unique_time_[profile_index].delta; + } + if (is_adjusted) { + if (new_start_min > task->StartMax()) { + task->demand()->SetMax(0); + } + } else { + task->start()->SetMin(new_start_min); + } + } + + typedef std::vector Profile; + + Profile profile_unique_time_; + Profile profile_non_unique_time_; + std::vector by_start_min_; + IntVar* const capacity_; + + DISALLOW_COPY_AND_ASSIGN(VariableCumulativeTimeTable); +}; +} // namespace + +Constraint* MakeIsBooleanSumInRange(Solver* const solver, + const std::vector& variables, + int64 range_min, int64 range_max, + IntVar* const target) { + return solver->RevAlloc( + new IsBooleanSumInRange(solver, variables, range_min, range_max, target)); +} + +Constraint* MakeBooleanSumInRange(Solver* const solver, + const std::vector& variables, + int64 range_min, int64 range_max) { + return solver->RevAlloc( + new BooleanSumInRange(solver, variables, range_min, range_max)); +} +Constraint* MakeBooleanSumOdd(Solver* const solver, + const std::vector& variables) { + return solver->RevAlloc(new BooleanSumOdd(solver, variables)); +} + +Constraint* MakeStrongScalProdEquality(Solver* const solver, + const std::vector& variables, + const std::vector& coefficients, + int64 rhs) { + const bool trace = FLAGS_cp_trace_search; + const bool propag = FLAGS_cp_trace_propagation; + FLAGS_cp_trace_search = false; + FLAGS_cp_trace_propagation = false; + const int size = variables.size(); + IntTupleSet tuples(size); + Solver s("build"); + std::vector copy_vars(size); + for (int i = 0; i < size; ++i) { + copy_vars[i] = s.MakeIntVar(variables[i]->Min(), variables[i]->Max()); + } + s.AddConstraint(s.MakeScalProdEquality(copy_vars, coefficients, rhs)); + s.NewSearch(s.MakePhase(copy_vars, Solver::CHOOSE_FIRST_UNBOUND, + Solver::ASSIGN_MIN_VALUE)); + while (s.NextSolution()) { + std::vector one_tuple(size); + for (int i = 0; i < size; ++i) { + one_tuple[i] = copy_vars[i]->Value(); + } + tuples.Insert(one_tuple); + } + s.EndSearch(); + FLAGS_cp_trace_search = trace; + FLAGS_cp_trace_propagation = propag; + return solver->MakeAllowedAssignments(variables, tuples); +} + +Constraint* MakeVariableCumulative(Solver* const solver, + const std::vector& starts, + const std::vector& durations, + const std::vector& usages, + IntVar* const capacity) { + std::vector tasks(starts.size()); + for (int i = 0; i < starts.size(); ++i) { + tasks[i] = solver->RevAlloc( + new VariableCumulativeTask(starts[i], durations[i], usages[i])); + } + return solver->RevAlloc( + new VariableCumulativeTimeTable(solver, tasks, capacity)); +} + +Constraint* MakeVariableOdd(Solver* const s, IntVar* const var) { + return s->RevAlloc(new VariableParity(s, var, true)); +} + +Constraint* MakeVariableEven(Solver* const s, IntVar* const var) { + return s->RevAlloc(new VariableParity(s, var, false)); +} +} // namespace operations_research diff --git a/src/flatzinc2/flatzinc_constraints.h b/src/flatzinc2/flatzinc_constraints.h new file mode 100644 index 0000000000..8d975de70c --- /dev/null +++ b/src/flatzinc2/flatzinc_constraints.h @@ -0,0 +1,46 @@ +// Copyright 2010-2013 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_FLATZINC2_FLATZINC_CONSTRAINTS_H_ +#define OR_TOOLS_FLATZINC2_FLATZINC_CONSTRAINTS_H_ + +#include "constraint_solver/constraint_solver.h" + +namespace operations_research { +Constraint* MakeStrongScalProdEquality(Solver* const solver, + const std::vector& variables, + const std::vector& coefficients, + int64 rhs); + +Constraint* MakeIsBooleanSumInRange(Solver* const solver, + const std::vector& variables, + int64 range_min, int64 range_max, + IntVar* const target); + +Constraint* MakeBooleanSumInRange(Solver* const solver, + const std::vector& variables, + int64 range_min, int64 range_max); + +Constraint* MakeBooleanSumOdd(Solver* const solver, + const std::vector& variables); + +Constraint* MakeVariableCumulative(Solver* const solver, + const std::vector& starts, + const std::vector& durations, + const std::vector& usages, + IntVar* const capacity); + +Constraint* MakeVariableOdd(Solver* const s, IntVar* const var); +Constraint* MakeVariableEven(Solver* const s, IntVar* const var); +} // namespace operations_research +#endif // OR_TOOLS_FLATZINC_FLATZINC_CONSTRAINTS_H_ diff --git a/src/flatzinc2/presolve.cc b/src/flatzinc2/presolve.cc index fb2de483bf..5350e8cfae 100644 --- a/src/flatzinc2/presolve.cc +++ b/src/flatzinc2/presolve.cc @@ -233,13 +233,13 @@ bool FzPresolver::PresolveIntDiv(FzConstraint* ct) { bool FzPresolver::PresolveArrayBoolOr(FzConstraint* ct) { if (!ct->presolve_propagation_done && ct->Arg(1).HasOneValue() && ct->Arg(1).Value() == 0) { - for (const FzIntegerVariable* const var : ct->arguments[0].variables) { + for (const FzIntegerVariable* const var : ct->Arg(0).variables) { if (!var->domain.Contains(0)) { return false; } } FZVLOG << "Propagate " << ct->DebugString() << FZENDL; - for (FzIntegerVariable* const var : ct->arguments[0].variables) { + for (FzIntegerVariable* const var : ct->Arg(0).variables) { var->domain.IntersectWithInterval(0, 0); } ct->presolve_propagation_done = true; @@ -252,13 +252,13 @@ bool FzPresolver::PresolveArrayBoolOr(FzConstraint* ct) { bool FzPresolver::PresolveArrayBoolAnd(FzConstraint* ct) { if (!ct->presolve_propagation_done && ct->Arg(1).HasOneValue() && ct->Arg(1).Value() == 1) { - for (const FzIntegerVariable* const var : ct->arguments[0].variables) { + for (const FzIntegerVariable* const var : ct->Arg(0).variables) { if (!var->domain.Contains(1)) { return false; } } FZVLOG << "Propagate " << ct->DebugString() << FZENDL; - for (FzIntegerVariable* const var : ct->arguments[0].variables) { + for (FzIntegerVariable* const var : ct->Arg(0).variables) { var->domain.IntersectWithInterval(1, 1); } ct->presolve_propagation_done = true; @@ -319,16 +319,16 @@ bool FzPresolver::PresolveArrayIntElement(FzConstraint* ct) { // Reverse the constraint to get the original >=. This is true for ==, !=, <=, // <, >, >= and their reified versions. bool FzPresolver::PresolveLinear(FzConstraint* ct) { - if (ct->arguments[0].values.empty()) { + if (ct->Arg(0).values.empty()) { return false; } - for (const int64 coef : ct->arguments[0].values) { + for (const int64 coef : ct->Arg(0).values) { if (coef > 0) { return false; } } FZVLOG << "Reverse " << ct->DebugString() << FZENDL; - for (int64& coef : ct->arguments[0].values) { + for (int64& coef : ct->MutableArg(0)->values) { coef *= -1; } ct->MutableArg(2)->values[0] *= -1; @@ -356,7 +356,7 @@ bool FzPresolver::PresolvePropagatePositiveLinear(FzConstraint* ct) { if (ct->presolve_propagation_done || rhs < 0) { return false; } - for (const int64 coef : ct->arguments[0].values) { + for (const int64 coef : ct->Arg(0).values) { if (coef < 0) { return false; } @@ -367,8 +367,8 @@ bool FzPresolver::PresolvePropagatePositiveLinear(FzConstraint* ct) { return false; } } - for (int i = 0; i < ct->arguments[0].values.size(); ++i) { - const int64 coef = ct->arguments[0].values[i]; + for (int i = 0; i < ct->Arg(0).values.size(); ++i) { + const int64 coef = ct->Arg(0).values[i]; if (coef > 0) { FzIntegerVariable* const var = ct->Arg(1).variables[i]; var->domain.IntersectWithInterval(0, rhs / coef); @@ -382,17 +382,29 @@ bool FzPresolver::PresolvePropagatePositiveLinear(FzConstraint* ct) { // constraint with an affine mapping between y, z and the new index. // This rule stores the mapping to reconstruct the 2d element constraint. bool FzPresolver::PresolveStoreMapping(FzConstraint* ct) { - if (ct->arguments[0].values.size() == 2 && + if (ct->Arg(0).values.size() == 2 && ct->Arg(1).variables[0] == ct->target_variable && - ct->arguments[0].values[0] == -1 && - !ContainsKey(affine_map_, ct->Arg(1).variables[0]) && + ct->Arg(0).values[0] == -1 && + !ContainsKey(affine_map_, ct->target_variable) && ct->strong_propagation) { - affine_map_[ct->Arg(1).variables[0]] = - AffineMapping(ct->Arg(1).variables[1], ct->arguments[0].values[1], + affine_map_[ct->target_variable] = + AffineMapping(ct->Arg(1).variables[1], ct->Arg(0).values[1], -ct->Arg(2).Value(), ct); FZVLOG << "Store affine mapping info for " << ct->DebugString() << FZENDL; return true; } + if (ct->Arg(0).values.size() == 3 && + ct->Arg(1).variables[0] && ct->target_variable && + ct->Arg(0).values[0] == -1 && + ct->Arg(0).values[2] == 1 && + !ContainsKey(flatten_map_, ct->target_variable) && + ct->strong_propagation) { + flatten_map_[ct->target_variable] = + FlatteningMapping(ct->Arg(1).variables[1], ct->Arg(0).values[1], + ct->Arg(1).variables[2], -ct->Arg(2).Value(), ct); + FZVLOG << "Store affine mapping info for " << ct->DebugString() << FZENDL; + return true; + } return false; } @@ -422,7 +434,7 @@ bool FzPresolver::PresolveSimplifyElement(FzConstraint* ct) { } // Rewrite constraint. FZVLOG << "Simplify " << ct->DebugString() << FZENDL; - ct->arguments[0].variables[0] = mapping.variable; + ct->MutableArg(0)->variables[0] = mapping.variable; // TODO(user): Encapsulate argument setters. ct->MutableArg(1)->values.swap(new_values); if (ct->Arg(1).values.size() == 1) { @@ -471,7 +483,7 @@ bool FzPresolver::PresolveOneConstraint(FzConstraint* ct) { ct->Arg(1).HasOneValue() && ct->Arg(1).Value() == 0 && ContainsKey(abs_map_, ct->Arg(0).Var())) { FZVLOG << "Remove abs() from " << ct->DebugString() << FZENDL; - ct->arguments[0].variables[0] = abs_map_[ct->Arg(0).Var()]; + ct->MutableArg(0)->variables[0] = abs_map_[ct->Arg(0).Var()]; changed = true; } if (id == "int_eq") changed |= PresolveIntEq(ct); @@ -666,7 +678,7 @@ void Regroup(FzConstraint* start, const std::vector& chain, // End of chain, reconstruct. FzIntegerVariable* const out = carry_over.back(); start->arguments.pop_back(); - start->arguments[0].variables[0] = out; + start->MutableArg(0)->variables[0] = out; start->MutableArg(1)->type = FzArgument::INT_VAR_REF_ARRAY; start->MutableArg(1)->variables = chain; const std::string old_type = start->type; @@ -701,12 +713,18 @@ void FzPresolver::CleanUpModelForTheCpSolver(FzModel* model) { for (FzConstraint* const ct : model->constraints()) { const std::string& id = ct->type; // Remove ignored annotations on int_lin_eq. - if (id == "int_lin_eq" && ct->arguments[0].values.size() > 2 && + if (id == "int_lin_eq" && ct->Arg(0).values.size() > 3 && ct->strong_propagation) { FZVLOG << "Remove strong_propagation from " << ct->DebugString() << FZENDL; ct->strong_propagation = false; } + if (id == "int_lin_eq" && ct->Arg(0).values.size() >= 3 && + ct->target_variable != nullptr) { + FZVLOG << "Remove target_variable from " << ct->DebugString() + << FZENDL; + ct->RemoveTargetVariable(); + } // Remove target variables from constraints passed to SAT. if (ct->target_variable != nullptr && (id == "array_bool_and" || id == "array_bool_or" || diff --git a/src/flatzinc2/presolve.h b/src/flatzinc2/presolve.h index 01ae43bc7a..d37754308e 100644 --- a/src/flatzinc2/presolve.h +++ b/src/flatzinc2/presolve.h @@ -36,6 +36,24 @@ struct AffineMapping { : variable(v), coefficient(c), offset(o), constraint(ct) {} }; +// This struct stores the flattening mapping of two variables onto +// one: it represents var1 * coefficient + var2 + offset. It also +// stores the constraint that defines this mapping. +struct FlatteningMapping { + FzIntegerVariable* variable1; + int64 coefficient; + FzIntegerVariable* variable2; + int64 offset; + FzConstraint* constraint; + + FlatteningMapping() + : variable1(nullptr), coefficient(0), variable2(nullptr), offset(0), constraint(nullptr) {} + FlatteningMapping(FzIntegerVariable* v1, int64 c, FzIntegerVariable* v2, + int64 o, FzConstraint* ct) + : variable1(v1), coefficient(c), variable2(v2), offset(o), + constraint(ct) {} +}; + // The FzPresolver "pre-solves" a FzModel by applying some iterative // transformations to it, which may simplify and/or reduce the model. class FzPresolver { @@ -105,6 +123,9 @@ class FzPresolver { // Stores affine_map_[x] = a * y + b. hash_map affine_map_; + + // Stores flatten_map_[z] = a * x + y + b. + hash_map flatten_map_; }; } // namespace operations_research