always-cautious/ortools-clone
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

improve primal integral computation

Browse Source
This commit is contained in:
Laurent Perron
2021-01-08 09:49:56 +01:00
parent 165c9fa938
commit 5964c5e2dd
6 changed files with 66 additions and 15 deletions
Download Patch File Download Diff File Expand all files Collapse all files

2
ortools/sat/cp_model.proto
View File

@@ -558,7 +558,7 @@ enum CpSolverStatus {
// detailed error by calling ValidateCpModel(model_proto).
MODEL_INVALID = 1;
// A feasible solution as been found. But the search was stopped before we
// A feasible solution has been found. But the search was stopped before we
// could prove optimality or before we enumerated all solutions of a
// feasibility problem (if asked).
FEASIBLE = 2;

2
ortools/sat/cp_model_presolve.cc
View File

@@ -1706,7 +1706,7 @@ bool CpModelPresolver::PropagateDomainsInLinear(int c, ConstraintProto* ct) {
// Can we perform some substitution?
//
// TODO(user): there is no guarante we will not miss some since we might
// TODO(user): there is no guarantee we will not miss some since we might
// not reprocess a constraint once other have been deleted.
// Skip affine constraint. It is more efficient to substitute them lazily

24
ortools/sat/cp_model_solver.cc
View File

@@ -2878,6 +2878,9 @@ CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) {
*model->GetOrCreate<SatParameters>());
SharedTimeLimit shared_time_limit(model->GetOrCreate<TimeLimit>());
#if defined(_MSC_VER)
// On windows, The final_response is optimized out in the return part, and is
// swapped out before the cleanup callback is called. A workaround is to
// create a unique ptr that will forbid this optimization.
std::unique_ptr<CpSolverResponse> final_response_ptr(new CpSolverResponse());
CpSolverResponse& final_response = *final_response_ptr.get();
#else
@@ -3093,8 +3096,7 @@ CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) {
const auto& observers = model->GetOrCreate<SolutionObservers>()->observers;
if (!observers.empty()) {
shared_response_manager.AddSolutionCallback(
[&model_proto, &observers, &wall_timer, &user_timer,
&postprocess_solution, &shared_time_limit](
[&model_proto, &observers, &postprocess_solution](
const CpSolverResponse& response_of_presolved_problem) {
CpSolverResponse response = response_of_presolved_problem;
postprocess_solution(&response);
@@ -3113,6 +3115,23 @@ CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) {
});
}
// If specified, we load the initial objective domain right away in the
// response manager. Note that the presolve will always fill it with the
// trivial min/max value if the user left it empty. This avoids to display
// [-infinity, infinity] for the initial objective search space.
if (new_cp_model_proto.has_objective()) {
const Domain domain = ReadDomainFromProto(new_cp_model_proto.objective());
if (!domain.IsEmpty()) {
shared_response_manager.UpdateInnerObjectiveBounds(
"initial domain", IntegerValue(domain.Min()),
IntegerValue(domain.Max()));
}
}
// Start counting the primal integral from the current determistic time and
// initial objective domain gap that we just filled.
shared_response_manager.UpdatePrimalIntegral();
#if !defined(__PORTABLE_PLATFORM__)
if (absl::GetFlag(FLAGS_cp_model_dump_models)) {
const std::string presolved_file = absl::StrCat(
@@ -3170,6 +3189,7 @@ CpSolverResponse SolveCpModel(const CpModelProto& model_proto, Model* model) {
LOG(INFO) << absl::StrFormat("*** starting to load the model at %.2fs",
wall_timer.Get());
}
shared_response_manager.SetUpdatePrimalIntegralOnEachChange(true);
LoadCpModel(new_cp_model_proto, &shared_response_manager, model);
shared_response_manager.LoadDebugSolution(model);
if (log_search) {

15
ortools/sat/samples/minimal_jobshop_sat.py
View File

@@ -60,8 +60,9 @@ def MinimalJobshopSat():
end_var = model.NewIntVar(0, horizon, 'end' + suffix)
interval_var = model.NewIntervalVar(start_var, duration, end_var,
'interval' + suffix)
all_tasks[job_id, task_id] = task_type(
start=start_var, end=end_var, interval=interval_var)
all_tasks[job_id, task_id] = task_type(start=start_var,
end=end_var,
interval=interval_var)
machine_to_intervals[machine].append(interval_var)
# [END variables]
@@ -101,11 +102,11 @@ def MinimalJobshopSat():
for task_id, task in enumerate(job):
machine = task[0]
assigned_jobs[machine].append(
assigned_task_type(
start=solver.Value(all_tasks[job_id, task_id].start),
job=job_id,
index=task_id,
duration=task[1]))
assigned_task_type(start=solver.Value(
all_tasks[job_id, task_id].start),
job=job_id,
index=task_id,
duration=task[1]))
# Create per machine output lines.
output = ''

27
ortools/sat/synchronization.cc
View File

@@ -132,13 +132,21 @@ void LogNewSatSolution(const std::string& event_or_solution_count,
} // namespace
void SharedResponseManager::SetUpdatePrimalIntegralOnEachChange(bool set) {
absl::MutexLock mutex_lock(&mutex_);
update_integral_on_each_change_ = set;
}
void SharedResponseManager::UpdatePrimalIntegral() {
absl::MutexLock mutex_lock(&mutex_);
UpdatePrimalIntegralInternal();
}
void SharedResponseManager::UpdatePrimalIntegralInternal() {
if (!model_proto_.has_objective()) return;
const double current_time = shared_time_limit_->GetElapsedDeterministicTime();
const double time_delta = current_time - last_primal_integral_time_stamp_;
last_primal_integral_time_stamp_ = current_time;
// We use the log of the absolute objective gap.
//
@@ -148,10 +156,14 @@ void SharedResponseManager::UpdatePrimalIntegral() {
const CpObjectiveProto& obj = model_proto_.objective();
const double factor =
obj.scaling_factor() != 0.0 ? std::abs(obj.scaling_factor()) : 1.0;
const double bounds_delta = std::log(
1 + factor * std::abs(static_cast<double>(inner_objective_upper_bound_) -
static_cast<double>(inner_objective_lower_bound_)));
const double bounds_delta = std::log(1 + factor * last_absolute_gap_);
primal_integral_ += time_delta * bounds_delta;
// Update with new value.
last_primal_integral_time_stamp_ = current_time;
last_absolute_gap_ =
std::max(0.0, static_cast<double>(inner_objective_upper_bound_) -
static_cast<double>(inner_objective_lower_bound_));
}
void SharedResponseManager::SetGapLimitsFromParameters(
@@ -163,6 +175,11 @@ void SharedResponseManager::SetGapLimitsFromParameters(
}
void SharedResponseManager::TestGapLimitsIfNeeded() {
// This is called on each internal limit change, so it is a good place to
// update the integral. Note that this is not called at the end of the search
// though.
if (update_integral_on_each_change_) UpdatePrimalIntegralInternal();
if (absolute_gap_limit_ == 0 && relative_gap_limit_ == 0) return;
if (best_solution_objective_value_ >= kMaxIntegerValue) return;
if (inner_objective_lower_bound_ <= kMinIntegerValue) return;
@@ -228,6 +245,7 @@ void SharedResponseManager::UpdateInnerObjectiveBounds(
} else {
best_response_.set_status(CpSolverStatus::INFEASIBLE);
}
if (update_integral_on_each_change_) UpdatePrimalIntegralInternal();
if (log_updates_) LogNewSatSolution("Done", wall_timer_.Get(), worker_info);
return;
}
@@ -264,6 +282,7 @@ void SharedResponseManager::NotifyThatImprovingProblemIsInfeasible(
// We just proved that the best solution cannot be improved uppon, so we
// have a new lower bound.
inner_objective_lower_bound_ = best_solution_objective_value_;
if (update_integral_on_each_change_) UpdatePrimalIntegralInternal();
} else {
CHECK_EQ(num_solutions_, 0);
best_response_.set_status(CpSolverStatus::INFEASIBLE);

11
ortools/sat/synchronization.h
View File

@@ -215,6 +215,9 @@ class SharedResponseManager {
// particular instance. Or to evaluate how efficient our LNS code is improving
// solution.
//
// Note: The integral will start counting on the first UpdatePrimalIntegral()
// call, since before the difference is assumed to be zero.
//
// Important: To report a proper deterministic integral, we only update it
// on UpdatePrimalIntegral() which should be called in the main subsolver
// synchronization loop.
@@ -225,6 +228,11 @@ class SharedResponseManager {
double PrimalIntegral() const;
void UpdatePrimalIntegral();
// Sets this to true to have the "real" but non-deterministic primal integral.
// If this is true, then there is no need to manually call
// UpdatePrimalIntegral() but it is not an issue to do so.
void SetUpdatePrimalIntegralOnEachChange(bool set);
// Updates the inner objective bounds.
void UpdateInnerObjectiveBounds(const std::string& worker_info,
IntegerValue lb, IntegerValue ub);
@@ -291,6 +299,7 @@ class SharedResponseManager {
ABSL_EXCLUSIVE_LOCKS_REQUIRED(mutex_);
void SetStatsFromModelInternal(Model* model)
ABSL_EXCLUSIVE_LOCKS_REQUIRED(mutex_);
void UpdatePrimalIntegralInternal() ABSL_EXCLUSIVE_LOCKS_REQUIRED(mutex_);
const bool log_updates_;
const bool enumerate_all_solutions_;
@@ -317,7 +326,9 @@ class SharedResponseManager {
IntegerValue synchronized_inner_objective_upper_bound_
ABSL_GUARDED_BY(mutex_) = IntegerValue(kint64max);
bool update_integral_on_each_change_ ABSL_GUARDED_BY(mutex_) = false;
double primal_integral_ ABSL_GUARDED_BY(mutex_) = 0.0;
double last_absolute_gap_ ABSL_GUARDED_BY(mutex_) = 0.0;
double last_primal_integral_time_stamp_ ABSL_GUARDED_BY(mutex_) = 0.0;
int next_callback_id_ ABSL_GUARDED_BY(mutex_) = 0;