1 // Copyright 2010-2021 Google LLC
2 // Licensed under the Apache License, Version 2.0 (the "License");
3 // you may not use this file except in compliance with the License.
4 // You may obtain a copy of the License at
5 //
6 // http://www.apache.org/licenses/LICENSE-2.0
7 //
8 // Unless required by applicable law or agreed to in writing, software
9 // distributed under the License is distributed on an "AS IS" BASIS,
10 // WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
11 // See the License for the specific language governing permissions and
12 // limitations under the License.
13
14 #include "ortools/linear_solver/model_validator.h"
15
16 #include <algorithm>
17 #include <cmath>
18 #include <limits>
19
20 #include "absl/container/flat_hash_map.h"
21 #include "absl/container/flat_hash_set.h"
22 #include "absl/status/status.h"
23 #include "absl/strings/match.h"
24 #include "absl/strings/str_cat.h"
25 #include "absl/strings/str_format.h"
26 #include "absl/types/optional.h"
27 #include "ortools/base/accurate_sum.h"
28 #include "ortools/base/commandlineflags.h"
29 #include "ortools/linear_solver/linear_solver.pb.h"
30 #include "ortools/port/file.h"
31 #include "ortools/port/proto_utils.h"
32 #include "ortools/util/fp_utils.h"
33 #include "ortools/util/lazy_mutable_copy.h"
34
35 ABSL_FLAG(
36 double, model_validator_infinity, 1e100,
37 "Anything above or equal to this magnitude will be considered infinity.");
38
39 namespace operations_research {
40 namespace {
41
IsNanOrAbsGreaterThanOrEqual(double value,double abs_value_threshold)42 bool IsNanOrAbsGreaterThanOrEqual(double value, double abs_value_threshold) {
43 return std::isnan(value) || std::abs(value) >= abs_value_threshold;
44 }
45
46 // Internal method to detect errors in bounds. The object passed as parameter
47 // must have "lower_bound" and "upper_bound" fields.
48 template <typename BoundedElement>
FindErrorInBounds(const BoundedElement & element,double abs_value_threshold,const bool accept_trivially_infeasible_bounds)49 std::string FindErrorInBounds(const BoundedElement& element,
50 double abs_value_threshold,
51 const bool accept_trivially_infeasible_bounds) {
52 if (std::isnan(element.lower_bound()) || std::isnan(element.upper_bound()) ||
53 element.lower_bound() >= abs_value_threshold ||
54 element.upper_bound() <= -abs_value_threshold ||
55 (!accept_trivially_infeasible_bounds &&
56 element.lower_bound() > element.upper_bound())) {
57 return absl::StrFormat("Infeasible bounds: [%f, %f]", element.lower_bound(),
58 element.upper_bound());
59 }
60 return "";
61 }
62
63 // Internal method to detect errors in a single variable.
FindErrorInMPVariable(const MPVariableProto & variable,double abs_value_threshold,const bool accept_trivially_infeasible_bounds)64 std::string FindErrorInMPVariable(
65 const MPVariableProto& variable, double abs_value_threshold,
66 const bool accept_trivially_infeasible_bounds) {
67 const std::string bound_error = FindErrorInBounds(
68 variable, abs_value_threshold, accept_trivially_infeasible_bounds);
69 if (!bound_error.empty()) return bound_error;
70
71 if (!accept_trivially_infeasible_bounds && variable.is_integer() &&
72 ceil(variable.lower_bound()) > floor(variable.upper_bound())) {
73 return absl::StrCat(
74 "Infeasible bounds for integer variable: [", (variable.lower_bound()),
75 ", ", (variable.upper_bound()), "]", " translate to the empty set");
76 }
77 if (IsNanOrAbsGreaterThanOrEqual(variable.objective_coefficient(),
78 abs_value_threshold)) {
79 return absl::StrCat("Invalid objective_coefficient: ",
80 (variable.objective_coefficient()));
81 }
82 return std::string();
83 }
84
85 // Returns an error message if 'var_indices' contains a duplicate index.
86 template <typename Iterable>
FindDuplicateVarIndex(const Iterable & var_indices,std::vector<bool> * var_mask)87 std::string FindDuplicateVarIndex(const Iterable& var_indices,
88 std::vector<bool>* var_mask) {
89 int duplicate_var_index = -1;
90 for (const int var_index : var_indices) {
91 if ((*var_mask)[var_index]) duplicate_var_index = var_index;
92 (*var_mask)[var_index] = true;
93 }
94 // Reset "var_mask" to all false, sparsely.
95 for (const int var_index : var_indices) {
96 (*var_mask)[var_index] = false;
97 }
98 if (duplicate_var_index >= 0) {
99 return absl::StrCat("var_index #", duplicate_var_index,
100 " appears several times");
101 }
102 return "";
103 }
104
105 // Internal method to detect errors in a single constraint.
106 // "var_mask" is a vector<bool> whose size is the number of variables in
107 // the model, and it will be all set to false before and after the call.
FindErrorInMPConstraint(const MPConstraintProto & constraint,std::vector<bool> * var_mask,double abs_value_threshold,const bool accept_trivially_infeasible_bounds)108 std::string FindErrorInMPConstraint(
109 const MPConstraintProto& constraint, std::vector<bool>* var_mask,
110 double abs_value_threshold, const bool accept_trivially_infeasible_bounds) {
111 const std::string bound_error = FindErrorInBounds(
112 constraint, abs_value_threshold, accept_trivially_infeasible_bounds);
113 if (!bound_error.empty()) return bound_error;
114
115 // TODO(user): clarify explicitly, at least in a comment, whether we want
116 // to accept empty constraints (i.e. without variables).
117
118 const int num_vars_in_model = var_mask->size();
119 const int num_vars_in_ct = constraint.var_index_size();
120 const int num_coeffs_in_ct = constraint.coefficient_size();
121 if (num_vars_in_ct != num_coeffs_in_ct) {
122 return absl::StrCat("var_index_size() != coefficient_size() (",
123 num_vars_in_ct, " VS ", num_coeffs_in_ct);
124 }
125 for (int i = 0; i < num_vars_in_ct; ++i) {
126 const int var_index = constraint.var_index(i);
127 if (var_index >= num_vars_in_model || var_index < 0) {
128 return absl::StrCat("var_index(", i, ")=", var_index,
129 " is out of bounds");
130 }
131 const double coeff = constraint.coefficient(i);
132 if (IsNanOrAbsGreaterThanOrEqual(coeff, abs_value_threshold)) {
133 return absl::StrCat("coefficient(", i, ")=", (coeff), " is invalid");
134 }
135 }
136
137 const std::string error =
138 FindDuplicateVarIndex(constraint.var_index(), var_mask);
139 if (!error.empty()) return error;
140
141 // We found no error, all is fine.
142 return std::string();
143 }
144
CroppedConstraintDebugString(const MPConstraintProto & constraint)145 std::string CroppedConstraintDebugString(const MPConstraintProto& constraint) {
146 const int kMaxPrintedVars = 10;
147
148 MPConstraintProto constraint_light = constraint;
149 std::string suffix_str;
150 if (constraint.var_index_size() > kMaxPrintedVars) {
151 constraint_light.mutable_var_index()->Truncate(kMaxPrintedVars);
152 absl::StrAppend(&suffix_str,
153 " (var_index cropped; size=", constraint.var_index_size(),
154 ").");
155 }
156 if (constraint.coefficient_size() > kMaxPrintedVars) {
157 constraint_light.mutable_coefficient()->Truncate(kMaxPrintedVars);
158 absl::StrAppend(&suffix_str, " (coefficient cropped; size=",
159 constraint.coefficient_size(), ").");
160 }
161 return absl::StrCat("Constraint proto: ",
162 ProtobufShortDebugString(constraint_light), suffix_str);
163 }
164
IsBoolean(const MPVariableProto & variable)165 bool IsBoolean(const MPVariableProto& variable) {
166 if (variable.lower_bound() < 0) return false;
167 if (variable.upper_bound() > 1) return false;
168 return variable.is_integer();
169 }
170
FindErrorInMPIndicatorConstraint(const MPModelProto & model,const MPIndicatorConstraint & indicator,std::vector<bool> * var_mask,double abs_value_threshold,bool accept_trivially_infeasible_bounds)171 std::string FindErrorInMPIndicatorConstraint(
172 const MPModelProto& model, const MPIndicatorConstraint& indicator,
173 std::vector<bool>* var_mask, double abs_value_threshold,
174 bool accept_trivially_infeasible_bounds) {
175 if (!indicator.has_var_index()) {
176 return "var_index is required.";
177 }
178 const int var_index = indicator.var_index();
179 if (var_index < 0 || var_index >= model.variable_size()) {
180 return absl::StrCat("var_index=", var_index, " is out of bounds.");
181 }
182 if (!IsBoolean(model.variable(var_index))) {
183 return absl::StrCat("var_index=", var_index, " is not Boolean.");
184 }
185 const int var_value = indicator.var_value();
186 if (var_value < 0 || var_value > 1) {
187 return absl::StrCat("var_value=", var_value, " must be 0 or 1.");
188 }
189 const MPConstraintProto& constraint = indicator.constraint();
190 std::string error =
191 FindErrorInMPConstraint(constraint, var_mask, abs_value_threshold,
192 accept_trivially_infeasible_bounds);
193 if (!error.empty()) {
194 // Constraint protos can be huge, theoretically. So we guard against
195 // that.
196 return absl::StrCat(error, " in constraint ",
197 CroppedConstraintDebugString(constraint));
198 }
199 return "";
200 }
201
FindErrorInMPSosConstraint(const MPModelProto & model,const MPSosConstraint & sos,std::vector<bool> * var_mask,double abs_value_threshold)202 std::string FindErrorInMPSosConstraint(const MPModelProto& model,
203 const MPSosConstraint& sos,
204 std::vector<bool>* var_mask,
205 double abs_value_threshold) {
206 if (sos.weight_size() != 0 && sos.weight_size() != sos.var_index_size()) {
207 return "weight_size() > 0 and var_index_size() != weight_size()";
208 }
209 for (const int var_index : sos.var_index()) {
210 if (var_index < 0 || var_index >= model.variable_size()) {
211 return absl::StrCat("var_index=", var_index, " is out of bounds.");
212 }
213 }
214 for (int i = 0; i < sos.weight_size(); ++i) {
215 if (IsNanOrAbsGreaterThanOrEqual(sos.weight(i), abs_value_threshold)) {
216 return absl::StrCat("Invalid weight: ", sos.weight(i));
217 }
218 if (i == 0) continue;
219 if (sos.weight(i - 1) >= sos.weight(i)) {
220 return "SOS weights must be strictly increasing";
221 }
222 }
223
224 const std::string error = FindDuplicateVarIndex(sos.var_index(), var_mask);
225 if (!error.empty()) return error;
226
227 return "";
228 }
229
FindErrorInMPQuadraticConstraint(const MPModelProto & model,const MPQuadraticConstraint & qcst,std::vector<bool> * var_mask,double abs_value_threshold,bool accept_trivially_infeasible_bounds)230 std::string FindErrorInMPQuadraticConstraint(
231 const MPModelProto& model, const MPQuadraticConstraint& qcst,
232 std::vector<bool>* var_mask, double abs_value_threshold,
233 bool accept_trivially_infeasible_bounds) {
234 const int num_vars = model.variable_size();
235
236 if (qcst.var_index_size() != qcst.coefficient_size()) {
237 return "var_index_size() != coefficient_size()";
238 }
239
240 const std::string bound_error = FindErrorInBounds(
241 qcst, abs_value_threshold, accept_trivially_infeasible_bounds);
242 if (!bound_error.empty()) return bound_error;
243
244 for (int i = 0; i < qcst.var_index_size(); ++i) {
245 if (qcst.var_index(i) < 0 || qcst.var_index(i) >= num_vars) {
246 return absl::StrCat("var_index(", i, ")=", qcst.var_index(i),
247 " is invalid.", " It must be in [0, ", num_vars, ")");
248 }
249
250 if (IsNanOrAbsGreaterThanOrEqual(qcst.coefficient(i),
251 abs_value_threshold)) {
252 return absl::StrCat("coefficient(", i, ")=", qcst.coefficient(i),
253 " is invalid");
254 }
255 }
256 const std::string duplicate_error =
257 FindDuplicateVarIndex(qcst.var_index(), var_mask);
258 if (!duplicate_error.empty()) return duplicate_error;
259
260 if (qcst.qvar1_index_size() != qcst.qvar2_index_size() ||
261 qcst.qvar1_index_size() != qcst.qcoefficient_size()) {
262 return "quadratic indices and coefficients must have the same size";
263 }
264 for (int i = 0; i < qcst.qvar1_index_size(); ++i) {
265 if (qcst.qvar1_index(i) >= num_vars || qcst.qvar1_index(i) < 0) {
266 return absl::StrCat("qvar1_index(", i, ")=", qcst.qvar1_index(i),
267 " is invalid.", " It must be in [0, ", num_vars, ")");
268 }
269
270 if (qcst.qvar2_index(i) >= num_vars || qcst.qvar2_index(i) < 0) {
271 return absl::StrCat("qvar2_index(", i, ")=", qcst.qvar2_index(i),
272 " is invalid.", " It must be in [0, ", num_vars, ")");
273 }
274
275 if (IsNanOrAbsGreaterThanOrEqual(qcst.qcoefficient(i),
276 abs_value_threshold)) {
277 return absl::StrCat("qcoefficient(", i, ")=", qcst.qcoefficient(i),
278 " is invalid");
279 }
280 }
281
282 return "";
283 }
284
FindErrorInMPAbsConstraint(const MPModelProto & model,const MPAbsConstraint & abs)285 std::string FindErrorInMPAbsConstraint(const MPModelProto& model,
286 const MPAbsConstraint& abs) {
287 if (!abs.has_var_index()) {
288 return "var_index is required.";
289 }
290 if (!abs.has_resultant_var_index()) {
291 return "resultant_var_index is required.";
292 }
293
294 const int num_vars = model.variable_size();
295 if (abs.var_index() < 0 || abs.var_index() >= num_vars) {
296 return absl::StrCat("var_index=", abs.var_index(), " is invalid.",
297 " It must be in [0, ", num_vars, ")");
298 }
299 if (abs.resultant_var_index() < 0 || abs.resultant_var_index() >= num_vars) {
300 return absl::StrCat("var_index=", abs.resultant_var_index(), " is invalid.",
301 " It must be in [0, ", num_vars, ")");
302 }
303 return "";
304 }
305
FindErrorInMPAndOrConstraint(const MPModelProto & model,const MPArrayConstraint & and_or)306 std::string FindErrorInMPAndOrConstraint(const MPModelProto& model,
307 const MPArrayConstraint& and_or) {
308 if (and_or.var_index_size() == 0) {
309 return "var_index cannot be empty.";
310 }
311 if (!and_or.has_resultant_var_index()) {
312 return "resultant_var_index is required.";
313 }
314
315 const int num_vars = model.variable_size();
316 for (int i = 0; i < and_or.var_index_size(); ++i) {
317 if (and_or.var_index(i) < 0 || and_or.var_index(i) >= num_vars) {
318 return absl::StrCat("var_index(", i, ")=", and_or.var_index(i),
319 " is invalid.", " It must be in [0, ", num_vars, ")");
320 }
321 if (!IsBoolean(model.variable(and_or.var_index(i)))) {
322 return absl::StrCat("var_index=", i, " is not Boolean.");
323 }
324 }
325 if (and_or.resultant_var_index() < 0 ||
326 and_or.resultant_var_index() >= num_vars) {
327 return absl::StrCat("resultant_var_index=", and_or.resultant_var_index(),
328 " is invalid.", " It must be in [0, ", num_vars, ")");
329 }
330 if (!IsBoolean(model.variable(and_or.resultant_var_index()))) {
331 return absl::StrCat("resultant_var_index is not Boolean.");
332 }
333 return "";
334 }
335
FindErrorInMPMinMaxConstraint(const MPModelProto & model,const MPArrayWithConstantConstraint & min_max,double abs_value_threshold)336 std::string FindErrorInMPMinMaxConstraint(
337 const MPModelProto& model, const MPArrayWithConstantConstraint& min_max,
338 double abs_value_threshold) {
339 if (min_max.var_index_size() == 0) {
340 return "var_index cannot be empty.";
341 }
342 if (!min_max.has_resultant_var_index()) {
343 return "resultant_var_index is required.";
344 }
345
346 if (IsNanOrAbsGreaterThanOrEqual(min_max.constant(), abs_value_threshold)) {
347 return absl::StrCat("Invalid constant: ", (min_max.constant()));
348 }
349
350 const int num_vars = model.variable_size();
351 for (int i = 0; i < min_max.var_index_size(); ++i) {
352 if (min_max.var_index(i) < 0 || min_max.var_index(i) >= num_vars) {
353 return absl::StrCat("var_index(", i, ")=", min_max.var_index(i),
354 " is invalid.", " It must be in [0, ", num_vars, ")");
355 }
356 }
357 if (min_max.resultant_var_index() < 0 ||
358 min_max.resultant_var_index() >= num_vars) {
359 return absl::StrCat("resultant_var_index=", min_max.resultant_var_index(),
360 " is invalid.", " It must be in [0, ", num_vars, ")");
361 }
362 return "";
363 }
364
FindErrorInQuadraticObjective(const MPQuadraticObjective & qobj,int num_vars,double abs_value_threshold)365 std::string FindErrorInQuadraticObjective(const MPQuadraticObjective& qobj,
366 int num_vars,
367 double abs_value_threshold) {
368 if (qobj.qvar1_index_size() != qobj.qvar2_index_size() ||
369 qobj.qvar1_index_size() != qobj.coefficient_size()) {
370 return "indices and coefficients must have the same size";
371 }
372
373 for (int i = 0; i < qobj.qvar1_index_size(); ++i) {
374 if (qobj.qvar1_index(i) >= num_vars || qobj.qvar1_index(i) < 0) {
375 return absl::StrCat("qvar1_index(", i, ")=", qobj.qvar1_index(i),
376 " is invalid.", " It must be in [0, ", num_vars, ")");
377 }
378
379 if (qobj.qvar2_index(i) >= num_vars || qobj.qvar2_index(i) < 0) {
380 return absl::StrCat("qvar2_index(", i, ")=", qobj.qvar2_index(i),
381 " is invalid.", " It must be in [0, ", num_vars, ")");
382 }
383
384 if (IsNanOrAbsGreaterThanOrEqual(qobj.coefficient(i),
385 abs_value_threshold)) {
386 return absl::StrCat("coefficient(", i, ")=", (qobj.coefficient(i)),
387 " is invalid");
388 }
389 }
390 return "";
391 }
392
FindErrorInSolutionHint(const PartialVariableAssignment & solution_hint,int num_vars,double abs_value_threshold)393 std::string FindErrorInSolutionHint(
394 const PartialVariableAssignment& solution_hint, int num_vars,
395 double abs_value_threshold) {
396 if (solution_hint.var_index_size() != solution_hint.var_value_size()) {
397 return absl::StrCat("var_index_size() != var_value_size() [",
398 solution_hint.var_index_size(), " VS ",
399 solution_hint.var_value_size());
400 }
401 std::vector<bool> var_in_hint(num_vars, false);
402 for (int i = 0; i < solution_hint.var_index_size(); ++i) {
403 const int var_index = solution_hint.var_index(i);
404 if (var_index >= num_vars || var_index < 0) {
405 return absl::StrCat("var_index(", i, ")=", var_index, " is invalid.",
406 " It must be in [0, ", num_vars, ")");
407 }
408 if (var_in_hint[var_index]) {
409 return absl::StrCat("Duplicate var_index = ", var_index);
410 }
411 var_in_hint[var_index] = true;
412 if (IsNanOrAbsGreaterThanOrEqual(solution_hint.var_value(i),
413 abs_value_threshold)) {
414 return absl::StrCat("var_value(", i, ")=", (solution_hint.var_value(i)),
415 " is invalid");
416 }
417 }
418 return std::string();
419 }
420
421 } // namespace
422
FindErrorInMPModelProto(const MPModelProto & model,double abs_value_threshold,const bool accept_trivially_infeasible_bounds)423 std::string FindErrorInMPModelProto(
424 const MPModelProto& model, double abs_value_threshold,
425 const bool accept_trivially_infeasible_bounds) {
426 // NOTE(user): Empty models are considered fine by this function, although
427 // it is not clear whether MPSolver::Solve() will always respond in the same
428 // way, depending on the solvers.
429 if (abs_value_threshold == 0.0) {
430 abs_value_threshold = absl::GetFlag(FLAGS_model_validator_infinity);
431 }
432
433 if (IsNanOrAbsGreaterThanOrEqual(model.objective_offset(),
434 abs_value_threshold)) {
435 return absl::StrCat("Invalid objective_offset: ",
436 (model.objective_offset()));
437 }
438 const int num_vars = model.variable_size();
439 const int num_cts = model.constraint_size();
440
441 // Validate variables.
442 std::string error;
443 for (int i = 0; i < num_vars; ++i) {
444 error = FindErrorInMPVariable(model.variable(i), abs_value_threshold,
445 accept_trivially_infeasible_bounds);
446 if (!error.empty()) {
447 return absl::StrCat("In variable #", i, ": ", error, ". Variable proto: ",
448 ProtobufShortDebugString(model.variable(i)));
449 }
450 }
451
452 // Validate constraints.
453 std::vector<bool> variable_appears(num_vars, false);
454 for (int i = 0; i < num_cts; ++i) {
455 const MPConstraintProto& constraint = model.constraint(i);
456 error = FindErrorInMPConstraint(constraint, &variable_appears,
457 abs_value_threshold,
458 accept_trivially_infeasible_bounds);
459 if (!error.empty()) {
460 // Constraint protos can be huge, theoretically. So we guard against that.
461 return absl::StrCat("In constraint #", i, ": ", error, ". ",
462 CroppedConstraintDebugString(constraint));
463 }
464 }
465
466 // Validate general constraints.
467 for (int i = 0; i < model.general_constraint_size(); ++i) {
468 const MPGeneralConstraintProto& gen_constraint =
469 model.general_constraint(i);
470 std::string error;
471 switch (gen_constraint.general_constraint_case()) {
472 case MPGeneralConstraintProto::kIndicatorConstraint:
473 error = FindErrorInMPIndicatorConstraint(
474 model, gen_constraint.indicator_constraint(), &variable_appears,
475 abs_value_threshold, accept_trivially_infeasible_bounds);
476 break;
477
478 case MPGeneralConstraintProto::kSosConstraint:
479 error =
480 FindErrorInMPSosConstraint(model, gen_constraint.sos_constraint(),
481 &variable_appears, abs_value_threshold);
482 break;
483
484 case MPGeneralConstraintProto::kQuadraticConstraint:
485 error = FindErrorInMPQuadraticConstraint(
486 model, gen_constraint.quadratic_constraint(), &variable_appears,
487 abs_value_threshold, accept_trivially_infeasible_bounds);
488 break;
489
490 case MPGeneralConstraintProto::kAbsConstraint:
491 error =
492 FindErrorInMPAbsConstraint(model, gen_constraint.abs_constraint());
493 break;
494
495 case MPGeneralConstraintProto::kAndConstraint:
496 error = FindErrorInMPAndOrConstraint(model,
497 gen_constraint.and_constraint());
498 break;
499
500 case MPGeneralConstraintProto::kOrConstraint:
501 error =
502 FindErrorInMPAndOrConstraint(model, gen_constraint.or_constraint());
503 break;
504
505 case MPGeneralConstraintProto::kMinConstraint:
506 error = FindErrorInMPMinMaxConstraint(
507 model, gen_constraint.min_constraint(), abs_value_threshold);
508 break;
509
510 case MPGeneralConstraintProto::kMaxConstraint:
511 error = FindErrorInMPMinMaxConstraint(
512 model, gen_constraint.max_constraint(), abs_value_threshold);
513 break;
514 default:
515 return absl::StrCat("Unknown general constraint type ",
516 gen_constraint.general_constraint_case());
517 }
518 if (!error.empty()) {
519 return absl::StrCat("In general constraint #", i, ": ", error);
520 }
521 }
522
523 // Validate objectives.
524 if (model.has_quadratic_objective()) {
525 error = FindErrorInQuadraticObjective(model.quadratic_objective(), num_vars,
526 abs_value_threshold);
527 if (!error.empty()) return absl::StrCat("In quadratic_objective: ", error);
528 }
529
530 // Validate the solution hint.
531 error = FindErrorInSolutionHint(model.solution_hint(), num_vars,
532 abs_value_threshold);
533 if (!error.empty()) {
534 return absl::StrCat("In solution_hint(): ", error);
535 }
536
537 return std::string();
538 }
539
540 absl::optional<LazyMutableCopy<MPModelProto>>
ExtractValidMPModelOrPopulateResponseStatus(const MPModelRequest & request,MPSolutionResponse * response)541 ExtractValidMPModelOrPopulateResponseStatus(const MPModelRequest& request,
542 MPSolutionResponse* response) {
543 CHECK(response != nullptr);
544
545 if (!request.has_model() && !request.has_model_delta()) {
546 response->set_status(MPSOLVER_OPTIMAL);
547 response->set_status_str("Requests without model are considered OPTIMAL");
548 return absl::nullopt;
549 }
550 if (request.has_model() && request.has_model_delta()) {
551 response->set_status(MPSOLVER_MODEL_INVALID);
552 response->set_status_str(
553 "Fields 'model' and 'model_delta' are mutually exclusive");
554 return absl::nullopt;
555 }
556
557 // Extract the baseline model.
558 LazyMutableCopy<MPModelProto> model(request.model());
559 if (request.has_model_delta()) {
560 // NOTE(user): This library needs to be portable, so we can't include
561 // ortools/base/file.h; see ../port/file.h.
562 std::string contents;
563 const absl::Status file_read_status = PortableFileGetContents(
564 request.model_delta().baseline_model_file_path(), &contents);
565 if (!file_read_status.ok()) {
566 response->set_status(MPSOLVER_MODEL_INVALID);
567 response->set_status_str(
568 "Error when reading model_delta.baseline_model_file_path: '" +
569 file_read_status.ToString());
570 return absl::nullopt;
571 }
572 if (!model.get_mutable()->ParseFromString(contents)) {
573 response->set_status(MPSOLVER_MODEL_INVALID);
574 response->set_status_str(
575 absl::StrFormat("The contents of baseline model file '%s' couldn't "
576 "be parsed as a raw serialized MPModelProto",
577 request.model_delta().baseline_model_file_path()));
578 return absl::nullopt;
579 }
580 }
581
582 // Validate the baseline model.
583 std::string error = FindErrorInMPModelProto(model.get());
584
585 // If the baseline is valid and we have a model delta, validate the delta,
586 // then apply it.
587 if (error.empty() && request.has_model_delta()) {
588 const MPModelDeltaProto& delta = request.model_delta();
589 error = FindErrorInMPModelDeltaProto(delta, model.get());
590 if (error.empty()) ApplyVerifiedMPModelDelta(delta, model.get_mutable());
591 }
592
593 // Deal with errors.
594 if (!error.empty()) {
595 if (request.enable_internal_solver_output()) {
596 LOG(ERROR) << absl::StrCat("Invalid model: ", error);
597 }
598 response->set_status(absl::StrContains(error, "Infeasible")
599 ? MPSOLVER_INFEASIBLE
600 : MPSOLVER_MODEL_INVALID);
601 response->set_status_str(error);
602 return absl::nullopt;
603 }
604
605 if (model.get().variable_size() == 0 && model.get().constraint_size() == 0 &&
606 model.get().general_constraint_size() == 0) {
607 response->set_status(MPSOLVER_OPTIMAL);
608 response->set_objective_value(model.get().objective_offset());
609 response->set_best_objective_bound(response->objective_value());
610 response->set_status_str(
611 "Requests without variables and constraints are considered OPTIMAL");
612 return absl::nullopt;
613 }
614
615 return std::move(model);
616 }
617
ExtractValidMPModelInPlaceOrPopulateResponseStatus(MPModelRequest * request,MPSolutionResponse * response)618 bool ExtractValidMPModelInPlaceOrPopulateResponseStatus(
619 MPModelRequest* request, MPSolutionResponse* response) {
620 absl::optional<LazyMutableCopy<MPModelProto>> lazy_copy =
621 ExtractValidMPModelOrPopulateResponseStatus(*request, response);
622 if (!lazy_copy) return false;
623 if (lazy_copy->was_copied()) {
624 lazy_copy->get_mutable()->Swap(request->mutable_model());
625 }
626 return true;
627 }
628
629 // TODO(user): Add a general FindFeasibilityErrorInSolution() and factor out the
630 // common code.
FindFeasibilityErrorInSolutionHint(const MPModelProto & model,double tolerance)631 std::string FindFeasibilityErrorInSolutionHint(const MPModelProto& model,
632 double tolerance) {
633 const int num_vars = model.variable_size();
634
635 // First, we validate the solution hint.
636 std::string error =
637 FindErrorInSolutionHint(model.solution_hint(), num_vars,
638 absl::GetFlag(FLAGS_model_validator_infinity));
639 if (!error.empty()) return absl::StrCat("Invalid solution_hint: ", error);
640
641 // Special error message for the empty case.
642 if (num_vars > 0 && model.solution_hint().var_index_size() == 0) {
643 return "Empty solution_hint.";
644 }
645
646 // To be feasible, the hint must not be partial.
647 if (model.solution_hint().var_index_size() != num_vars) {
648 return absl::StrCat("Partial solution_hint: only ",
649 model.solution_hint().var_index_size(), " out of the ",
650 num_vars, " problem variables are set.");
651 }
652
653 // All the values must be exactly in the variable bounds.
654 std::vector<double> var_value(num_vars);
655 for (int i = 0; i < model.solution_hint().var_index_size(); ++i) {
656 const int var_index = model.solution_hint().var_index(i);
657 const double value = model.solution_hint().var_value(i);
658 var_value[var_index] = value;
659 const double lb = model.variable(var_index).lower_bound();
660 const double ub = model.variable(var_index).upper_bound();
661 if (!IsSmallerWithinTolerance(value, ub, tolerance) ||
662 !IsSmallerWithinTolerance(lb, value, tolerance)) {
663 return absl::StrCat("Variable '", model.variable(var_index).name(),
664 "' is set to ", (value),
665 " which is not in the variable bounds [", (lb), ", ",
666 (ub), "] modulo a tolerance of ", (tolerance), ".");
667 }
668 }
669
670 // All the constraints must be satisfiable.
671 for (int cst_index = 0; cst_index < model.constraint_size(); ++cst_index) {
672 const MPConstraintProto& constraint = model.constraint(cst_index);
673 AccurateSum<double> activity;
674 for (int j = 0; j < constraint.var_index_size(); ++j) {
675 activity.Add(constraint.coefficient(j) *
676 var_value[constraint.var_index(j)]);
677 }
678 const double lb = model.constraint(cst_index).lower_bound();
679 const double ub = model.constraint(cst_index).upper_bound();
680 if (!IsSmallerWithinTolerance(activity.Value(), ub, tolerance) ||
681 !IsSmallerWithinTolerance(lb, activity.Value(), tolerance)) {
682 return absl::StrCat(
683 "Constraint '", model.constraint(cst_index).name(), "' has activity ",
684 (activity.Value()), " which is not in the constraint bounds [", (lb),
685 ", ", (ub), "] modulo a tolerance of ", (tolerance), ".");
686 }
687 }
688
689 return "";
690 }
691
FindErrorInMPModelDeltaProto(const MPModelDeltaProto & delta,const MPModelProto & model)692 std::string FindErrorInMPModelDeltaProto(const MPModelDeltaProto& delta,
693 const MPModelProto& model) {
694 const double abs_value_threshold =
695 absl::GetFlag(FLAGS_model_validator_infinity);
696 int num_vars = model.variable_size();
697 // Validate delta variables.
698 std::string error;
699 absl::flat_hash_set<int> new_var_indices;
700 int max_var_index = num_vars - 1;
701 MPVariableProto tmp_var_proto;
702 for (const auto& pair : delta.variable_overrides()) {
703 const int var_index = pair.first;
704 const MPVariableProto& var_override_proto = pair.second;
705 if (var_index < 0) {
706 error = "Invalid key";
707 } else if (var_index >= num_vars) {
708 max_var_index = std::max(max_var_index, var_index);
709 new_var_indices.insert(var_index);
710 error =
711 FindErrorInMPVariable(var_override_proto, abs_value_threshold,
712 /*accept_trivially_infeasible_bounds=*/false);
713 } else {
714 tmp_var_proto = model.variable(var_index);
715 // NOTE(user): It is OK for the override proto to be empty, i.e. be a
716 // non-override.
717 tmp_var_proto.MergeFrom(var_override_proto);
718 error =
719 FindErrorInMPVariable(tmp_var_proto, abs_value_threshold,
720 /*accept_trivially_infeasible_bounds=*/false);
721 }
722 if (!error.empty()) {
723 return absl::StrFormat(
724 "variable_overrides with key (eg. var index) = %d: %s", var_index,
725 error);
726 }
727 }
728 if (max_var_index != num_vars + new_var_indices.size() - 1) {
729 return absl::StrFormat(
730 "The added and existing variable indices do not form a dense integer "
731 "interval: oldmax=%d, max=%d, num added=%d",
732 num_vars - 1, max_var_index, new_var_indices.size());
733 }
734 // Now we "officially" add the new variables to "num_vars".
735 num_vars += new_var_indices.size();
736
737 // Validate delta constraints. We can avoid going over the full
738 // var_index/coefficient of the original constraint, since the overrides are
739 // self-sufficient (i.e. the override var_index/coefficients are valid iff
740 // they would be valid in a standalone, new constraint). So we use a partial
741 // proto merger to avoid those in the baseline constraint.
742 std::vector<bool> variable_appears(num_vars, false);
743 MPConstraintProto tmp_constraint_proto;
744 const int num_constraints = model.constraint_size();
745 absl::flat_hash_set<int> new_ct_indices;
746 int max_ct_index = num_constraints - 1;
747 for (const auto& pair : delta.constraint_overrides()) {
748 const int ct_index = pair.first;
749 const MPConstraintProto& constraint_override_proto = pair.second;
750 if (ct_index < 0) {
751 error = "Invalid constraint index";
752 } else if (ct_index >= num_constraints) {
753 max_ct_index = std::max(max_ct_index, ct_index);
754 new_ct_indices.insert(ct_index);
755 error = FindErrorInMPConstraint(
756 constraint_override_proto, &variable_appears, abs_value_threshold,
757 /*accept_trivially_infeasible_bounds=*/false);
758 } else {
759 // NOTE(user): We don't need to do the merging of var_index/coefficient:
760 // that part of the merged constraint will be valid iff the override is
761 // valid as a standalone var_index/coefficient map.
762 // So we simply validate a reduced version of the actual "merged"
763 // constraint, by removing the var_index/coefficient of the baseline.
764 // Benefit: the complexity is O(|constraint override|) even if the
765 // baseline constraint was huge.
766 tmp_constraint_proto.Clear();
767 MergeMPConstraintProtoExceptTerms(model.constraint(ct_index),
768 &tmp_constraint_proto);
769 tmp_constraint_proto.MergeFrom(constraint_override_proto);
770 error = FindErrorInMPConstraint(
771 tmp_constraint_proto, &variable_appears, abs_value_threshold,
772 /*accept_trivially_infeasible_bounds=*/false);
773 }
774 if (!error.empty()) {
775 return absl::StrFormat(
776 "constraint_overrides with key (eg. constraint index) = %d: %s",
777 ct_index, error);
778 }
779 }
780 if (max_ct_index != num_constraints + new_ct_indices.size() - 1) {
781 return absl::StrFormat(
782 "The added and existing constraint indices do not form a dense integer "
783 "interval: oldmax=%d, max=%d, num added=%d",
784 num_constraints - 1, max_ct_index, new_ct_indices.size());
785 }
786
787 return "";
788 }
789
MergeMPConstraintProtoExceptTerms(const MPConstraintProto & from,MPConstraintProto * to)790 void MergeMPConstraintProtoExceptTerms(const MPConstraintProto& from,
791 MPConstraintProto* to) {
792 #define COPY_FIELD_IF_PRESENT(field) \
793 if (from.has_##field()) to->set_##field(from.field())
794 COPY_FIELD_IF_PRESENT(lower_bound);
795 COPY_FIELD_IF_PRESENT(upper_bound);
796 COPY_FIELD_IF_PRESENT(name);
797 COPY_FIELD_IF_PRESENT(is_lazy);
798 #undef COPY_FIELD_IF_PRESENT
799 }
800
801 namespace {
PruneZeroTermsInMpConstraint(MPConstraintProto * ct)802 void PruneZeroTermsInMpConstraint(MPConstraintProto* ct) {
803 // Optimize the fast path (when no term is pruned) by doing a first quick scan
804 // until the first zero.
805 int first_zero = 0;
806 while (first_zero < ct->var_index_size() &&
807 ct->coefficient(first_zero) != 0.0) {
808 ++first_zero;
809 }
810 int num_kept = first_zero;
811 for (int i = first_zero; i < ct->var_index_size(); ++i) {
812 if (ct->coefficient(i) == 0.0) continue;
813 if (num_kept != i) {
814 ct->set_var_index(num_kept, ct->var_index(i));
815 ct->set_coefficient(num_kept, ct->coefficient(i));
816 }
817 ++num_kept;
818 }
819 ct->mutable_var_index()->Truncate(num_kept);
820 ct->mutable_coefficient()->Truncate(num_kept);
821 }
822
823 // Adds default entries to a repeated message field until it has the wanted
824 // size. We don't use google::protobuf::util::Resize() because it's not
825 // compatible with 'light' protos.
826 template <class T>
ExtendRepeatedPtrFieldToSize(const int size,T * repeated_messages)827 void ExtendRepeatedPtrFieldToSize(const int size, T* repeated_messages) {
828 DCHECK_GE(size, repeated_messages->size());
829 while (repeated_messages->size() < size) repeated_messages->Add();
830 }
831 } // namespace
832
ApplyVerifiedMPModelDelta(const MPModelDeltaProto & delta,MPModelProto * model)833 void ApplyVerifiedMPModelDelta(const MPModelDeltaProto& delta,
834 MPModelProto* model) {
835 // Apply the delta to the variables: first, resize the variable array.
836 int max_var_index = -1;
837 for (const auto& p : delta.variable_overrides()) {
838 max_var_index = std::max(max_var_index, p.first);
839 }
840 if (max_var_index >= model->variable_size()) {
841 ExtendRepeatedPtrFieldToSize(max_var_index + 1, model->mutable_variable());
842 }
843 // Then, apply the variable overrides.
844 for (const auto& p : delta.variable_overrides()) {
845 model->mutable_variable(p.first)->MergeFrom(p.second);
846 }
847
848 // Apply the delta to the constraints: first, resize the constraint array.
849 int max_ct_index = -1;
850 for (const auto& p : delta.constraint_overrides()) {
851 max_ct_index = std::max(max_ct_index, p.first);
852 }
853 const int old_num_constraints = model->constraint_size();
854 if (max_ct_index >= old_num_constraints) {
855 ExtendRepeatedPtrFieldToSize(max_ct_index + 1, model->mutable_constraint());
856 }
857 // Then, apply the constraint overrides.
858 for (const auto& p : delta.constraint_overrides()) {
859 const MPConstraintProto& override_ct = p.second;
860 MPConstraintProto* baseline = model->mutable_constraint(p.first);
861 // Fast path for added constraints.
862 if (p.first >= old_num_constraints) {
863 *baseline = override_ct;
864 continue;
865 }
866 MergeMPConstraintProtoExceptTerms(/*from=*/override_ct, /*to=*/baseline);
867 // Special case: the override is neutralized.
868 if (override_ct.has_lower_bound() &&
869 override_ct.lower_bound() <=
870 -absl::GetFlag(FLAGS_model_validator_infinity) &&
871 override_ct.has_upper_bound() &&
872 override_ct.upper_bound() >=
873 absl::GetFlag(FLAGS_model_validator_infinity)) {
874 baseline->clear_var_index();
875 baseline->clear_coefficient();
876 continue;
877 }
878 // Otherwise we have to apply the term overrides. We can't do that in less
879 // than O(|baseline| + |override_ct|) because the baseline doesn't have a
880 // lookup-friendly data structure. But we still try to do it as efficiently
881 // as possible. In particular, we only use O(|override_ct|) extra memory.
882 absl::flat_hash_map<int, double> term_overrides;
883 term_overrides.reserve(override_ct.var_index_size());
884 for (int i = 0; i < override_ct.var_index_size(); ++i) {
885 term_overrides[override_ct.var_index(i)] = override_ct.coefficient(i);
886 }
887 for (int i = 0; i < baseline->var_index_size(); ++i) {
888 auto it = term_overrides.find(baseline->var_index(i));
889 if (it == term_overrides.end()) continue;
890 baseline->set_coefficient(i, it->second);
891 it->second = 0.0; // To mark this term override as 'has been applied'.
892 }
893 PruneZeroTermsInMpConstraint(baseline);
894 // Add the term overrides which haven't been used: those are added terms.
895 for (const auto& p : term_overrides) {
896 if (p.second != 0.0) {
897 baseline->add_var_index(p.first);
898 baseline->add_coefficient(p.second);
899 }
900 }
901 }
902 }
903
904 } // namespace operations_research
905