From 425901d7d01aa052e4ebf5d436547c3c7b8a3be0 Mon Sep 17 00:00:00 2001 From: Laurent Perron Date: Tue, 5 Jan 2021 22:13:51 +0100 Subject: [PATCH] new course scheduling example; sync c++ examples with internal code --- examples/cpp/course_scheduling.cc | 843 +++++++++++++++++++++++++ examples/cpp/course_scheduling.h | 86 +++ examples/cpp/course_scheduling.proto | 191 ++++++ examples/cpp/course_scheduling_run.cc | 108 ++++ examples/cpp/cvrp_disjoint_tw.cc | 13 +- examples/cpp/cvrptw.cc | 15 +- examples/cpp/jobshop_sat.cc | 656 ++++++++++++++----- examples/cpp/multi_knapsack_sat.cc | 3 +- examples/cpp/sat_cnf_reader.h | 2 +- examples/cpp/sat_runner.cc | 36 +- examples/cpp/shift_minimization_sat.cc | 15 +- examples/cpp/sports_scheduling_sat.cc | 38 +- examples/cpp/weighted_tardiness_sat.cc | 6 +- 13 files changed, 1807 insertions(+), 205 deletions(-) create mode 100644 examples/cpp/course_scheduling.cc create mode 100644 examples/cpp/course_scheduling.h create mode 100644 examples/cpp/course_scheduling.proto create mode 100644 examples/cpp/course_scheduling_run.cc diff --git a/examples/cpp/course_scheduling.cc b/examples/cpp/course_scheduling.cc new file mode 100644 index 0000000000..1d0f0ec001 --- /dev/null +++ b/examples/cpp/course_scheduling.cc @@ -0,0 +1,843 @@ +// Copyright 2010-2018 Google LLC +// 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 "examples/cpp/course_scheduling.h" + +#include +#include + +#include "absl/container/flat_hash_set.h" +#include "absl/status/status.h" +#include "absl/strings/numbers.h" +#include "absl/strings/str_cat.h" +#include "absl/strings/str_format.h" +#include "absl/strings/str_split.h" +#include "examples/cpp/course_scheduling.pb.h" +#include "ortools/base/logging.h" +#include "ortools/base/mathutil.h" +#include "ortools/base/status_macros.h" +#include "ortools/linear_solver/linear_solver.h" + +namespace operations_research { + +absl::Status CourseSchedulingSolver::ValidateModelAndLoadClasses( + const CourseSchedulingModel& model) { + time_slot_count_ = model.days_count() * model.daily_time_slot_count(); + room_count_ = model.rooms_size(); + solve_for_rooms_ = room_count_ > 0; + // If there are no rooms given, we have to give room_count at least one in + // order for the loops creating the solver variables and constraints to work. + if (!solve_for_rooms_) { + room_count_ = 1; + } + + // Validate the information given for each course. + for (const Course& course : model.courses()) { + if (course.consecutive_slots_count() != 1 && + course.consecutive_slots_count() != 2) { + return absl::InvalidArgumentError(absl::StrFormat( + "The course titled %s has %d consecutive time slots specified when " + "it can only have 1 or 2.", + course.display_name(), course.consecutive_slots_count())); + } + + if (course.teacher_section_counts_size() != course.teacher_indices_size()) { + return absl::InvalidArgumentError( + absl::StrFormat("The course titled %s should have the same number of " + "teacher indices and section numbers.", + course.display_name())); + } + + for (const int teacher_index : course.teacher_indices()) { + if (teacher_index >= model.teachers_size()) { + return absl::InvalidArgumentError(absl::StrFormat( + "The course titled %s has teacher %d assigned to it but there are " + "only %d teachers.", + course.display_name(), teacher_index, model.teachers_size())); + } + } + + for (const int room_index : course.room_indices()) { + if (room_index >= model.rooms_size()) { + return absl::InvalidArgumentError( + absl::StrFormat("The course titled %s is slotted for room index %d " + "but there are only %d rooms.", + course.display_name(), room_index, room_count_)); + } + } + } + + // Validate the information given for each teacher and create hash sets of the + // restricted indices for each teacher. + teacher_to_restricted_slots_ = + std::vector>(model.teachers_size()); + for (int teacher_index = 0; teacher_index < model.teachers_size(); + ++teacher_index) { + for (const int restricted_slot : + model.teachers(teacher_index).restricted_time_slots()) { + if (restricted_slot >= time_slot_count_) { + return absl::InvalidArgumentError( + absl::StrFormat("Teacher with name %s has restricted time slot %d " + "but there are only %d time slots.", + model.teachers(teacher_index).display_name(), + restricted_slot, time_slot_count_)); + } + teacher_to_restricted_slots_[teacher_index].insert(restricted_slot); + } + } + + // Since each course can have multiple sections (classes), we need to + // "flatten" out each course so that each of its sections gets a unique index. + // This vector stores the information for calculating each unique index. The + // first value is the unique index the course's sections begin at. The second + // value is the total number of sections of that course. + course_to_classes_ = std::vector>(model.courses_size()); + // For each teacher, store the class unique indices that they teach. + teacher_to_classes_ = + std::vector>(model.teachers_size()); + // Store course indices of courses that have a single section. + absl::flat_hash_set singleton_courses; + int flattened_course_index = 0; + for (int course_index = 0; course_index < model.courses_size(); + ++course_index) { + const Course& course = model.courses(course_index); + int total_section_count = 0; + for (int teacher = 0; teacher < course.teacher_indices_size(); ++teacher) { + for (int section = 0; section < course.teacher_section_counts(teacher); + ++section) { + teacher_to_classes_[course.teacher_indices(teacher)].insert( + flattened_course_index); + course_to_classes_[course_index].push_back(flattened_course_index); + ++flattened_course_index; + } + total_section_count += course.teacher_section_counts(teacher); + } + if (total_section_count == 1) { + singleton_courses.insert(course_index); + } + } + class_count_ = flattened_course_index; + + // Validate the information given for each student. Store each student's + // course pairs. + for (const Student& student : model.students()) { + for (const int course_index : student.course_indices()) { + if (course_index >= model.courses_size()) { + return absl::InvalidArgumentError(absl::StrFormat( + "Student with name %s has course index %d but there are only %d " + "courses.", + student.display_name(), course_index, model.courses_size())); + } + } + InsertSortedPairs(std::vector(student.course_indices().begin(), + student.course_indices().end()), + &course_conflicts_); + } + + LOG(INFO) << "Number of days: " << model.days_count(); + LOG(INFO) << "Number of daily time slots: " << model.daily_time_slot_count(); + LOG(INFO) << "Total number of time slots: " << time_slot_count_; + LOG(INFO) << "Number of courses: " << model.courses_size(); + LOG(INFO) << "Total number of classes: " << class_count_; + LOG(INFO) << "Number of teachers: " << model.teachers_size(); + LOG(INFO) << "Number of students: " << model.students_size(); + if (solve_for_rooms_) { + LOG(INFO) << "Number of rooms: " << model.rooms_size(); + } + + return absl::OkStatus(); +} + +CourseSchedulingResult CourseSchedulingSolver::Solve( + const CourseSchedulingModel& model) { + CourseSchedulingResult result; + + const auto validation_status = ValidateModelAndLoadClasses(model); + if (!validation_status.ok()) { + result.set_solver_status( + CourseSchedulingResultStatus::SOLVER_MODEL_INVALID); + result.set_message(validation_status.message()); + return result; + } + + ConflictPairs class_conflicts; + result = SolveModel(model, class_conflicts); + + if (result.solver_status() != SOLVER_FEASIBLE && + result.solver_status() != SOLVER_OPTIMAL) { + return result; + } + + const auto verifier_status = VerifyCourseSchedulingResult(model, result); + if (!verifier_status.ok()) { + result.set_solver_status(CourseSchedulingResultStatus::ABNORMAL); + result.set_message(verifier_status.message()); + } + + return result; +} + +CourseSchedulingResult CourseSchedulingSolver::SolveModel( + const CourseSchedulingModel& model, const ConflictPairs& class_conflicts) { + CourseSchedulingResult result; + + result = ScheduleCourses(class_conflicts, model); + if (result.solver_status() != SOLVER_FEASIBLE && + result.solver_status() != SOLVER_OPTIMAL) { + if (result.solver_status() == SOLVER_INFEASIBLE) { + result.set_message("The problem is infeasible with the given courses."); + } + return result; + } + + ConflictPairs class_conflicts_to_try = AssignStudents(model, &result); + + if (class_conflicts_to_try.empty()) return result; + + std::vector> conflicts(class_conflicts_to_try.begin(), + class_conflicts_to_try.end()); + + const int conflicts_count = conflicts.size(); + const int conflicts_log = conflicts_count == 1 ? 1 : log2(conflicts_count); + for (int i = 0; i < conflicts_log; ++i) { + const int divisions = MathUtil::IPow(2, i); + for (int j = 0; j < divisions; ++j) { + const int start = std::floor(conflicts_count * j / divisions); + const int end = std::floor(conflicts_count * (j + 1) / divisions); + + ConflictPairs new_class_conflicts = class_conflicts; + if (end >= conflicts_count) { + new_class_conflicts.insert(conflicts.begin() + start, conflicts.end()); + } else { + new_class_conflicts.insert(conflicts.begin() + start, + conflicts.begin() + end); + } + + result = SolveModel(model, new_class_conflicts); + if (result.solver_status() == SOLVER_FEASIBLE || + result.solver_status() == SOLVER_OPTIMAL) { + return result; + } + } + } + + return result; +} + +std::vector CourseSchedulingSolver::GetRoomIndices(const Course& course) { + if (solve_for_rooms_) { + return std::vector(course.room_indices().begin(), + course.room_indices().end()); + } + + return {0}; +} + +void CourseSchedulingSolver::InsertSortedPairs(const std::vector& list, + ConflictPairs* pairs) { + for (int first = 1; first < list.size(); ++first) { + for (int second = first; second < list.size(); ++second) { + pairs->insert(std::minmax(list[first - 1], list[second])); + } + } +} + +std::vector> +CourseSchedulingSolver::GetClassesByTimeSlot( + const CourseSchedulingResult* result) { + std::vector> time_slot_to_classes(time_slot_count_); + + for (const ClassAssignment& class_assignment : result->class_assignments()) { + const int course_index = class_assignment.course_index(); + const int section_number = class_assignment.section_number(); + for (const int time_slot : class_assignment.time_slots()) { + time_slot_to_classes[time_slot].insert( + course_to_classes_[course_index][section_number]); + } + } + + return time_slot_to_classes; +} + +void CourseSchedulingSolver::AddVariableIfNonNull(double coeff, + const MPVariable* var, + MPConstraint* ct) { + if (var == nullptr) return; + + ct->SetCoefficient(var, coeff); +} + +CourseSchedulingResult CourseSchedulingSolver::ScheduleCourses( + const ConflictPairs& class_conflicts, const CourseSchedulingModel& model) { + LOG(INFO) << "Starting schedule courses solver with " + << class_conflicts.size() << " class conflicts."; + MPSolver mip_solver("CourseSchedulingMIP", + MPSolver::SCIP_MIXED_INTEGER_PROGRAMMING); + const double kInfinity = std::numeric_limits::infinity(); + + // Create binary variables x(n,t,r) where x(n,t,r) = 1 indicates that course n + // is scheduled for time slot t in course r. Variables are only created if + // attendees of class n are available for time slot t and if the course can be + // placed into room r. + std::vector>> variables( + class_count_, + std::vector>( + time_slot_count_, std::vector(room_count_, nullptr))); + for (int proto_index = 0; proto_index < model.courses_size(); ++proto_index) { + const Course& course = model.courses(proto_index); + const int course_index = course_to_classes_[proto_index][0]; + int total_section_index = 0; + for (int i = 0; i < course.teacher_section_counts_size(); ++i) { + const int teacher = course.teacher_indices(i); + const int section_count = course.teacher_section_counts(i); + for (int section = 0; section < section_count; ++section) { + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + if (teacher_to_restricted_slots_[teacher].contains(time_slot)) { + continue; + } + + for (const int room : GetRoomIndices(course)) { + variables[course_index + total_section_index][time_slot][room] = + mip_solver.MakeBoolVar(absl::StrFormat( + "x_%d_%d_%d", course_index + total_section_index, time_slot, + room)); + } + } + ++total_section_index; + } + } + } + + std::vector> intermediate_vars( + class_count_, std::vector(time_slot_count_)); + for (int class_index = 0; class_index < class_count_; ++class_index) { + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, 0); + for (int room = 0; room < room_count_; ++room) { + if (variables[class_index][time_slot][room] == nullptr) continue; + + ct->SetCoefficient(variables[class_index][time_slot][room], 1); + } + if (!ct->terms().empty()) { + intermediate_vars[class_index][time_slot] = mip_solver.MakeBoolVar( + absl::StrFormat("intermediate_%d_%d", class_index, time_slot)); + ct->SetCoefficient(intermediate_vars[class_index][time_slot], -1); + } + } + } + + // 1. Each course meets no more than its number of consecutive time slots a + // day. + for (int day = 0; day < model.days_count(); ++day) { + for (int course = 0; course < model.courses_size(); ++course) { + const int consecutive_slot_count = + model.courses(course).consecutive_slots_count(); + for (const int class_index : course_to_classes_[course]) { + MPConstraint* const ct = + mip_solver.MakeRowConstraint(0, consecutive_slot_count); + for (int daily_time_slot = 0; + daily_time_slot < model.daily_time_slot_count(); + ++daily_time_slot) { + AddVariableIfNonNull( + 1, + intermediate_vars[class_index] + [(day * model.daily_time_slot_count()) + + daily_time_slot], + ct); + } + } + } + } + + // 2. Each course must meet the given number of times * its number of + // consecutive time slots. + for (int course = 0; course < model.courses_size(); ++course) { + const int meeting_count = model.courses(course).meetings_count(); + const int consecutive_slot_count = + model.courses(course).consecutive_slots_count(); + for (const int class_index : course_to_classes_[course]) { + MPConstraint* const ct = + mip_solver.MakeRowConstraint(meeting_count * consecutive_slot_count, + meeting_count * consecutive_slot_count); + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + AddVariableIfNonNull(1, intermediate_vars[class_index][time_slot], ct); + } + } + } + + // 3. Teachers are scheduled for no more than one course per time slot. + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + for (int teacher = 0; teacher < model.teachers_size(); ++teacher) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, 1); + for (const int class_index : teacher_to_classes_[teacher]) { + AddVariableIfNonNull(1, intermediate_vars[class_index][time_slot], ct); + } + } + } + + // 4. Each room can only be occupied by one course for each time slot. + if (solve_for_rooms_) { + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + for (int room = 0; room < room_count_; ++room) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, 1); + for (int class_index = 0; class_index < class_count_; ++class_index) { + AddVariableIfNonNull(1, variables[class_index][time_slot][room], ct); + } + } + } + } + + // 5. Ensure each class is scheduled for the correct number of consecutive + // time slots. + for (int course = 0; course < model.courses_size(); ++course) { + const int consecutive_slot_count = + model.courses(course).consecutive_slots_count(); + if (consecutive_slot_count == 1) continue; + for (const int class_index : course_to_classes_[course]) { + for (int day = 0; day < model.days_count(); ++day) { + for (int room = 0; room < room_count_; ++room) { + // If only the previous time slot is chosen, force the current time + // slot to be chosen. + for (int daily_time_slot = 0; + daily_time_slot < model.daily_time_slot_count(); + ++daily_time_slot) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, kInfinity); + const int time_slot_offset = + day * model.daily_time_slot_count() + daily_time_slot; + + if (daily_time_slot > 0) { + AddVariableIfNonNull( + 1, variables[class_index][time_slot_offset - 1][room], ct); + } + AddVariableIfNonNull( + -0.5, variables[class_index][time_slot_offset][room], ct); + if (daily_time_slot < model.daily_time_slot_count() - 1) { + AddVariableIfNonNull( + 0.5, variables[class_index][time_slot_offset + 1][room], ct); + } + } + } + } + } + } + + // 6. Ensure that there are at least two sections of each course_conflicts_ + // pair that are scheduled for different time slots. + for (const auto& conflict_pair : course_conflicts_) { + const int course_1 = conflict_pair.first; + const int course_2 = conflict_pair.second; + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + const int bound = course_to_classes_[course_1].size() + + course_to_classes_[course_2].size() - 1; + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, bound); + for (const int class_1 : course_to_classes_[course_1]) { + AddVariableIfNonNull(1, intermediate_vars[class_1][time_slot], ct); + } + for (const int class_2 : course_to_classes_[course_2]) { + AddVariableIfNonNull(1, intermediate_vars[class_2][time_slot], ct); + } + } + } + + // 7. Ensure that conflicting class pairs are not scheduled for the same time + // slot. + for (const auto& conflict_pair : class_conflicts) { + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, 1); + AddVariableIfNonNull(1, intermediate_vars[conflict_pair.first][time_slot], + ct); + AddVariableIfNonNull( + 1, intermediate_vars[conflict_pair.second][time_slot], ct); + } + } + + MPSolver::ResultStatus status = mip_solver.Solve(); + + CourseSchedulingResult result; + result.set_solver_status(MipStatusToCourseSchedulingResultStatus(status)); + if (status != MPSolver::OPTIMAL && status != MPSolver::FEASIBLE) { + if (status == MPSolver::UNBOUNDED) { + result.set_message( + "MIP solver returned UNBOUNDED: the model is solved but the solution " + "is infinity"); + } else if (status == MPSolver::ABNORMAL) { + result.set_message( + "MIP solver returned ABNORMAL: some error occurred while solving"); + } + return result; + } + + for (int course = 0; course < model.courses_size(); ++course) { + for (int section = 0; section < course_to_classes_[course].size(); + ++section) { + ClassAssignment class_assignment; + class_assignment.set_course_index(course); + class_assignment.set_section_number(section); + const int class_index = course_to_classes_[course][section]; + + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + for (int room = 0; room < room_count_; ++room) { + MPVariable* x_i = variables[class_index][time_slot][room]; + if (x_i != nullptr && x_i->solution_value() == 1) { + if (solve_for_rooms_) { + class_assignment.add_room_indices(room); + } + class_assignment.add_time_slots(time_slot); + } + } + } + *result.add_class_assignments() = class_assignment; + } + } + return result; +} + +CourseSchedulingSolver::ConflictPairs CourseSchedulingSolver::AssignStudents( + const CourseSchedulingModel& model, CourseSchedulingResult* result) { + LOG(INFO) << "Starting assign students solver."; + MPSolver mip_solver("AssignStudentsMIP", + MPSolver::SCIP_MIXED_INTEGER_PROGRAMMING); + + // Create binary variables y(s,n) where y(s,n) = 1 indicates that student s is + // enrolled in class n. Variables are created for a student and each section + // of a course that they are signed up to take. + std::vector> variables( + model.students_size(), std::vector(class_count_, nullptr)); + for (int student_index = 0; student_index < model.students_size(); + ++student_index) { + const Student& student = model.students(student_index); + for (const int course_index : student.course_indices()) { + for (const int class_index : course_to_classes_[course_index]) { + variables[student_index][class_index] = mip_solver.MakeBoolVar( + absl::StrFormat("y_%d_%d", student_index, class_index)); + } + } + } + + // 1. Each student must be assigned to exactly one section for each course + // they are signed up to take. + for (int student_index = 0; student_index < model.students_size(); + ++student_index) { + const Student& student = model.students(student_index); + for (const int course_index : student.course_indices()) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(1, 1); + for (const int class_index : course_to_classes_[course_index]) { + AddVariableIfNonNull(1, variables[student_index][class_index], ct); + } + } + } + + // 2. Ensure that the minimum and maximum capacities for each class are met. + for (int course_index = 0; course_index < model.courses_size(); + ++course_index) { + const Course& course = model.courses(course_index); + const int min_cap = course.min_capacity(); + const int max_cap = course.max_capacity(); + for (const int class_index : course_to_classes_[course_index]) { + MPConstraint* const ct = mip_solver.MakeRowConstraint(min_cap, max_cap); + for (int student = 0; student < model.students_size(); ++student) { + AddVariableIfNonNull(1, variables[student][class_index], ct); + } + } + } + + // 3. Each student should be assigned to one class per time slot. This is a + // soft constraint -- if violated, then the variable infeasibility_var(s,t) + // will be greater than 0 for that student s and time slot t. + std::vector> infeasibility_vars( + model.students_size(), + std::vector(time_slot_count_, nullptr)); + std::vector> time_slot_to_classes = + GetClassesByTimeSlot(result); + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + for (int student_index = 0; student_index < model.students_size(); + ++student_index) { + infeasibility_vars[student_index][time_slot] = mip_solver.MakeIntVar( + 0, class_count_, + absl::StrFormat("f_%d_%d", student_index, time_slot)); + const Student& student = model.students(student_index); + + MPConstraint* const ct = mip_solver.MakeRowConstraint(0, 1); + ct->SetCoefficient(infeasibility_vars[student_index][time_slot], -1); + for (const int course_index : student.course_indices()) { + for (const int class_index : course_to_classes_[course_index]) { + if (!time_slot_to_classes[time_slot].contains(class_index)) continue; + + AddVariableIfNonNull(1, variables[student_index][class_index], ct); + } + } + } + } + + // Minimize the infeasibility vars. If the objective is 0, then we have found + // a feasible solution for the course schedule. If the objective is greater + // than 0, then some students were assigned to multiple courses for the same + // time slot and we need to find a new schedule for the courses. + MPObjective* const objective = mip_solver.MutableObjective(); + for (int student_index = 0; student_index < model.students_size(); + ++student_index) { + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + objective->SetCoefficient(infeasibility_vars[student_index][time_slot], + 1); + } + } + + mip_solver.SetSolverSpecificParametersAsString("limits/gap=0.01"); + MPSolver::ResultStatus status = mip_solver.Solve(); + ConflictPairs class_conflict_pairs; + + // This model will only be infeasible if the minimum or maximum capacities + // are violated. + if (status != MPSolver::OPTIMAL && status != MPSolver::FEASIBLE) { + result->set_solver_status(MipStatusToCourseSchedulingResultStatus(status)); + result->clear_class_assignments(); + if (status == MPSolver::INFEASIBLE) { + result->set_message( + "Check the minimum or maximum capacity constraints for your " + "classes."); + } + + return class_conflict_pairs; + } + + LOG(INFO) << "Finished assign students solver with " << objective->Value() + << " schedule violations."; + if (objective->Value() > 0) { + for (int time_slot = 0; time_slot < time_slot_count_; ++time_slot) { + for (int student_index = 0; student_index < model.students_size(); + ++student_index) { + std::vector conflicting_classes; + MPVariable* const f_i = infeasibility_vars[student_index][time_slot]; + if (f_i != nullptr && f_i->solution_value() == 0) continue; + + for (const int class_index : time_slot_to_classes[time_slot]) { + MPVariable* const y_i = variables[student_index][class_index]; + if (y_i != nullptr && y_i->solution_value() == 1) { + conflicting_classes.push_back(class_index); + } + } + InsertSortedPairs(conflicting_classes, &class_conflict_pairs); + } + } + return class_conflict_pairs; + } + + for (int student_index = 0; student_index < model.students_size(); + ++student_index) { + StudentAssignment student_assignment; + student_assignment.set_student_index(student_index); + + const Student& student = model.students(student_index); + for (const int course_index : student.course_indices()) { + for (int section_index = 0; + section_index < course_to_classes_[course_index].size(); + ++section_index) { + int class_index = course_to_classes_[course_index][section_index]; + MPVariable* const y_i = variables[student_index][class_index]; + + if (y_i->solution_value() == 1) { + student_assignment.add_course_indices(course_index); + student_assignment.add_section_indices(section_index); + } + } + } + *result->add_student_assignments() = student_assignment; + } + + return class_conflict_pairs; +} + +CourseSchedulingResultStatus +CourseSchedulingSolver::MipStatusToCourseSchedulingResultStatus( + MPSolver::ResultStatus mip_status) { + switch (mip_status) { + case MPSolver::ResultStatus::OPTIMAL: + return SOLVER_OPTIMAL; + case MPSolver::ResultStatus::FEASIBLE: + return SOLVER_FEASIBLE; + case MPSolver::ResultStatus::INFEASIBLE: + return SOLVER_INFEASIBLE; + case MPSolver::ResultStatus::UNBOUNDED: + case MPSolver::ResultStatus::MODEL_INVALID: + return SOLVER_MODEL_INVALID; + case MPSolver::ResultStatus::NOT_SOLVED: + return SOLVER_NOT_SOLVED; + case MPSolver::ResultStatus::ABNORMAL: + return ABNORMAL; + default: + return COURSE_SCHEDULING_RESULT_STATUS_UNSPECIFIED; + } +} + +absl::Status CourseSchedulingSolver::VerifyCourseSchedulingResult( + const CourseSchedulingModel& model, const CourseSchedulingResult& result) { + std::vector> class_to_time_slots(class_count_); + std::vector> room_to_time_slots(model.rooms_size()); + for (const ClassAssignment& class_assignment : result.class_assignments()) { + const int course_index = class_assignment.course_index(); + const int meetings_count = model.courses(course_index).meetings_count(); + const int consecutive_time_slots = + model.courses(course_index).consecutive_slots_count(); + + // Verify that each class meets the correct number of times. + if (class_assignment.time_slots_size() != + meetings_count * consecutive_time_slots) { + return absl::InvalidArgumentError(absl::StrFormat( + "Verification failed: The course titled %s and section number %d " + "meets %d times when it should meet %d times.", + model.courses(course_index).display_name(), + class_assignment.section_number(), class_assignment.time_slots_size(), + meetings_count * consecutive_time_slots)); + } + + const int class_index = + course_to_classes_[course_index][class_assignment.section_number()]; + std::vector> day_to_time_slots(model.days_count()); + for (const int time_slot : class_assignment.time_slots()) { + class_to_time_slots[class_index].insert(time_slot); + day_to_time_slots[std::floor(time_slot / model.daily_time_slot_count())] + .push_back(time_slot); + } + + // Verify that a class meets no more than its consecutive time slot count + // per day. If a class needs 2 consecutive time slots, verify that it is + // scheduled accordingly. + for (int day = 0; day < model.days_count(); ++day) { + const int day_count = day_to_time_slots[day].size(); + if (day_count != 0 && day_count != consecutive_time_slots) { + return absl::InvalidArgumentError( + absl::StrFormat("Verification failed: The course titled %s does " + "not meet the correct number of times in " + "day %d.", + model.courses(course_index).display_name(), day)); + } + if (day_count != 2) continue; + + const int first_slot = day_to_time_slots[day][0]; + const int second_slot = day_to_time_slots[day][1]; + if (std::abs(first_slot - second_slot) != 1) { + return absl::InvalidArgumentError( + absl::StrFormat("Verification failed: The course titled %s is not " + "scheduled for consecutive time slots " + "in day %d.", + model.courses(course_index).display_name(), day)); + } + } + + // Verify that their is no more than 1 class per room for each time slot. + if (solve_for_rooms_) { + for (int i = 0; i < class_assignment.room_indices_size(); ++i) { + const int room = class_assignment.room_indices(i); + const int time_slot = class_assignment.time_slots(i); + if (room_to_time_slots[room].contains(time_slot)) { + return absl::InvalidArgumentError( + absl::StrFormat("Verification failed: Multiple classes have " + "been assigned to room %s during time slot %d.", + model.rooms(room).display_name(), time_slot)); + } + room_to_time_slots[room].insert(time_slot); + } + } + } + + // Verify that each teacher is assigned to no more than one class per time + // slot and that each teacher is not assigned to their restricted time slots. + for (int teacher = 0; teacher < model.teachers_size(); ++teacher) { + const auto& class_list = teacher_to_classes_[teacher]; + absl::flat_hash_set teacher_time_slots; + for (const int class_index : class_list) { + for (const int time_slot : class_to_time_slots[class_index]) { + if (teacher_to_restricted_slots_[teacher].contains(time_slot)) { + return absl::InvalidArgumentError(absl::StrFormat( + "Verification failed: Teacher with name %s has been assigned to " + "restricted time slot %d.", + model.teachers(teacher).display_name(), time_slot)); + } + if (teacher_time_slots.contains(time_slot)) { + return absl::InvalidArgumentError(absl::StrFormat( + "Verification failed: Teacher with name %s has been assigned to " + "multiple classes at time slot %d.", + model.teachers(teacher).display_name(), time_slot)); + } + teacher_time_slots.insert(time_slot); + } + } + } + + std::vector class_student_count(class_count_); + for (const StudentAssignment& student_assignment : + result.student_assignments()) { + const int student_index = student_assignment.student_index(); + + // Verify that each student is assigned to the correct courses. + std::vector enrolled_courses = + std::vector(model.students(student_index).course_indices().begin(), + model.students(student_index).course_indices().end()); + std::vector assigned_courses = + std::vector(student_assignment.course_indices().begin(), + student_assignment.course_indices().end()); + std::sort(enrolled_courses.begin(), enrolled_courses.end()); + std::sort(assigned_courses.begin(), assigned_courses.end()); + if (enrolled_courses != assigned_courses) { + return absl::InvalidArgumentError( + absl::StrFormat("Verification failed: Student with name %s has not " + "been assigned the correct courses.", + model.students(student_index).display_name())); + } + + // Verify that each student is assigned to no more than one class per time + // slot. + absl::flat_hash_set student_time_slots; + for (int i = 0; i < student_assignment.course_indices_size(); ++i) { + const int course_index = student_assignment.course_indices(i); + const int section = student_assignment.section_indices(i); + const int class_index = course_to_classes_[course_index][section]; + ++class_student_count[class_index]; + + for (const int time_slot : class_to_time_slots[class_index]) { + if (student_time_slots.contains(time_slot)) { + return absl::InvalidArgumentError(absl::StrFormat( + "Verification failed: Student with name %s has been assigned to " + "multiple classes at time slot %d.", + model.students(student_index).display_name(), time_slot)); + } + student_time_slots.insert(time_slot); + } + } + } + + // Verify size of each class is within the minimum and maximum capacities. + for (int course = 0; course < model.courses_size(); ++course) { + const int min_cap = model.courses(course).min_capacity(); + const int max_cap = model.courses(course).max_capacity(); + for (const int class_index : course_to_classes_[course]) { + const int class_size = class_student_count[class_index]; + if (class_size < min_cap) { + return absl::InvalidArgumentError(absl::StrFormat( + "Verification failed: The course titled %s has %d students when it " + "should have at least %d students.", + model.courses(course).display_name(), class_size, min_cap)); + } + if (class_size > max_cap) { + return absl::InvalidArgumentError(absl::StrFormat( + "Verification failed: The course titled %s has %d students when it " + "should have no more than %d students.", + model.courses(course).display_name(), class_size, max_cap)); + } + } + } + + return absl::OkStatus(); +} + +} // namespace operations_research diff --git a/examples/cpp/course_scheduling.h b/examples/cpp/course_scheduling.h new file mode 100644 index 0000000000..8a35e1e9c2 --- /dev/null +++ b/examples/cpp/course_scheduling.h @@ -0,0 +1,86 @@ +// Copyright 2010-2018 Google LLC +// 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_EXAMPLES_COURSE_SCHEDULING_H_ +#define OR_TOOLS_EXAMPLES_COURSE_SCHEDULING_H_ + +#include +#include + +#include "absl/container/flat_hash_set.h" +#include "absl/status/status.h" +#include "absl/strings/str_format.h" +#include "examples/cpp/course_scheduling.pb.h" +#include "ortools/linear_solver/linear_solver.h" +#include "ortools/sat/cp_model.pb.h" + +namespace operations_research { +class CourseSchedulingSolver { + public: + CourseSchedulingSolver() : solve_for_rooms_(false) {} + virtual ~CourseSchedulingSolver() {} + + using ConflictPairs = absl::flat_hash_set>; + + CourseSchedulingResult Solve(const CourseSchedulingModel& model); + + protected: + virtual absl::Status ValidateModelAndLoadClasses( + const CourseSchedulingModel& model); + + virtual CourseSchedulingResult SolveModel( + const CourseSchedulingModel& model, const ConflictPairs& class_conflicts); + + virtual absl::Status VerifyCourseSchedulingResult( + const CourseSchedulingModel& model, const CourseSchedulingResult& result); + + private: + CourseSchedulingResult ScheduleCourses(const ConflictPairs& class_conflicts, + const CourseSchedulingModel& model); + + // This method will modify the CourseSchedulingResult returned from + // ScheduleCoursesMip, which is why the result is passed in as a pointer. + ConflictPairs AssignStudents(const CourseSchedulingModel& model, + CourseSchedulingResult* result); + + int GetTeacherIndex(int course_index, int section); + + void InsertSortedPairs(const std::vector& list, ConflictPairs* pairs); + + bool ShouldCreateVariable(int course_index, int section, int time_slot, + int room); + + std::vector GetRoomIndices(const Course& course); + + std::vector> GetClassesByTimeSlot( + const CourseSchedulingResult* result); + + void AddVariableIfNonNull(double coeff, const MPVariable* var, + MPConstraint* ct); + + CourseSchedulingResultStatus MipStatusToCourseSchedulingResultStatus( + MPSolver::ResultStatus mip_status); + + bool solve_for_rooms_; + int class_count_; + int time_slot_count_; + int room_count_; + ConflictPairs course_conflicts_; + std::vector> teacher_to_classes_; + std::vector> teacher_to_restricted_slots_; + std::vector> course_to_classes_; +}; + +} // namespace operations_research + +#endif // OR_TOOLS_EXAMPLES_COURSE_SCHEDULING_H_ diff --git a/examples/cpp/course_scheduling.proto b/examples/cpp/course_scheduling.proto new file mode 100644 index 0000000000..3ad0d688d4 --- /dev/null +++ b/examples/cpp/course_scheduling.proto @@ -0,0 +1,191 @@ +// Copyright 2010-2018 Google LLC +// 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. + +syntax = "proto3"; + +package operations_research; + +import "google/api/field_behavior.proto"; + +// Information required to create a schedule for a school system. +message CourseSchedulingModel { + // Schedule name, used only for logging purposes. + string display_name = 1; + + // The number of days in a schedule rotation. If the school system uses a + // block schedule, this value should be 1. + int32 days_count = 2 [(google.api.field_behavior) = REQUIRED]; + + // The number of time slots each day in a schedule rotation. If the school + // system uses a block schedule, this value is the number of blocks. + int32 daily_time_slot_count = 3 [(google.api.field_behavior) = REQUIRED]; + + // List of courses that need to be scheduled. + repeated Course courses = 4; + + // List of teachers. + repeated Teacher teachers = 5; + + // List of students that need to be assigned to these courses. + repeated Student students = 6; + + // List of rooms that the courses can be assigned to. + repeated Room rooms = 7; +} + +// Holds the solution to the course scheduling problem. +message CourseSchedulingResult { + // Human readable message about the solver or given model. + string message = 1; + + // Status of the solver. + CourseSchedulingResultStatus solver_status = 2 + [(google.api.field_behavior) = REQUIRED]; + + // List of the time slot and room assignments for each section of a course. + repeated ClassAssignment class_assignments = 3; + + // List of course and section assignments for each student. + repeated StudentAssignment student_assignments = 4; +} + +message ClassAssignment { + // Index of the course in the CourseSchedulingModel. + int32 course_index = 1 [(google.api.field_behavior) = REQUIRED]; + + // Specific section of the course in the CourseSchedulingModel. + int32 section_number = 2 [(google.api.field_behavior) = REQUIRED]; + + // Time slots that this class has been assigned to in the + // CourseSchedulingModel. + repeated int32 time_slots = 3 [(google.api.field_behavior) = REQUIRED]; + + // Indices of the rooms that the class is assigned to in the + // CourseSchedulingModel. If this is not empty, then the number of indices + // must match the number of time slots. + repeated int32 room_indices = 4; +} + +message StudentAssignment { + // Index of the student in the CourseSchedulingModel. + int32 student_index = 1 [(google.api.field_behavior) = REQUIRED]; + + // Course indices in the CourseSchedulingModel that this student has been + // assigned to. The number of indices must match the number of section + // indices. + repeated int32 course_indices = 2 [(google.api.field_behavior) = REQUIRED]; + + // Section indices for each Course in the CourseSchedulingModel this this + // student has been assigned to. The number of indices must match the number + // of course indices. + repeated int32 section_indices = 3 [(google.api.field_behavior) = REQUIRED]; +} + +message Course { + // Course name, used only for logging purposes. + string display_name = 1; + + // The number of times each section of this course needs to meet during a + // schedule rotation. Each section of the course meets no more than once a + // day. If the school system uses a block schedule, then this value should + // be 1. + int32 meetings_count = 2 [(google.api.field_behavior) = REQUIRED]; + + // The maximum number of students for this course. This value can be equal to + // +Infinity to encode a course has no maximum capacity. + int32 max_capacity = 3 [(google.api.field_behavior) = REQUIRED]; + + // The minimum number of students for this course. + int32 min_capacity = 4; + + // The number of consecutive time slots that each section of this course needs + // to be scheduled for. This value can only be 1 or 2. If the value is 2, then + // 2 consecutive time slots in a day counts as 1 meeting time for the section. + int32 consecutive_slots_count = 5 [(google.api.field_behavior) = REQUIRED]; + + // List of indices for the teachers of this course. We are assuming that each + // teacher teaches separately. Must have the same number of elements as the + // number of sections list. + repeated int32 teacher_indices = 6; + + // The number of sections each teacher teaches of this course. Must have the + // same number of elements as the teacher index list. + repeated int32 teacher_section_counts = 7; + + // List of the possible rooms that this course can be assigned to. This can + // be empty. + repeated int32 room_indices = 8; +} + +message Teacher { + // Teacher name, used only for logging purposes. + string display_name = 1; + + // List of time slots that the teacher cannot be scheduled for. These time + // slot values index to the accumulative number of time slots starting at 0. + // For example, if a schedule rotation has 5 days and 8 time slots per day, + // and a teacher cannot be scheduled for the last time slot of the fourth + // day, the number here would be 31. + repeated int32 restricted_time_slots = 2; +} + +message Student { + // Student name, used only for logging purposes. + string display_name = 1; + + // List of course indices that the student needs to be enrolled in. + repeated int32 course_indices = 2; +} + +message Room { + // Room name, used only for logging purposes. + string display_name = 1; + + // Maximum number of students that can fit into this room. + int32 capacity = 2 [(google.api.field_behavior) = REQUIRED]; +} + +// Status returned by the solver. +enum CourseSchedulingResultStatus { + COURSE_SCHEDULING_RESULT_STATUS_UNSPECIFIED = 0; + + // The solver had enough time to find some solution that satisfies all + // constraints, but it did not prove optimality (which means it may or may + // not have reached the optimal). + // + // This can happen for large LP models (linear programming), and is a frequent + // response for time-limited MIPs (mixed integer programming). This is also + // what the CP (constraint programming) solver will return if there is no + // objective specified. + SOLVER_FEASIBLE = 1; + + // The solver found the proven optimal solution. + SOLVER_OPTIMAL = 2; + + // The model does not have any solution, according to the solver (which + // "proved" it, with the caveat that numerical proofs aren't actual proofs), + // or based on trivial considerations (eg. a variable whose lower bound is + // strictly greater than its upper bound). + SOLVER_INFEASIBLE = 3; + + // Model errors. These are always deterministic and repeatable. + // They should be accompanied with a string description of the error. + SOLVER_MODEL_INVALID = 4; + + // The model has not been solved in the given time or the solver was not able + // to solve the model given. + SOLVER_NOT_SOLVED = 5; + + // An error (either numerical or from a bug in the code) occurred. + ABNORMAL = 6; +} diff --git a/examples/cpp/course_scheduling_run.cc b/examples/cpp/course_scheduling_run.cc new file mode 100644 index 0000000000..fab20ca582 --- /dev/null +++ b/examples/cpp/course_scheduling_run.cc @@ -0,0 +1,108 @@ +// Copyright 2010-2018 Google LLC +// 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. + +// This file implements the main function for the Course Scheduling solver. It +// reads the problem specification from an input file specified via command-line +// flags, and prints the time slots for each course. See go/or-course-scheduling +// for more information. +// +// Example usage: +// ./course_scheduling_run --input_file=testdata/my_input_proto.textproto + +#include + +#include "examples/cpp/course_scheduling.h" +#include "examples/cpp/course_scheduling.pb.h" +#include "ortools/base/commandlineflags.h" +#include "ortools/base/timer.h" + +ABSL_FLAG(std::string, input, "", + "Input file containing a CourseSchedulingModel in text format."); + +namespace operations_research { + +void Main() { + CourseSchedulingModel input; + const auto proto_status = + file::GetTextProto(absl::GetFlag(FLAGS_input), &input, file::Defaults()); + if (!proto_status.ok()) { + LOG(ERROR) << proto_status.message(); + return; + } + + CourseSchedulingSolver solver; + WallTimer timer; + timer.Start(); + const CourseSchedulingResult result = solver.Solve(input); + timer.Stop(); + + LOG(INFO) << "Solver result status: " + << CourseSchedulingResultStatus_Name(result.solver_status()) << ". " + << result.message(); + + for (const ClassAssignment& class_assignment : result.class_assignments()) { + const int course_index = class_assignment.course_index(); + const int section_number = class_assignment.section_number(); + + int teacher_index = 0; + const Course& course = input.courses(course_index); + int sections = 0; + for (int section_index = 0; + section_index < course.teacher_section_counts_size(); + ++section_index) { + sections += course.teacher_section_counts(section_index); + if (section_number < sections) { + teacher_index = course.teacher_indices(section_index); + break; + } + } + + LOG(INFO) << course.display_name(); + LOG(INFO) << " Section: " << section_number; + LOG(INFO) << " Teacher: " << input.teachers(teacher_index).display_name(); + for (int i = 0; i < class_assignment.time_slots_size(); ++i) { + if (input.rooms_size() > 0) { + LOG(INFO) + << " Scheduled for time slot " << class_assignment.time_slots(i) + << " in room " + << input.rooms(class_assignment.room_indices(i)).display_name(); + } else { + LOG(INFO) << " Scheduled for time slot " + << class_assignment.time_slots(i); + } + } + } + + for (const StudentAssignment& student_assignment : + result.student_assignments()) { + const int student_index = student_assignment.student_index(); + + LOG(INFO) << input.students(student_index).display_name(); + for (int i = 0; i < student_assignment.course_indices_size(); ++i) { + LOG(INFO) + << " " + << input.courses(student_assignment.course_indices(i)).display_name() + << " " << student_assignment.section_indices(i); + } + } + + LOG(INFO) << "Solved model in " << timer.GetDuration(); +} + +} // namespace operations_research + +int main(int argc, char** argv) { + absl::ParseCommandLine(argc, argv); + operations_research::Main(); + return EXIT_SUCCESS; +} diff --git a/examples/cpp/cvrp_disjoint_tw.cc b/examples/cpp/cvrp_disjoint_tw.cc index 9d7acf5994..dbc4465241 100644 --- a/examples/cpp/cvrp_disjoint_tw.cc +++ b/examples/cpp/cvrp_disjoint_tw.cc @@ -24,6 +24,7 @@ #include +#include "absl/random/random.h" #include "examples/cpp/cvrptw_lib.h" #include "google/protobuf/text_format.h" #include "ortools/base/commandlineflags.h" @@ -34,7 +35,6 @@ #include "ortools/constraint_solver/routing_parameters.h" #include "ortools/constraint_solver/routing_parameters.pb.h" -using operations_research::ACMRandom; using operations_research::Assignment; using operations_research::DefaultRoutingSearchParameters; using operations_research::GetSeed; @@ -73,8 +73,7 @@ int main(int argc, char** argv) { << "Specify a non-null vehicle fleet size."; // VRP of size absl::GetFlag(FLAGS_vrp_size). // Nodes are indexed from 0 to absl::GetFlag(FLAGS_vrp_orders), the starts and - // ends of - // the routes are at node 0. + // ends of the routes are at node 0. const RoutingIndexManager::NodeIndex kDepot(0); RoutingIndexManager manager(absl::GetFlag(FLAGS_vrp_orders) + 1, absl::GetFlag(FLAGS_vrp_vehicles), kDepot); @@ -109,8 +108,8 @@ int main(int argc, char** argv) { routing.RegisterTransitCallback([&demand, &manager](int64 i, int64 j) { return demand.Demand(manager.IndexToNode(i), manager.IndexToNode(j)); }), - kNullCapacitySlack, kVehicleCapacity, /*fix_start_cumul_to_zero=*/true, - kCapacity); + kNullCapacitySlack, kVehicleCapacity, + /*fix_start_cumul_to_zero=*/true, kCapacity); // Adding time dimension constraints. const int64 kTimePerDemandUnit = 300; @@ -132,12 +131,12 @@ int main(int argc, char** argv) { // Adding disjoint time windows. Solver* solver = routing.solver(); - ACMRandom randomizer( + std::mt19937 randomizer( GetSeed(absl::GetFlag(FLAGS_vrp_use_deterministic_random_seed))); for (int order = 1; order < manager.num_nodes(); ++order) { std::vector forbid_points(2 * absl::GetFlag(FLAGS_vrp_windows), 0); for (int i = 0; i < forbid_points.size(); ++i) { - forbid_points[i] = randomizer.Uniform(kHorizon); + forbid_points[i] = absl::Uniform(randomizer, 0, kHorizon); } std::sort(forbid_points.begin(), forbid_points.end()); std::vector forbid_starts(1, 0); diff --git a/examples/cpp/cvrptw.cc b/examples/cpp/cvrptw.cc index 14ff29fc20..3037c82684 100644 --- a/examples/cpp/cvrptw.cc +++ b/examples/cpp/cvrptw.cc @@ -23,18 +23,17 @@ #include +#include "absl/random/random.h" #include "examples/cpp/cvrptw_lib.h" #include "google/protobuf/text_format.h" #include "ortools/base/commandlineflags.h" #include "ortools/base/integral_types.h" #include "ortools/base/logging.h" -#include "ortools/base/random.h" #include "ortools/constraint_solver/routing.h" #include "ortools/constraint_solver/routing_index_manager.h" #include "ortools/constraint_solver/routing_parameters.h" #include "ortools/constraint_solver/routing_parameters.pb.h" -using operations_research::ACMRandom; using operations_research::Assignment; using operations_research::DefaultRoutingSearchParameters; using operations_research::GetSeed; @@ -71,8 +70,7 @@ int main(int argc, char** argv) { << "Specify a non-null vehicle fleet size."; // VRP of size absl::GetFlag(FLAGS_vrp_size). // Nodes are indexed from 0 to absl::GetFlag(FLAGS_vrp_orders), the starts and - // ends of - // the routes are at node 0. + // ends of the routes are at node 0. const RoutingIndexManager::NodeIndex kDepot(0); RoutingIndexManager manager(absl::GetFlag(FLAGS_vrp_orders) + 1, absl::GetFlag(FLAGS_vrp_vehicles), kDepot); @@ -107,8 +105,8 @@ int main(int argc, char** argv) { routing.RegisterTransitCallback([&demand, &manager](int64 i, int64 j) { return demand.Demand(manager.IndexToNode(i), manager.IndexToNode(j)); }), - kNullCapacitySlack, kVehicleCapacity, /*fix_start_cumul_to_zero=*/true, - kCapacity); + kNullCapacitySlack, kVehicleCapacity, + /*fix_start_cumul_to_zero=*/true, kCapacity); // Adding time dimension constraints. const int64 kTimePerDemandUnit = 300; @@ -129,11 +127,12 @@ int main(int argc, char** argv) { const RoutingDimension& time_dimension = routing.GetDimensionOrDie(kTime); // Adding time windows. - ACMRandom randomizer( + std::mt19937 randomizer( GetSeed(absl::GetFlag(FLAGS_vrp_use_deterministic_random_seed))); const int64 kTWDuration = 5 * 3600; for (int order = 1; order < manager.num_nodes(); ++order) { - const int64 start = randomizer.Uniform(kHorizon - kTWDuration); + const int64 start = + absl::Uniform(randomizer, 0, kHorizon - kTWDuration); time_dimension.CumulVar(order)->SetRange(start, start + kTWDuration); } diff --git a/examples/cpp/jobshop_sat.cc b/examples/cpp/jobshop_sat.cc index bcce5f105b..5b6335c107 100644 --- a/examples/cpp/jobshop_sat.cc +++ b/examples/cpp/jobshop_sat.cc @@ -15,7 +15,9 @@ #include #include +#include "absl/flags/flag.h" #include "absl/strings/match.h" +#include "absl/strings/str_join.h" #include "google/protobuf/text_format.h" #include "google/protobuf/wrappers.pb.h" #include "ortools/base/commandlineflags.h" @@ -23,6 +25,7 @@ #include "ortools/base/timer.h" #include "ortools/data/jobshop_scheduling.pb.h" #include "ortools/data/jobshop_scheduling_parser.h" +#include "ortools/graph/connected_components.h" #include "ortools/sat/cp_model.h" #include "ortools/sat/cp_model.pb.h" #include "ortools/sat/model.h" @@ -32,8 +35,20 @@ ABSL_FLAG(std::string, params, "", "Sat parameters in text proto format."); ABSL_FLAG(bool, use_optional_variables, true, "Whether we use optional variables for bounds of an optional " "interval or not."); +ABSL_FLAG(bool, use_interval_makespan, true, + "Whether we encode the makespan using an interval or not."); +ABSL_FLAG(bool, use_expanded_precedences, false, + "Whether we add precedences between alternative tasks within the " + "same job."); +ABSL_FLAG( + bool, use_cumulative_relaxation, true, + "Whether we regroup multiple machines to create a cumulative relaxation."); +ABSL_FLAG( + int, job_suffix_relaxation_length, 5, + "The maximum length of the suffix of a job used in the linear relaxation."); ABSL_FLAG(bool, display_model, false, "Display jobshop proto before solving."); ABSL_FLAG(bool, display_sat_model, false, "Display sat proto before solving."); +ABSL_FLAG(int, horizon, -1, "Override horizon computation."); using operations_research::data::jssp::Job; using operations_research::data::jssp::JobPrecedence; @@ -88,193 +103,262 @@ int64 ComputeHorizon(const JsspInputProblem& problem) { // TODO(user): Uses transitions. } -// Solve a JobShop scheduling problem using SAT. -void Solve(const JsspInputProblem& problem) { - if (absl::GetFlag(FLAGS_display_model)) { - LOG(INFO) << problem.DebugString(); - } +// A job is a sequence of tasks. For each task, we store the main interval, as +// well as its start, size, and end variables. +struct JobTaskData { + IntervalVar interval; + IntVar start; + IntVar duration; + IntVar end; +}; - CpModelBuilder cp_model; +// Each task in a job can have multiple alternative ways of being performed. +// This structure stores the start, end, and presence variables attached to one +// alternative for a given task. +struct AlternativeTaskData { + IntervalVar interval; + IntVar start; + IntVar end; + BoolVar presence; +}; +// Create the job structure as a chain of tasks. Fills in the job_to_tasks +// vector. +void CreateJobs(const JsspInputProblem& problem, int64 horizon, + std::vector>& job_to_tasks, + bool& has_variable_duration_tasks, CpModelBuilder& cp_model) { const int num_jobs = problem.jobs_size(); - const int num_machines = problem.machines_size(); - const int64 horizon = ComputeHorizon(problem); - - std::vector starts; - std::vector ends; - - const Domain all_horizon(0, horizon); - - const IntVar makespan = cp_model.NewIntVar(all_horizon); - - std::vector > machine_to_intervals(num_machines); - std::vector > machine_to_jobs(num_machines); - std::vector > machine_to_starts(num_machines); - std::vector > machine_to_ends(num_machines); - std::vector > machine_to_presences(num_machines); - std::vector job_starts(num_jobs); - std::vector job_ends(num_jobs); - std::vector task_starts; - int64 objective_offset = 0; - std::vector objective_vars; - std::vector objective_coeffs; for (int j = 0; j < num_jobs; ++j) { const Job& job = problem.jobs(j); - IntVar previous_end; + const int num_tasks_in_job = job.tasks_size(); + std::vector& task_data = job_to_tasks[j]; + const int64 hard_start = job.has_earliest_start() ? job.earliest_start().value() : 0L; const int64 hard_end = job.has_latest_end() ? job.latest_end().value() : horizon; - for (int t = 0; t < job.tasks_size(); ++t) { + for (int t = 0; t < num_tasks_in_job; ++t) { const Task& task = job.tasks(t); const int num_alternatives = task.machine_size(); CHECK_EQ(num_alternatives, task.duration_size()); - // Add the "main" task interval. It will englobe all the alternative ones - // if there is many, or be a normal task otherwise. + // Add the "main" task interval. + std::vector durations; int64 min_duration = task.duration(0); int64 max_duration = task.duration(0); - for (int i = 1; i < num_alternatives; ++i) { - min_duration = std::min(min_duration, task.duration(i)); - max_duration = std::max(max_duration, task.duration(i)); + durations.push_back(task.duration(0)); + for (int a = 1; a < num_alternatives; ++a) { + min_duration = std::min(min_duration, task.duration(a)); + max_duration = std::max(max_duration, task.duration(a)); + durations.push_back(task.duration(a)); } + + if (min_duration != max_duration) has_variable_duration_tasks = true; + const IntVar start = cp_model.NewIntVar(Domain(hard_start, hard_end)); - const IntVar duration = - cp_model.NewIntVar(Domain(min_duration, max_duration)); + const IntVar duration = cp_model.NewIntVar(Domain::FromValues(durations)); const IntVar end = cp_model.NewIntVar(Domain(hard_start, hard_end)); const IntervalVar interval = cp_model.NewIntervalVar(start, duration, end); - // Store starts and ends of jobs for precedences. - if (t == 0) { - job_starts[j] = start; - } - if (t == job.tasks_size() - 1) { - job_ends[j] = end; - } - task_starts.push_back(start); + // Fill in job_to_tasks. + task_data.push_back({interval, start, duration, end}); // Chain the task belonging to the same job. if (t > 0) { - cp_model.AddLessOrEqual(previous_end, start); + cp_model.AddLessOrEqual(task_data[t - 1].end, task_data[t].start); + } + } + } +} + +// For each task of each jobs, create the alternative tasks and link them to the +// main task of the job. +void CreateAlternativeTasks( + const JsspInputProblem& problem, + const std::vector>& job_to_tasks, int64 horizon, + std::vector>>& + job_task_to_alternatives, + CpModelBuilder& cp_model) { + const int num_jobs = problem.jobs_size(); + const BoolVar true_var = cp_model.TrueVar(); + + for (int j = 0; j < num_jobs; ++j) { + const Job& job = problem.jobs(j); + const int num_tasks_in_job = job.tasks_size(); + job_task_to_alternatives[j].resize(num_tasks_in_job); + const std::vector& tasks = job_to_tasks[j]; + + const int64 hard_start = + job.has_earliest_start() ? job.earliest_start().value() : 0L; + const int64 hard_end = + job.has_latest_end() ? job.latest_end().value() : horizon; + + for (int t = 0; t < num_tasks_in_job; ++t) { + const Task& task = job.tasks(t); + const int num_alternatives = task.machine_size(); + CHECK_EQ(num_alternatives, task.duration_size()); + std::vector& alt_data = + job_task_to_alternatives[j][t]; + + absl::flat_hash_map> duration_supports; + duration_supports[task.duration(0)].push_back(0); + for (int a = 1; a < num_alternatives; ++a) { + duration_supports[task.duration(a)].push_back(a); } - previous_end = end; if (num_alternatives == 1) { - const int m = task.machine(0); - machine_to_intervals[m].push_back(interval); - machine_to_jobs[m].push_back(j); - machine_to_starts[m].push_back(start); - machine_to_ends[m].push_back(end); - machine_to_presences[m].push_back(cp_model.TrueVar()); - if (task.cost_size() > 0) { - objective_offset += task.cost(0); - } + alt_data.push_back( + {tasks[t].interval, tasks[t].start, tasks[t].end, true_var}); } else { - std::vector presences; for (int a = 0; a < num_alternatives; ++a) { - const BoolVar presence = cp_model.NewBoolVar(); + const BoolVar local_presence = cp_model.NewBoolVar(); const IntVar local_start = absl::GetFlag(FLAGS_use_optional_variables) ? cp_model.NewIntVar(Domain(hard_start, hard_end)) - : start; + : tasks[t].start; const IntVar local_duration = cp_model.NewConstant(task.duration(a)); const IntVar local_end = absl::GetFlag(FLAGS_use_optional_variables) ? cp_model.NewIntVar(Domain(hard_start, hard_end)) - : end; + : tasks[t].end; const IntervalVar local_interval = cp_model.NewOptionalIntervalVar( - local_start, local_duration, local_end, presence); + local_start, local_duration, local_end, local_presence); // Link local and global variables. if (absl::GetFlag(FLAGS_use_optional_variables)) { - cp_model.AddEquality(start, local_start).OnlyEnforceIf(presence); - cp_model.AddEquality(end, local_end).OnlyEnforceIf(presence); + cp_model.AddEquality(tasks[t].start, local_start) + .OnlyEnforceIf(local_presence); + cp_model.AddEquality(tasks[t].end, local_end) + .OnlyEnforceIf(local_presence); // TODO(user): Experiment with the following implication. - cp_model.AddEquality(duration, local_duration) - .OnlyEnforceIf(presence); + cp_model.AddEquality(tasks[t].duration, task.duration(a)) + .OnlyEnforceIf(local_presence); } - // Record relevant variables for later use. - const int m = task.machine(a); - machine_to_intervals[m].push_back(local_interval); - machine_to_jobs[m].push_back(j); - machine_to_starts[m].push_back(local_start); - machine_to_ends[m].push_back(local_end); - machine_to_presences[m].push_back(presence); - - // Add cost if present. - if (task.cost_size() > 0) { - objective_vars.push_back(presence); - objective_coeffs.push_back(task.cost(a)); - } - // Collect presence variables. - presences.push_back(presence); + alt_data.push_back( + {local_interval, local_start, local_end, local_presence}); } // Exactly one alternative interval is present. - cp_model.AddEquality(LinearExpr::BooleanSum(presences), 1); + std::vector interval_presences; + for (const AlternativeTaskData& alternative : alt_data) { + interval_presences.push_back(alternative.presence); + } + cp_model.AddEquality(LinearExpr::BooleanSum(interval_presences), 1); + + // Implement supporting literals for the duration of the main interval. + if (duration_supports.size() > 1) { // duration is not fixed. + for (const auto& duration_alternative_indices : duration_supports) { + const int64 value = duration_alternative_indices.first; + const BoolVar duration_eq_value = cp_model.NewBoolVar(); + + // duration_eq_value <=> duration == value. + cp_model.AddEquality(tasks[t].duration, value) + .OnlyEnforceIf(duration_eq_value); + cp_model.AddNotEqual(tasks[t].duration, value) + .OnlyEnforceIf(duration_eq_value.Not()); + // Implement the support part. If all literals pointing to the same + // duration are false, then the duration cannot take this value. + std::vector support_clause; + for (const int a : duration_alternative_indices.second) { + support_clause.push_back(alt_data[a].presence); + } + support_clause.push_back(duration_eq_value.Not()); + cp_model.AddBoolOr(support_clause); + } + } + } + + // Chain the alternative tasks belonging to the same job. + if (t > 0 && absl::GetFlag(FLAGS_use_expanded_precedences)) { + const std::vector& prev_data = + job_task_to_alternatives[j][t - 1]; + const std::vector& curr_data = + job_task_to_alternatives[j][t]; + for (int p = 0; p < prev_data.size(); ++p) { + const IntVar previous_end = prev_data[p].end; + const BoolVar previous_presence = prev_data[p].presence; + for (int c = 0; c < curr_data.size(); ++c) { + const IntVar current_start = curr_data[c].start; + const BoolVar current_presence = curr_data[c].presence; + cp_model.AddLessOrEqual(previous_end, current_start) + .OnlyEnforceIf({previous_presence, current_presence}); + } + } } } + } +} - // The makespan will be greater than the end of each job. - if (problem.makespan_cost_per_time_unit() != 0L) { - cp_model.AddLessOrEqual(previous_end, makespan); - } +// Tasks or alternative tasks are added to machines one by one. +// This structure records the characteristics of each task added on a machine. +// This information is indexed on each vector by the order of addition. +struct MachineTaskData { + IntervalVar interval; + int job; + IntVar start; + int64 duration; + IntVar end; + BoolVar presence; +}; +void CreateMachines( + const JsspInputProblem& problem, + const std::vector>>& + job_task_to_alternatives, + IntervalVar makespan_interval, CpModelBuilder& cp_model) { + const int num_jobs = problem.jobs_size(); + const int num_machines = problem.machines_size(); + std::vector> machine_to_tasks(num_machines); - const int64 lateness_penalty = job.lateness_cost_per_time_unit(); - // Lateness cost. - if (lateness_penalty != 0L) { - const int64 due_date = job.late_due_date(); - if (due_date == 0) { - objective_vars.push_back(previous_end); - objective_coeffs.push_back(lateness_penalty); - } else { - const IntVar shifted_var = - cp_model.NewIntVar(Domain(-due_date, horizon - due_date)); - cp_model.AddEquality(shifted_var, - LinearExpr(previous_end).AddConstant(-due_date)); - const IntVar lateness_var = cp_model.NewIntVar(all_horizon); - cp_model.AddMaxEquality(lateness_var, - {cp_model.NewConstant(0), shifted_var}); - objective_vars.push_back(lateness_var); - objective_coeffs.push_back(lateness_penalty); - } - } - const int64 earliness_penalty = job.earliness_cost_per_time_unit(); - // Earliness cost. - if (earliness_penalty != 0L) { - const int64 due_date = job.early_due_date(); - if (due_date > 0) { - const IntVar shifted_var = - cp_model.NewIntVar(Domain(due_date - horizon, due_date)); - cp_model.AddEquality(LinearExpr::Sum({shifted_var, previous_end}), - due_date); - const IntVar earliness_var = cp_model.NewIntVar(all_horizon); - cp_model.AddMaxEquality(earliness_var, - {cp_model.NewConstant(0), shifted_var}); - objective_vars.push_back(earliness_var); - objective_coeffs.push_back(earliness_penalty); + // Fills in the machine data vector. + for (int j = 0; j < num_jobs; ++j) { + const Job& job = problem.jobs(j); + const int num_tasks_in_job = job.tasks_size(); + + for (int t = 0; t < num_tasks_in_job; ++t) { + const Task& task = job.tasks(t); + const int num_alternatives = task.machine_size(); + CHECK_EQ(num_alternatives, task.duration_size()); + const std::vector& alt_data = + job_task_to_alternatives[j][t]; + + for (int a = 0; a < num_alternatives; ++a) { + // Record relevant variables for later use. + machine_to_tasks[task.machine(a)].push_back( + {alt_data[a].interval, j, alt_data[a].start, task.duration(a), + alt_data[a].end, alt_data[a].presence}); } } } // Add one no_overlap constraint per machine. for (int m = 0; m < num_machines; ++m) { - cp_model.AddNoOverlap(machine_to_intervals[m]); + std::vector intervals; + for (const MachineTaskData& task : machine_to_tasks[m]) { + intervals.push_back(task.interval); + } + if (absl::GetFlag(FLAGS_use_interval_makespan) && + problem.makespan_cost_per_time_unit() != 0L) { + intervals.push_back(makespan_interval); + } + cp_model.AddNoOverlap(intervals); + } + // Add transition times if needed. + for (int m = 0; m < num_machines; ++m) { if (problem.machines(m).has_transition_time_matrix()) { + const int num_intervals = machine_to_tasks[m].size(); const TransitionTimeMatrix& transitions = problem.machines(m).transition_time_matrix(); - const int num_intervals = machine_to_intervals[m].size(); // Create circuit constraint on a machine. // Node 0 and num_intervals + 1 are source and sink. CircuitConstraint circuit = cp_model.AddCircuitConstraint(); for (int i = 0; i < num_intervals; ++i) { - const int job_i = machine_to_jobs[m][i]; + const int job_i = machine_to_tasks[m][i].job; // Source to nodes. circuit.AddArc(0, i + 1, cp_model.NewBoolVar()); // Node to sink. @@ -282,14 +366,14 @@ void Solve(const JsspInputProblem& problem) { // Node to node. for (int j = 0; j < num_intervals; ++j) { if (i == j) { - circuit.AddArc(i + 1, i + 1, Not(machine_to_presences[m][i])); + circuit.AddArc(i + 1, i + 1, Not(machine_to_tasks[m][i].presence)); } else { - const int job_j = machine_to_jobs[m][j]; + const int job_j = machine_to_tasks[m][i].job; const int64 transition = transitions.transition_time(job_i * num_jobs + job_j); const BoolVar lit = cp_model.NewBoolVar(); - const IntVar start = machine_to_starts[m][j]; - const IntVar end = machine_to_ends[m][i]; + const IntVar start = machine_to_tasks[m][j].start; + const IntVar end = machine_to_tasks[m][i].end; circuit.AddArc(i + 1, j + 1, lit); // Push the new start with an extra transition. cp_model @@ -300,46 +384,326 @@ void Solve(const JsspInputProblem& problem) { } } } +} - // Add job precedences. - for (const JobPrecedence& precedence : problem.precedences()) { - const IntVar start = job_starts[precedence.second_job_index()]; - const IntVar end = job_ends[precedence.first_job_index()]; - cp_model.AddLessOrEqual(LinearExpr(end).AddConstant(precedence.min_delay()), - start); +// Collect all objective terms and add them to the model. +void CreateObjective( + const JsspInputProblem& problem, + const std::vector>& job_to_tasks, + const std::vector>>& + job_task_to_alternatives, + int64 horizon, IntVar makespan, CpModelBuilder& cp_model) { + int64 objective_offset = 0; + std::vector objective_vars; + std::vector objective_coeffs; + + const int num_jobs = problem.jobs_size(); + for (int j = 0; j < num_jobs; ++j) { + const Job& job = problem.jobs(j); + const int num_tasks_in_job = job.tasks_size(); + + // Add the cost associated to each task. + for (int t = 0; t < num_tasks_in_job; ++t) { + const Task& task = job.tasks(t); + const int num_alternatives = task.machine_size(); + + for (int a = 0; a < num_alternatives; ++a) { + // Add cost if present. + if (task.cost_size() > 0) { + objective_vars.push_back(job_task_to_alternatives[j][t][a].presence); + objective_coeffs.push_back(task.cost(a)); + } + } + } + + // Job lateness cost. + const int64 lateness_penalty = job.lateness_cost_per_time_unit(); + if (lateness_penalty != 0L) { + const int64 due_date = job.late_due_date(); + const IntVar job_end = job_to_tasks[j].back().end; + if (due_date == 0) { + objective_vars.push_back(job_end); + objective_coeffs.push_back(lateness_penalty); + } else { + const IntVar lateness_var = cp_model.NewIntVar(Domain(0, horizon)); + cp_model.AddLinMaxEquality( + lateness_var, + {LinearExpr(0), LinearExpr(job_end).AddConstant(-due_date)}); + objective_vars.push_back(lateness_var); + objective_coeffs.push_back(lateness_penalty); + } + } + + // Job earliness cost. + const int64 earliness_penalty = job.earliness_cost_per_time_unit(); + if (earliness_penalty != 0L) { + const int64 due_date = job.early_due_date(); + const IntVar job_end = job_to_tasks[j].back().end; + + if (due_date > 0) { + const IntVar earliness_var = cp_model.NewIntVar(Domain(0, horizon)); + cp_model.AddLinMaxEquality( + earliness_var, + {LinearExpr(0), + LinearExpr::Term(job_end, -1).AddConstant(due_date)}); + objective_vars.push_back(earliness_var); + objective_coeffs.push_back(earliness_penalty); + } + } } - // Add objective. + // Makespan objective. if (problem.makespan_cost_per_time_unit() != 0L) { objective_coeffs.push_back(problem.makespan_cost_per_time_unit()); objective_vars.push_back(makespan); } + + // Add the objective to the model. cp_model.Minimize(LinearExpr::ScalProd(objective_vars, objective_coeffs) .AddConstant(objective_offset)); if (problem.has_scaling_factor()) { cp_model.ScaleObjectiveBy(problem.scaling_factor().value()); } +} + +// This is a relaxation of the problem where we only consider the main tasks, +// and not the alternate copies. +void AddCumulativeRelaxation( + const JsspInputProblem& problem, + const std::vector>& job_to_tasks, + IntervalVar makespan_interval, CpModelBuilder& cp_model) { + const int num_jobs = problem.jobs_size(); + const int num_machines = problem.machines_size(); + + // Build a graph where two machines are connected if they appear in the same + // set of alternate machines for a given task. + std::vector> neighbors(num_machines); + for (int j = 0; j < num_jobs; ++j) { + const Job& job = problem.jobs(j); + const int num_tasks_in_job = job.tasks_size(); + for (int t = 0; t < num_tasks_in_job; ++t) { + const Task& task = job.tasks(t); + for (int a = 1; a < task.machine_size(); ++a) { + neighbors[task.machine(0)].insert(task.machine(a)); + } + } + } + + // Search for connected components in the above graph. + std::vector components = + util::GetConnectedComponents(num_machines, neighbors); + absl::flat_hash_map> machines_per_component; + for (int c = 0; c < components.size(); ++c) { + machines_per_component[components[c]].push_back(c); + } + + const IntVar one = cp_model.NewConstant(1); + for (const auto& it : machines_per_component) { + // Ignore the trivial cases. + if (it.second.size() < 2 || it.second.size() == num_machines) continue; + + LOG(INFO) << "Found machine connected component: [" + << absl::StrJoin(it.second, ", ") << "]"; + absl::flat_hash_set component(it.second.begin(), it.second.end()); + const IntVar capacity = cp_model.NewConstant(component.size()); + int num_intervals_in_cumulative = 0; + CumulativeConstraint cumul = cp_model.AddCumulative(capacity); + for (int j = 0; j < num_jobs; ++j) { + const Job& job = problem.jobs(j); + const int num_tasks_in_job = job.tasks_size(); + for (int t = 0; t < num_tasks_in_job; ++t) { + const Task& task = job.tasks(t); + for (const int m : task.machine()) { + if (component.contains(m)) { + cumul.AddDemand(job_to_tasks[j][t].interval, one); + num_intervals_in_cumulative++; + break; + } + } + } + } + if (absl::GetFlag(FLAGS_use_interval_makespan)) { + cumul.AddDemand(makespan_interval, capacity); + } + LOG(INFO) << " - created cumulative with " << num_intervals_in_cumulative + << " intervals"; + } +} + +// There are two linear redundant constraints. +// +// The first one states that the sum of durations of all tasks is a lower bound +// of the makespan * number of machines. +// +// The second one takes a suffix of one job chain, and states that the start of +// the suffix + the sum of all task durations in the suffix is a lower bound of +// the makespan. +void AddMakespanRedundantConstraints( + const JsspInputProblem& problem, + const std::vector>& job_to_tasks, IntVar makespan, + bool has_variable_duration_tasks, CpModelBuilder& cp_model) { + const int num_jobs = problem.jobs_size(); + const int num_machines = problem.machines_size(); + + // Global energetic reasoning. + std::vector all_task_durations; + for (const std::vector& tasks : job_to_tasks) { + for (const JobTaskData& task : tasks) { + all_task_durations.push_back(task.duration); + } + } + cp_model.AddLessOrEqual(LinearExpr::Sum(all_task_durations), + LinearExpr::Term(makespan, num_machines)); + + // Suffix linear equations. + if (has_variable_duration_tasks) { + for (int j = 0; j < num_jobs; ++j) { + const int job_length = job_to_tasks[j].size(); + const int start_suffix = std::max( + 0, job_length - absl::GetFlag(FLAGS_job_suffix_relaxation_length)); + for (int first_t = start_suffix; first_t + 1 < job_length; ++first_t) { + std::vector terms = {job_to_tasks[j][first_t].start}; + for (int t = first_t; t < job_length; ++t) { + terms.push_back(job_to_tasks[j][t].duration); + } + cp_model.AddLessOrEqual(LinearExpr::Sum(terms), makespan); + } + } + } +} + +// Solve a JobShop scheduling problem using CP-SAT. +void Solve(const JsspInputProblem& problem) { + if (absl::GetFlag(FLAGS_display_model)) { + LOG(INFO) << problem.DebugString(); + } + + CpModelBuilder cp_model; + + // Compute an over estimate of the horizon. + const int64 horizon = absl::GetFlag(FLAGS_horizon) != -1 + ? absl::GetFlag(FLAGS_horizon) + : ComputeHorizon(problem); + + // Create the main job structure. + const int num_jobs = problem.jobs_size(); + std::vector> job_to_tasks(num_jobs); + bool has_variable_duration_tasks = false; + CreateJobs(problem, horizon, job_to_tasks, has_variable_duration_tasks, + cp_model); + + // For each task of each jobs, create the alternative copies if needed and + // fill in the AlternativeTaskData vector. + std::vector>> + job_task_to_alternatives(num_jobs); + CreateAlternativeTasks(problem, job_to_tasks, horizon, + job_task_to_alternatives, cp_model); + + // Create the makespan variable and interval. + // If this flag is true, we will add to each no overlap constraint a special + // "makespan interval" that must necessarily be last by construction. This + // gives us a better lower bound on the makespan because this way we known + // that it must be after all other intervals in each no-overlap constraint. + // + // Otherwise, we will just add precence constraints between the last task of + // each job and the makespan variable. Alternatively, we could have added a + // precedence relation between all tasks and the makespan for a similar + // propagation thanks to our "precedence" propagator in the dijsunctive but + // that was slower than the interval trick when I tried. + const IntVar makespan = cp_model.NewIntVar(Domain(0, horizon)); + IntervalVar makespan_interval; + if (absl::GetFlag(FLAGS_use_interval_makespan)) { + makespan_interval = cp_model.NewIntervalVar( + /*start=*/makespan, + /*size=*/cp_model.NewIntVar(Domain(1, horizon)), + /*end=*/cp_model.NewIntVar(Domain(horizon + 1))); + } else if (problem.makespan_cost_per_time_unit() != 0L) { + for (int j = 0; j < num_jobs; ++j) { + // The makespan will be greater than the end of each job. + // This is not needed if we add the makespan "interval" to each + // disjunctive. + cp_model.AddLessOrEqual(job_to_tasks[j].back().end, makespan); + } + } + + // Machine constraints. + CreateMachines(problem, job_task_to_alternatives, makespan_interval, + cp_model); + + // Try to detect connected components of alternative machines. + // If this is happens, we can add a cumulative constraint as a relaxation of + // all no_ovelap constraints on the set of alternative machines. + if (absl::GetFlag(FLAGS_use_cumulative_relaxation)) { + AddCumulativeRelaxation(problem, job_to_tasks, makespan_interval, cp_model); + } + + // Various redundant constraints. They are mostly here to improve the LP + // relaxation. + if (problem.makespan_cost_per_time_unit() != 0L) { + AddMakespanRedundantConstraints(problem, job_to_tasks, makespan, + has_variable_duration_tasks, cp_model); + } + + // Add job precedences. + for (const JobPrecedence& precedence : problem.precedences()) { + const IntVar start = + job_to_tasks[precedence.second_job_index()].front().start; + const IntVar end = job_to_tasks[precedence.first_job_index()].back().end; + cp_model.AddLessOrEqual(LinearExpr(end).AddConstant(precedence.min_delay()), + start); + } + + // Objective. + CreateObjective(problem, job_to_tasks, job_task_to_alternatives, horizon, + makespan, cp_model); // Decision strategy. - cp_model.AddDecisionStrategy(task_starts, + std::vector all_task_starts; + for (const std::vector& job : job_to_tasks) { + for (const JobTaskData& task : job) { + all_task_starts.push_back(task.start); + } + } + cp_model.AddDecisionStrategy(all_task_starts, DecisionStrategyProto::CHOOSE_LOWEST_MIN, DecisionStrategyProto::SELECT_MIN_VALUE); - LOG(INFO) << "#machines:" << num_machines; + // Display problem statistics. + int num_tasks = 0; + int num_tasks_with_variable_duration = 0; + int num_tasks_with_alternatives = 0; + for (const std::vector& job : job_to_tasks) { + num_tasks += job.size(); + for (const JobTaskData& task : job) { + if (task.duration.Proto().domain_size() != 2 || + task.duration.Proto().domain(0) != task.duration.Proto().domain(1)) { + num_tasks_with_variable_duration++; + } + } + } + for (const std::vector>& + task_to_alternatives : job_task_to_alternatives) { + for (const std::vector& alternatives : + task_to_alternatives) { + if (alternatives.size() > 1) num_tasks_with_alternatives++; + } + } + + LOG(INFO) << "#machines:" << problem.machines_size(); LOG(INFO) << "#jobs:" << num_jobs; LOG(INFO) << "horizon:" << horizon; + LOG(INFO) << "#tasks: " << num_tasks; + LOG(INFO) << "#tasks with alternative: " << num_tasks_with_alternatives; + LOG(INFO) << "#tasks with variable duration: " + << num_tasks_with_variable_duration; if (absl::GetFlag(FLAGS_display_sat_model)) { LOG(INFO) << cp_model.Proto().DebugString(); } - LOG(INFO) << CpModelStats(cp_model.Proto()); - Model model; model.Add(NewSatParameters(absl::GetFlag(FLAGS_params))); - const CpSolverResponse response = SolveCpModel(cp_model.Build(), &model); - LOG(INFO) << CpSolverResponseStats(response); // Abort if we don't have any solution. if (response.status() != CpSolverStatus::OPTIMAL && @@ -350,18 +714,20 @@ void Solve(const JsspInputProblem& problem) { int64 final_cost = 0; if (problem.makespan_cost_per_time_unit() != 0) { int64 makespan = 0; - for (IntVar v : job_ends) { - makespan = std::max(makespan, SolutionIntegerValue(response, v)); + for (const std::vector& tasks : job_to_tasks) { + const IntVar job_end = tasks.back().end; + makespan = std::max(makespan, SolutionIntegerValue(response, job_end)); } final_cost += makespan * problem.makespan_cost_per_time_unit(); } - for (int i = 0; i < job_ends.size(); ++i) { - const int64 early_due_date = problem.jobs(i).early_due_date(); - const int64 late_due_date = problem.jobs(i).late_due_date(); - const int64 early_penalty = problem.jobs(i).earliness_cost_per_time_unit(); - const int64 late_penalty = problem.jobs(i).lateness_cost_per_time_unit(); - const int64 end = SolutionIntegerValue(response, job_ends[i]); + for (int j = 0; j < num_jobs; ++j) { + const int64 early_due_date = problem.jobs(j).early_due_date(); + const int64 late_due_date = problem.jobs(j).late_due_date(); + const int64 early_penalty = problem.jobs(j).earliness_cost_per_time_unit(); + const int64 late_penalty = problem.jobs(j).lateness_cost_per_time_unit(); + const int64 end = + SolutionIntegerValue(response, job_to_tasks[j].back().end); if (end < early_due_date && early_penalty != 0) { final_cost += (early_due_date - end) * early_penalty; } diff --git a/examples/cpp/multi_knapsack_sat.cc b/examples/cpp/multi_knapsack_sat.cc index 594c4f995c..0db7139abb 100644 --- a/examples/cpp/multi_knapsack_sat.cc +++ b/examples/cpp/multi_knapsack_sat.cc @@ -21,6 +21,7 @@ #include +#include "absl/flags/flag.h" #include "ortools/base/commandlineflags.h" #include "ortools/base/logging.h" #include "ortools/sat/cp_model.h" @@ -49,7 +50,7 @@ void MultiKnapsackSat(int scaling, const std::string& params) { const int num_items = scaling * kNumItems; const int num_bins = scaling; - std::vector > items_in_bins(num_bins); + std::vector> items_in_bins(num_bins); for (int b = 0; b < num_bins; ++b) { for (int i = 0; i < num_items; ++i) { items_in_bins[b].push_back(builder.NewBoolVar()); diff --git a/examples/cpp/sat_cnf_reader.h b/examples/cpp/sat_cnf_reader.h index 66415b2406..889071f8c8 100644 --- a/examples/cpp/sat_cnf_reader.h +++ b/examples/cpp/sat_cnf_reader.h @@ -186,7 +186,7 @@ class SatCnfReader { // Since the problem name is not stored in the cnf format, we infer it from // the file name. static std::string ExtractProblemName(const std::string& filename) { - const int found = filename.find_last_of("/"); + const int found = filename.find_last_of('/'); const std::string problem_name = found != std::string::npos ? filename.substr(found + 1) : filename; return problem_name; diff --git a/examples/cpp/sat_runner.cc b/examples/cpp/sat_runner.cc index 09aafddc79..19fa4b763c 100644 --- a/examples/cpp/sat_runner.cc +++ b/examples/cpp/sat_runner.cc @@ -18,6 +18,7 @@ #include #include +#include "absl/flags/flag.h" #include "absl/memory/memory.h" #include "absl/status/status.h" #include "absl/strings/match.h" @@ -239,6 +240,7 @@ int Run() { // The SAT competition requires a particular exit code and since we don't // really use it for any other purpose, we comply. + if (response.status() == CpSolverStatus::OPTIMAL) return 10; if (response.status() == CpSolverStatus::FEASIBLE) return 10; if (response.status() == CpSolverStatus::INFEASIBLE) return 20; return EXIT_SUCCESS; @@ -300,7 +302,7 @@ int Run() { CHECK(!absl::GetFlag(FLAGS_reduce_memory_usage)) << "incompatible"; CHECK(!absl::GetFlag(FLAGS_presolve)) << "incompatible"; LOG(INFO) << "Finding symmetries of the problem."; - std::vector > generators; + std::vector> generators; FindLinearBooleanProblemSymmetries(problem, &generators); std::unique_ptr propagator(new SymmetryPropagator); for (int i = 0; i < generators.size(); ++i) { @@ -366,7 +368,7 @@ int Run() { if (result == SatSolver::FEASIBLE) { if (absl::GetFlag(FLAGS_fu_malik) || absl::GetFlag(FLAGS_linear_scan) || absl::GetFlag(FLAGS_wpm1) || absl::GetFlag(FLAGS_core_enc)) { - printf("s OPTIMUM FOUND\n"); + absl::PrintF("s OPTIMUM FOUND\n"); CHECK(!solution.empty()); const Coefficient objective = ComputeObjectiveValue(problem, solution); scaled_best_bound = AddOffsetAndScaleObjectiveValue(problem, objective); @@ -377,13 +379,13 @@ int Run() { problem = original_problem; } } else { - printf("s SATISFIABLE\n"); + absl::PrintF("s SATISFIABLE\n"); } // Check and output the solution. CHECK(IsAssignmentValid(problem, solution)); if (absl::GetFlag(FLAGS_output_cnf_solution)) { - printf("v %s\n", SolutionString(problem, solution).c_str()); + absl::PrintF("v %s\n", SolutionString(problem, solution)); } if (!absl::GetFlag(FLAGS_output).empty()) { CHECK(!absl::GetFlag(FLAGS_reduce_memory_usage)) << "incompatible"; @@ -400,36 +402,32 @@ int Run() { } } if (result == SatSolver::INFEASIBLE) { - printf("s UNSATISFIABLE\n"); + absl::PrintF("s UNSATISFIABLE\n"); } // Print status. - printf("c status: %s\n", SatStatusString(result).c_str()); + absl::PrintF("c status: %s\n", SatStatusString(result)); // Print objective value. if (solution.empty()) { - printf("c objective: na\n"); - printf("c best bound: na\n"); + absl::PrintF("c objective: na\n"); + absl::PrintF("c best bound: na\n"); } else { const Coefficient objective = ComputeObjectiveValue(problem, solution); - printf("c objective: %.16g\n", - AddOffsetAndScaleObjectiveValue(problem, objective)); - printf("c best bound: %.16g\n", scaled_best_bound); + absl::PrintF("c objective: %.16g\n", + AddOffsetAndScaleObjectiveValue(problem, objective)); + absl::PrintF("c best bound: %.16g\n", scaled_best_bound); } // Print final statistics. - printf("c booleans: %d\n", solver->NumVariables()); + absl::PrintF("c booleans: %d\n", solver->NumVariables()); absl::PrintF("c conflicts: %d\n", solver->num_failures()); absl::PrintF("c branches: %d\n", solver->num_branches()); absl::PrintF("c propagations: %d\n", solver->num_propagations()); - printf("c walltime: %f\n", wall_timer.Get()); - printf("c usertime: %f\n", user_timer.Get()); - printf("c deterministic_time: %f\n", solver->deterministic_time()); + absl::PrintF("c walltime: %f\n", wall_timer.Get()); + absl::PrintF("c usertime: %f\n", user_timer.Get()); + absl::PrintF("c deterministic_time: %f\n", solver->deterministic_time()); - // The SAT competition requires a particular exit code and since we don't - // really use it for any other purpose, we comply. - if (result == SatSolver::FEASIBLE) return 10; - if (result == SatSolver::INFEASIBLE) return 20; return EXIT_SUCCESS; } diff --git a/examples/cpp/shift_minimization_sat.cc b/examples/cpp/shift_minimization_sat.cc index 622d9b2410..d29059233e 100644 --- a/examples/cpp/shift_minimization_sat.cc +++ b/examples/cpp/shift_minimization_sat.cc @@ -31,6 +31,7 @@ #include #include +#include "absl/flags/flag.h" #include "absl/strings/numbers.h" #include "absl/strings/str_split.h" #include "ortools/base/commandlineflags.h" @@ -63,10 +64,10 @@ class ShiftMinimizationParser { num_workers_read_(0) {} const std::vector& jobs() const { return jobs_; } - const std::vector >& possible_jobs_per_worker() const { + const std::vector>& possible_jobs_per_worker() const { return possible_jobs_per_worker_; } - const std::vector >& possible_assignments_per_job() + const std::vector>& possible_assignments_per_job() const { return possible_assignments_per_job_; } @@ -160,8 +161,8 @@ class ShiftMinimizationParser { } std::vector jobs_; - std::vector > possible_jobs_per_worker_; - std::vector > possible_assignments_per_job_; + std::vector> possible_jobs_per_worker_; + std::vector> possible_assignments_per_job_; LoadStatus load_status_; int declared_num_jobs_; int declared_num_workers_; @@ -186,8 +187,8 @@ void LoadAndSolve(const std::string& file_name) { const int num_jobs = jobs.size(); std::vector active_workers(num_workers); - std::vector > worker_job_vars(num_workers); - std::vector > possible_workers_per_job(num_jobs); + std::vector> worker_job_vars(num_workers); + std::vector> possible_workers_per_job(num_jobs); for (int w = 0; w < num_workers; ++w) { // Status variables for workers, are they active or not? @@ -235,7 +236,7 @@ void LoadAndSolve(const std::string& file_name) { // then the number of active workers on these jobs is equal to the number of // active jobs. std::set time_points; - std::set > visited_job_lists; + std::set> visited_job_lists; for (int j = 0; j < num_jobs; ++j) { time_points.insert(parser.jobs()[j].start); diff --git a/examples/cpp/sports_scheduling_sat.cc b/examples/cpp/sports_scheduling_sat.cc index 0203b7d8e4..a5123c488b 100644 --- a/examples/cpp/sports_scheduling_sat.cc +++ b/examples/cpp/sports_scheduling_sat.cc @@ -26,7 +26,7 @@ // - If team A meets team B, the reverse match cannot happen less that 6 weeks // after. // -// We model this problem with three matrices of variables, each with +// In the opponent model, we use three matrices of variables, each with // num_teams rows and 2*(num_teams - 1) columns: the var at position [i][j] // corresponds to the match of team #i at day #j. There are // 2*(num_teams - 1) columns because each team meets num_teams - 1 @@ -38,32 +38,38 @@ // - The 'signed_opponent' var [i][j] is the 'opponent' var [i][j] + // num_teams * the 'home_away' var [i][j]. // -// This aggregated variable will be useful to state constraints of the model -// and to do search on it. +// In the fixture model, we have a cube of Boolean variables fixtures. +// fixtures[d][i][j] is true if team i plays team j at home on day d. +// We also introduces a variable at_home[d][i] which is true if team i +// plays any opponent at home on day d. #include "absl/strings/str_cat.h" +#include "absl/strings/str_format.h" +#include "absl/strings/str_join.h" #include "ortools/base/commandlineflags.h" -#include "ortools/base/integral_types.h" #include "ortools/base/logging.h" #include "ortools/sat/cp_model.h" +#include "ortools/sat/cp_model.pb.h" +#include "ortools/sat/model.h" // Problem main flags. ABSL_FLAG(int, num_teams, 10, "Number of teams in the problem."); ABSL_FLAG(std::string, params, "", "Sat parameters."); +ABSL_FLAG(int, model, 1, "1 = opponent model, 2 = fixture model"); namespace operations_research { namespace sat { -void FirstModel(int num_teams) { +void OpponentModel(int num_teams) { const int num_days = 2 * num_teams - 2; const int kNoRematch = 6; CpModelBuilder builder; // Calendar variables. - std::vector > opponents(num_teams); - std::vector > home_aways(num_teams); - std::vector > signed_opponents(num_teams); + std::vector> opponents(num_teams); + std::vector> home_aways(num_teams); + std::vector> signed_opponents(num_teams); for (int t = 0; t < num_teams; ++t) { for (int d = 0; d < num_days; ++d) { @@ -185,7 +191,7 @@ void FirstModel(int num_teams) { } } -void SecondModel(int num_teams) { +void FixtureModel(int num_teams) { const int num_days = 2 * num_teams - 2; // const int kNoRematch = 6; const int matches_per_day = num_teams - 1; @@ -193,7 +199,7 @@ void SecondModel(int num_teams) { CpModelBuilder builder; // Does team i receive team j at home on day d? - std::vector > > fixtures(num_days); + std::vector>> fixtures(num_days); for (int d = 0; d < num_days; ++d) { fixtures[d].resize(num_teams); for (int i = 0; i < num_teams; ++i) { @@ -209,14 +215,14 @@ void SecondModel(int num_teams) { } // Is team t at home on day d? - std::vector > at_home(num_days); + std::vector> at_home(num_days); for (int d = 0; d < num_days; ++d) { for (int t = 0; t < num_teams; ++t) { at_home[d].push_back(builder.NewBoolVar()); } } - // Each day, Team t plays either at home or away. + // Each day, Team t plays another team, either at home or away. for (int d = 0; d < num_days; ++d) { for (int team = 0; team < num_teams; ++team) { std::vector possible_opponents; @@ -317,11 +323,15 @@ static const char kUsage[] = "There is no output besides the debug LOGs of the solver."; int main(int argc, char** argv) { - absl::SetProgramUsageMessage(kUsage); + absl::SetFlag(&FLAGS_logtostderr, true); absl::ParseCommandLine(argc, argv); CHECK_EQ(0, absl::GetFlag(FLAGS_num_teams) % 2) << "The number of teams must be even"; CHECK_GE(absl::GetFlag(FLAGS_num_teams), 2) << "At least 2 teams"; - operations_research::sat::SecondModel(absl::GetFlag(FLAGS_num_teams)); + if (absl::GetFlag(FLAGS_model) == 1) { + operations_research::sat::OpponentModel(absl::GetFlag(FLAGS_num_teams)); + } else { + operations_research::sat::FixtureModel(absl::GetFlag(FLAGS_num_teams)); + } return EXIT_SUCCESS; } diff --git a/examples/cpp/weighted_tardiness_sat.cc b/examples/cpp/weighted_tardiness_sat.cc index 22980b96aa..0a98b40646 100644 --- a/examples/cpp/weighted_tardiness_sat.cc +++ b/examples/cpp/weighted_tardiness_sat.cc @@ -16,6 +16,7 @@ #include #include +#include "absl/flags/flag.h" #include "absl/strings/match.h" #include "absl/strings/numbers.h" #include "absl/strings/str_join.h" @@ -118,8 +119,7 @@ void Solve(const std::vector& durations, // TODO(user): We can't set an objective upper bound with the current cp_model // interface, so we can't use heuristic or absl::GetFlag(FLAGS_upper_bound) - // here. The best is - // probably to provide a "solution hint" instead. + // here. The best is probably to provide a "solution hint" instead. // // Set a known upper bound (or use the flag). This has a bigger impact than // can be expected at first: @@ -176,7 +176,7 @@ void Solve(const std::vector& durations, for (int i = 0; i < num_tasks; ++i) { const int64 end = SolutionIntegerMin(r, task_ends[i]); CHECK_EQ(end, SolutionIntegerMax(r, task_ends[i])); - objective += weights[i] * std::max(0ll, end - due_dates[i]); + objective += weights[i] * std::max(int64{0}, end - due_dates[i]); } LOG(INFO) << "Cost " << objective;