this repo has no description
dekker.one
1#include <minizinc/solvers/nl/nl_components.hh>
2#include <minizinc/solvers/nl/nl_file.hh>
3
4namespace MiniZinc {
5
6/* *** *** *** NLBound *** *** *** */
7
8// Constructors
9NLBound::NLBound(Bound tag, double lb, double ub) : tag(tag), lb(lb), ub(ub) {}
10
11NLBound NLBound::make_bounded(double lb, double ub) {
12 if (lb == ub) {
13 return make_equal(lb);
14 } else {
15 assert(lb < ub);
16 return NLBound(LB_UB, lb, ub);
17 }
18}
19
20NLBound NLBound::make_ub_bounded(double ub) { return NLBound(UB, 0, ub); }
21
22NLBound NLBound::make_lb_bounded(double lb) { return NLBound(LB, lb, 0); }
23
24NLBound NLBound::make_nobound() { return NLBound(NONE, 0, 0); }
25
26NLBound NLBound::make_equal(double val) { return NLBound(EQ, val, val); }
27
28/** Update the lower bound */
29void NLBound::update_lb(double new_lb) {
30 switch (tag) {
31 // LB <= var <= UB. Same tag
32 case NLBound::LB_UB: {
33 assert(new_lb <= ub);
34 if (new_lb > lb) {
35 lb = new_lb;
36 }
37 break;
38 }
39 // var <= UB. Update tag
40 case NLBound::UB: {
41 assert(new_lb <= ub);
42 tag = LB_UB;
43 lb = new_lb;
44 break;
45 }
46 // LB <= var. Same tag
47 case NLBound::LB: {
48 if (new_lb > lb) {
49 lb = new_lb;
50 }
51 break;
52 }
53 // No bound. Update tag
54 case NLBound::NONE: {
55 tag = LB;
56 lb = new_lb;
57 break;
58 }
59 // LB = var = UB. Should not happen
60 case NLBound::EQ: {
61 should_not_happen(
62 "Updating a bound already set to \"equals\". We only allow tightening update.");
63 }
64 }
65}
66
67/** Update the upper bound */
68void NLBound::update_ub(double new_ub) {
69 switch (tag) {
70 // LB <= var <= UB. Same tag
71 case NLBound::LB_UB: {
72 assert(lb <= new_ub);
73 if (new_ub < ub) {
74 ub = new_ub;
75 }
76 break;
77 }
78 // var <= UB. Same tag
79 case NLBound::UB: {
80 if (new_ub < ub) {
81 ub = new_ub;
82 }
83 break;
84 }
85 // LB <= var. Update tag
86 case NLBound::LB: {
87 assert(lb <= new_ub);
88 tag = LB_UB;
89 ub = new_ub;
90 break;
91 }
92 // No bound. Update tag
93 case NLBound::NONE: {
94 tag = UB;
95 ub = new_ub;
96 break;
97 }
98 // LB = var = UB. Should not happen
99 case NLBound::EQ: {
100 should_not_happen(
101 "Updating a bound already set to \"equals\". We only allow tightening update.");
102 }
103 }
104}
105
106/** Update equal bound */
107void NLBound::update_eq(double new_eq) {
108 switch (tag) {
109 // LB <= var <= UB. Update tag
110 case NLBound::LB_UB: {
111 assert(lb <= new_eq && new_eq <= ub);
112 tag = EQ;
113 lb = new_eq;
114 ub = new_eq;
115 }
116 // var <= UB. Update tag
117 case NLBound::UB: {
118 assert(new_eq <= ub);
119 tag = EQ;
120 lb = new_eq;
121 ub = new_eq;
122 }
123 // LB <= var. Update tag
124 case NLBound::LB: {
125 assert(lb <= new_eq);
126 tag = EQ;
127 lb = new_eq;
128 ub = new_eq;
129 }
130 // No bound. Update tag
131 case NLBound::NONE: {
132 tag = EQ;
133 lb = new_eq;
134 ub = new_eq;
135 }
136 // LB = var = UB. Can happen only if we "update" with the same bound.
137 case NLBound::EQ: {
138 assert(lb == new_eq);
139 }
140 }
141}
142
143/** Printing with a name. */
144ostream& NLBound::print_on(ostream& os, string vname) const {
145 switch (tag) {
146 case LB_UB: {
147 os << "0 " << lb << " " << ub << " # " << lb << " =< " << vname << " =< " << ub;
148 break;
149 }
150 case UB: {
151 os << "1 " << ub << " # " << vname << " =< " << ub;
152 break;
153 }
154 case LB: {
155 os << "2 " << lb << " # " << lb << " =< " << vname;
156 break;
157 }
158 case NONE: {
159 os << "3"
160 << " # No constraint";
161 break;
162 }
163 case EQ: {
164 os << "4 " << lb << " # " << vname << " = " << lb;
165 break;
166 }
167 }
168 return os;
169}
170
171/** Default printing */
172ostream& NLBound::print_on(ostream& os) const {
173 print_on(os, "body");
174 return os;
175}
176
177/* *** *** *** NLVar *** *** *** */
178
179NLVar NLVar::copy_with_bound(NLBound bound) const {
180 NLVar v = NLVar(*this); // copy constructor
181 v.bound = bound;
182 return v;
183}
184
185/* *** *** *** NLArray *** *** *** */
186
187/* *** *** *** NLToken *** *** *** */
188
189const char* NLToken::get_name(OpCode oc) {
190 switch (oc) {
191 case OpCode::OPPLUS:
192 return "OPPLUS";
193 case OpCode::OPMINUS:
194 return "OPMINUS";
195 case OpCode::OPMULT:
196 return "OPMULT";
197 case OpCode::OPDIV:
198 return "OPDIV";
199 case OpCode::OPREM:
200 return "OPREM";
201 case OpCode::OPPOW:
202 return "OPPOW";
203 case OpCode::OPLESS:
204 return "OPLESS";
205 case OpCode::FLOOR:
206 return "FLOOR";
207 case OpCode::CEIL:
208 return "CEIL";
209 case OpCode::ABS:
210 return "ABS";
211 case OpCode::OPUMINUS:
212 return "OPUMINUS";
213 case OpCode::OPOR:
214 return "OPOR";
215 case OpCode::OPAND:
216 return "OPAND";
217 case OpCode::LT:
218 return "LT";
219 case OpCode::LE:
220 return "LE";
221 case OpCode::EQ:
222 return "EQ";
223 case OpCode::GE:
224 return "GE";
225 case OpCode::GT:
226 return "GT";
227 case OpCode::NE:
228 return "NE";
229 case OpCode::OPNOT:
230 return "OPNOT";
231 case OpCode::OPIFnl:
232 return "OPIFnl";
233 case OpCode::OP_tanh:
234 return "OP_tanh";
235 case OpCode::OP_tan:
236 return "OP_tan";
237 case OpCode::OP_sqrt:
238 return "OP_sqrt";
239 case OpCode::OP_sinh:
240 return "OP_sinh";
241 case OpCode::OP_sin:
242 return "OP_sin";
243 case OpCode::OP_log10:
244 return "OP_log10";
245 case OpCode::OP_log:
246 return "OP_log";
247 case OpCode::OP_exp:
248 return "OP_exp";
249 case OpCode::OP_cosh:
250 return "OP_cosh";
251 case OpCode::OP_cos:
252 return "OP_cos";
253 case OpCode::OP_atanh:
254 return "OP_atanh";
255 case OpCode::OP_atan2:
256 return "OP_atan2";
257 case OpCode::OP_atan:
258 return "OP_atan";
259 case OpCode::OP_asinh:
260 return "OP_asinh";
261 case OpCode::OP_asin:
262 return "OP_asin";
263 case OpCode::OP_acosh:
264 return "OP_acosh";
265 case OpCode::OP_acos:
266 return "OP_acos";
267 case OpCode::OPintDIV:
268 return "OPintDIV";
269 case OpCode::OPprecision:
270 return "OPprecision";
271 case OpCode::OPround:
272 return "OPround";
273 case OpCode::OPtrunc:
274 return "OPtrunc";
275 case OpCode::OPATLEAST:
276 return "OPATLEAST";
277 case OpCode::OPATMOST:
278 return "OPATMOST";
279 case OpCode::OPPLTERM:
280 return "OPPLTERM";
281 case OpCode::OPIFSYM:
282 return "OPIFSYM";
283 case OpCode::OPEXACTLY:
284 return "OPEXACTLY";
285 case OpCode::OPNOTATLEAST:
286 return "OPNOTATLEAST";
287 case OpCode::OPNOTATMOST:
288 return "OPNOTATMOST";
289 case OpCode::OPNOTEXACTLY:
290 return "OPNOTEXACTLY";
291 case OpCode::OPIMPELSE:
292 return "OPIMPELSE";
293 case OpCode::OP_IFF:
294 return "OP_IFF";
295 case OpCode::OPSOMESAME:
296 return "OPSOMESAME";
297 case OpCode::OP1POW:
298 return "OP1POW";
299 case OpCode::OP2POW:
300 return "OP2POW";
301 case OpCode::OPCPOW:
302 return "OPCPOW";
303 case OpCode::OPFUNCALL:
304 return "OPFUNCALL";
305 case OpCode::OPNUM:
306 return "OPNUM";
307 case OpCode::OPHOL:
308 return "OPHOL";
309 case OpCode::OPVARVAL:
310 return "OPVARVAL";
311 case OpCode::N_OPS:
312 return "N_OPS";
313 default:
314 assert(false);
315 return NULL;
316 }
317};
318
319const char* NLToken::get_name(MOpCode moc) {
320 switch (moc) {
321 case MOpCode::MINLIST:
322 return "MINLIST";
323 case MOpCode::MAXLIST:
324 return "MAXLIST";
325 case MOpCode::OPSUMLIST:
326 return "OPSUMLIST";
327 case MOpCode::OPCOUNT:
328 return "OPCOUNT";
329 case MOpCode::OPNUMBEROF:
330 return "OPNUMBEROF";
331 case MOpCode::OPNUMBEROFs:
332 return "OPNUMBEROFs";
333 case MOpCode::ANDLIST:
334 return "ANDLIST";
335 case MOpCode::ORLIST:
336 return "ORLIST";
337 case MOpCode::OPALLDIFF:
338 return "OPALLDIFF";
339 default:
340 assert(false);
341 return NULL;
342 }
343}
344
345NLToken NLToken::n(double value) {
346 NLToken tok;
347 tok.kind = Kind::NUMERIC;
348 tok.numeric_value = value;
349 return tok;
350}
351
352NLToken NLToken::v(string vname) {
353 NLToken tok;
354 tok.kind = Kind::VARIABLE;
355 tok.str = vname;
356 return tok;
357}
358
359NLToken NLToken::o(OpCode opc) {
360 NLToken tok;
361 tok.kind = Kind::OP;
362 tok.oc = opc;
363 return tok;
364}
365
366NLToken NLToken::mo(MOpCode mopc, int nb) {
367 NLToken tok;
368 tok.kind = Kind::MOP;
369 tok.moc = mopc;
370 tok.nb_args = nb;
371 return tok;
372}
373
374bool NLToken::is_variable() const { return kind == VARIABLE; }
375
376bool NLToken::is_constant() const { return kind == NUMERIC; }
377
378ostream& NLToken::print_on(ostream& os, const NLFile& nl_file) const {
379 switch (kind) {
380 case Kind::NUMERIC: {
381 os << "n" << numeric_value;
382 break;
383 }
384
385 case Kind::VARIABLE: {
386 os << "v" << nl_file.variable_indexes.at(str) << " # " << str;
387 break;
388 }
389
390 case Kind::STRING: {
391 should_not_happen("NL string token (Kind::STRING) not implemented");
392 }
393
394 case Kind::FUNCALL: {
395 should_not_happen("NL function call token (Kind::FUNCALL) not implemented");
396 }
397
398 case Kind::OP: {
399 os << "o" << oc << " # " << get_name(oc);
400 break;
401 }
402
403 case Kind::MOP: {
404 os << "o" << moc << " # " << get_name(moc) << " " << nb_args << endl;
405 os << nb_args;
406 break;
407 }
408
409 default:
410 should_not_happen("Unknown token kind: " << kind);
411 }
412
413 return os;
414}
415
416/* *** *** *** NLAlgCons *** *** *** */
417
418/** Method to build the var_coeff vector. */
419void NLAlgCons::set_jacobian(const vector<string>& vnames, const vector<double>& coeffs,
420 NLFile* nl_file) {
421 assert(vnames.size() == coeffs.size());
422 for (int i = 0; i < vnames.size(); ++i) {
423 string vn = vnames[i];
424 nl_file->variables.at(vn).jacobian_count++;
425 jacobian.push_back(pair<string, double>(vn, coeffs[i]));
426 }
427}
428
429/** A constraint is considered linear if the expression_graph is empty. */
430bool NLAlgCons::is_linear() const { return expression_graph.empty(); }
431
432/** Printing. */
433ostream& NLAlgCons::print_on(ostream& os, const NLFile& nl_file) const {
434 int idx = nl_file.constraint_indexes.at(name);
435
436 // Print the 'C' segment: if no expression graph, print "n0".
437 os << "C" << idx << " # Non linear part of " << name << endl;
438 if (expression_graph.empty()) {
439 os << "n0 # No non linear part coded as the value '0'" << endl;
440 } else {
441 for (auto t : expression_graph) {
442 t.print_on(os, nl_file);
443 os << endl;
444 }
445 }
446
447 // Print the 'J' segment if present.
448 if (!jacobian.empty()) {
449 os << "J" << idx << " " << jacobian.size() << " # Linear part of " << name << endl;
450 for (auto& vn_coef : jacobian) {
451 os << nl_file.variable_indexes.at(vn_coef.first) << " " << vn_coef.second << " # "
452 << vn_coef.first << endl;
453 }
454 }
455
456 return os;
457}
458
459/* *** *** *** NLLogicalCons *** *** *** */
460
461/** Printing. */
462ostream& NLLogicalCons::print_on(ostream& os, const NLFile& nl_file) const {
463 os << "L" << index << " # Logical constraint " << name << endl;
464 for (auto t : expression_graph) {
465 t.print_on(os, nl_file);
466 os << endl;
467 }
468
469 return os;
470}
471
472/* *** *** *** NLHeader *** *** *** */
473
474/** Printing the header.
475 * The header is composed of then lines that we describe as we proceed.
476 * A '#' starts a comment until the end of the line. However, it cannot be a line on its own!*/
477ostream& NLHeader::print_on(ostream& os, const NLFile& nl_file) const {
478 // 1st line:
479 // 'g': file will be in text format
480 // other numbers: as given in the doc (no other explanation...)
481 os << "g3 1 1 0" << endl;
482
483 // 2nd line:
484 os << nl_file.variables.size() << " " // Total number of variables
485 << nl_file.constraints.size()
486 << " " // Total number of algebraic constraint (including 'range' and 'eq')
487 << 1 << " " // Always 1 objective
488 << nl_file.nb_alg_cons_range << " " // Number of algebraic range constraints
489 << nl_file.nb_alg_cons_eq << " " // Number of algebraic eq constraints
490 << nl_file.logical_constraints.size() << " " // Number of logical constraints
491 << "# Total nb of: variables, algebraic constraints, objectives, ranges, eqs, logical "
492 "constraints"
493 << endl;
494
495 // 3rd line: Nonlinear and complementary information
496 os << nl_file.cnames_nl_general.size() << " " // Non linear constraints
497 << (nl_file.objective.is_linear() ? 0 : 1) << " " // Non linear objective
498 << "# Nb of nonlinear constraints, nonlinar objectives." << endl;
499 /* This was found in the online source of the ASL parser, but is not produce in our ampl tests.
500 * If needed, should be put on the same line
501 << nb_complementarity_linear_conditions << " "
502 << nb_complementarity_nonlinear_conditions << " "
503 << nb_complementarity_items_double_inequalities << " "
504 << nb_complemented_vars_non_zero_lowerbound << " "
505 << "Nb of complementary: linear & nonlinear conditions, double inequalities, vars with non-0 lower
506 bound."
507 << endl;
508 */
509
510 // 4th line: Network constraints
511 os << nl_file.cnames_nl_network.size() << " " // Number of nonlinear network constraints
512 << nl_file.cnames_lin_network.size() << " " // Number of linear network constraints
513 << "# Nb of network constraints: nonlinear, linear." << endl;
514
515 // 5th line: nonlinear variables:
516 os << nl_file.nlvc() << " " // Nb of nonlinear vars in constraints
517 << nl_file.nlvo() << " " // Nb of nonlinear vars in objectives
518 << nl_file.nlvb() << " " // Nb of nonlinear vars in both
519 << "# Nb of non linear vars in: constraints, objectives, both." << endl;
520
521 // 6th line:
522 os << nl_file.nwv() << " " // Nb of linear network vars
523 << "0"
524 << " " // Nb of functions. Not Implemented
525 << "0 1 "
526 << "# Nb of: linear network vars, functions. Floating point arithmetic mode (TEXT == 0). "
527 "Flag: if 1, add .sol suffixe."
528 << endl;
529
530 // 7th line: discrete variables
531 os << nl_file.nbv() << " " // Nb of linear binary vars
532 << nl_file.niv() << " " // Nb of linear integer vars
533 << nl_file.nlvbi() << " " // Nb of nonlinear integer vars in both
534 << nl_file.nlvci() << " " // Nb of nonlinear integer vars in constraints only
535 << nl_file.nlvoi() << " " // Nb of nonlinear integer vars in objectives only
536 << "# Nb of linear vars: binary, integer (non binary). "
537 << "Nb of nonlinear integer vars in: both, constraints only, objectives only." << endl;
538
539 // 8th line: non zeros
540 os << nl_file.jacobian_count() << " " // Nb of nonzero in jacobian
541 << nl_file.objective.gradient_count() << " " // Nb of nonzero in gradient
542 << "# Nb of non zeros in: jacobian, objective gradients." << endl;
543
544 // 9th line: name length. Our tests always produce 0...
545 os << "0"
546 << " " // Max constraint name length (??)
547 << "0"
548 << " " // Max var name length (??)
549 << "# Longest name among: contraints' name, vars' name." << endl;
550
551 // 10th line: common expressions. Not enough infor for now...
552 os << "0"
553 << " " // Nb common exprs in both
554 << "0"
555 << " " // Nb common exprs in constraints only
556 << "0"
557 << " " // Nb common exprs in objectives only
558 << "0"
559 << " " // Nb common exprs in single constraint only
560 << "0"
561 << " " // Nb common exprs in single objectives only
562 << "# Nb of common expressions in: both, constraints only, objectives only, single "
563 "constraint, single objective.";
564
565 return os;
566}
567
568/* *** *** *** NLObjective *** *** *** */
569
570/** Gradient count. */
571int NLObjective::gradient_count() const { return gradient.size(); }
572
573/** Set the gradient. */
574void NLObjective::set_gradient(const vector<string>& vnames, const vector<double>& coeffs) {
575 assert(vnames.size() == coeffs.size());
576 for (int i = 0; i < vnames.size(); ++i) {
577 string vn = vnames[i];
578 gradient.push_back(pair<string, double>(vn, coeffs[i]));
579 }
580}
581
582/** A objective is considered as linear if its expression graph is non empty. */
583bool NLObjective::is_defined() const { return minmax != UNDEF; }
584
585bool NLObjective::is_linear() const { return expression_graph.empty(); }
586
587bool NLObjective::is_optimisation() const { return minmax >= MINIMIZE; }
588
589/** Printing. */
590ostream& NLObjective::print_on(ostream& os, const NLFile& nl_file) const {
591 if (minmax != UNDEF) {
592 if (minmax == SATISFY) {
593 os << "O0 0 # Satisfy objectif implemented as 'minimize 0'" << endl;
594 os << "n0" << endl;
595 } else {
596 os << "O0 " << minmax << " # Objectif (0: minimize, 1: maximize)" << endl;
597 if (expression_graph.empty()) {
598 os << "n0 # No expression graph" << endl;
599 } else {
600 for (auto& tok : expression_graph) {
601 tok.print_on(os, nl_file) << endl;
602 }
603 }
604 // Print gradient
605 if (!gradient.empty()) {
606 os << "G0 " << gradient.size() << " # Objective Linear part" << endl;
607 for (auto& vn_coef : gradient) {
608 os << nl_file.variable_indexes.at(vn_coef.first) << " " << vn_coef.second << " # "
609 << vn_coef.first << endl;
610 }
611 }
612 }
613 }
614
615 return os;
616}
617} // namespace MiniZinc