|
1 /* -*- mode: C++; indent-tabs-mode: nil; -*- |
|
2 * |
|
3 * This file is a part of LEMON, a generic C++ optimization library. |
|
4 * |
|
5 * Copyright (C) 2003-2008 |
|
6 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport |
|
7 * (Egervary Research Group on Combinatorial Optimization, EGRES). |
|
8 * |
|
9 * Permission to use, modify and distribute this software is granted |
|
10 * provided that this copyright notice appears in all copies. For |
|
11 * precise terms see the accompanying LICENSE file. |
|
12 * |
|
13 * This software is provided "AS IS" with no warranty of any kind, |
|
14 * express or implied, and with no claim as to its suitability for any |
|
15 * purpose. |
|
16 * |
|
17 */ |
|
18 |
|
19 #ifndef LEMON_LP_BASE_H |
|
20 #define LEMON_LP_BASE_H |
|
21 |
|
22 #include<iostream> |
|
23 #include<vector> |
|
24 #include<map> |
|
25 #include<limits> |
|
26 #include<lemon/math.h> |
|
27 |
|
28 #include<lemon/core.h> |
|
29 #include<lemon/bits/lp_id.h> |
|
30 |
|
31 ///\file |
|
32 ///\brief The interface of the LP solver interface. |
|
33 ///\ingroup lp_group |
|
34 namespace lemon { |
|
35 |
|
36 /// Function to decide whether a floating point value is finite or not. |
|
37 |
|
38 /// Retruns true if the argument is not infinity, minus infinity or NaN. |
|
39 /// It does the same as the isfinite() function defined by C99. |
|
40 template <typename T> |
|
41 bool isFinite(T value) |
|
42 { |
|
43 typedef std::numeric_limits<T> Lim; |
|
44 if ((Lim::has_infinity && (value == Lim::infinity() || value == |
|
45 -Lim::infinity())) || |
|
46 ((Lim::has_quiet_NaN || Lim::has_signaling_NaN) && value != value)) |
|
47 { |
|
48 return false; |
|
49 } |
|
50 return true; |
|
51 } |
|
52 |
|
53 ///Common base class for LP solvers |
|
54 |
|
55 ///\todo Much more docs |
|
56 ///\ingroup lp_group |
|
57 class LpSolverBase { |
|
58 |
|
59 protected: |
|
60 |
|
61 _lp_bits::LpId rows; |
|
62 _lp_bits::LpId cols; |
|
63 |
|
64 public: |
|
65 |
|
66 ///Possible outcomes of an LP solving procedure |
|
67 enum SolveExitStatus { |
|
68 ///This means that the problem has been successfully solved: either |
|
69 ///an optimal solution has been found or infeasibility/unboundedness |
|
70 ///has been proved. |
|
71 SOLVED = 0, |
|
72 ///Any other case (including the case when some user specified |
|
73 ///limit has been exceeded) |
|
74 UNSOLVED = 1 |
|
75 }; |
|
76 |
|
77 ///\e |
|
78 enum SolutionStatus { |
|
79 ///Feasible solution hasn't been found (but may exist). |
|
80 |
|
81 ///\todo NOTFOUND might be a better name. |
|
82 /// |
|
83 UNDEFINED = 0, |
|
84 ///The problem has no feasible solution |
|
85 INFEASIBLE = 1, |
|
86 ///Feasible solution found |
|
87 FEASIBLE = 2, |
|
88 ///Optimal solution exists and found |
|
89 OPTIMAL = 3, |
|
90 ///The cost function is unbounded |
|
91 |
|
92 ///\todo Give a feasible solution and an infinite ray (and the |
|
93 ///corresponding bases) |
|
94 INFINITE = 4 |
|
95 }; |
|
96 |
|
97 ///\e The type of the investigated LP problem |
|
98 enum ProblemTypes { |
|
99 ///Primal-dual feasible |
|
100 PRIMAL_DUAL_FEASIBLE = 0, |
|
101 ///Primal feasible dual infeasible |
|
102 PRIMAL_FEASIBLE_DUAL_INFEASIBLE = 1, |
|
103 ///Primal infeasible dual feasible |
|
104 PRIMAL_INFEASIBLE_DUAL_FEASIBLE = 2, |
|
105 ///Primal-dual infeasible |
|
106 PRIMAL_DUAL_INFEASIBLE = 3, |
|
107 ///Could not determine so far |
|
108 UNKNOWN = 4 |
|
109 }; |
|
110 |
|
111 ///The floating point type used by the solver |
|
112 typedef double Value; |
|
113 ///The infinity constant |
|
114 static const Value INF; |
|
115 ///The not a number constant |
|
116 static const Value NaN; |
|
117 |
|
118 static inline bool isNaN(const Value& v) { return v!=v; } |
|
119 |
|
120 friend class Col; |
|
121 friend class ColIt; |
|
122 friend class Row; |
|
123 |
|
124 ///Refer to a column of the LP. |
|
125 |
|
126 ///This type is used to refer to a column of the LP. |
|
127 /// |
|
128 ///Its value remains valid and correct even after the addition or erase of |
|
129 ///other columns. |
|
130 /// |
|
131 ///\todo Document what can one do with a Col (INVALID, comparing, |
|
132 ///it is similar to Node/Edge) |
|
133 class Col { |
|
134 protected: |
|
135 int id; |
|
136 friend class LpSolverBase; |
|
137 friend class MipSolverBase; |
|
138 explicit Col(int _id) : id(_id) {} |
|
139 public: |
|
140 typedef Value ExprValue; |
|
141 typedef True LpSolverCol; |
|
142 Col() {} |
|
143 Col(const Invalid&) : id(-1) {} |
|
144 bool operator< (Col c) const {return id< c.id;} |
|
145 bool operator> (Col c) const {return id> c.id;} |
|
146 bool operator==(Col c) const {return id==c.id;} |
|
147 bool operator!=(Col c) const {return id!=c.id;} |
|
148 }; |
|
149 |
|
150 class ColIt : public Col { |
|
151 const LpSolverBase *_lp; |
|
152 public: |
|
153 ColIt() {} |
|
154 ColIt(const LpSolverBase &lp) : _lp(&lp) |
|
155 { |
|
156 _lp->cols.firstFix(id); |
|
157 } |
|
158 ColIt(const Invalid&) : Col(INVALID) {} |
|
159 ColIt &operator++() |
|
160 { |
|
161 _lp->cols.nextFix(id); |
|
162 return *this; |
|
163 } |
|
164 }; |
|
165 |
|
166 static int id(const Col& col) { return col.id; } |
|
167 |
|
168 |
|
169 ///Refer to a row of the LP. |
|
170 |
|
171 ///This type is used to refer to a row of the LP. |
|
172 /// |
|
173 ///Its value remains valid and correct even after the addition or erase of |
|
174 ///other rows. |
|
175 /// |
|
176 ///\todo Document what can one do with a Row (INVALID, comparing, |
|
177 ///it is similar to Node/Edge) |
|
178 class Row { |
|
179 protected: |
|
180 int id; |
|
181 friend class LpSolverBase; |
|
182 explicit Row(int _id) : id(_id) {} |
|
183 public: |
|
184 typedef Value ExprValue; |
|
185 typedef True LpSolverRow; |
|
186 Row() {} |
|
187 Row(const Invalid&) : id(-1) {} |
|
188 |
|
189 bool operator< (Row c) const {return id< c.id;} |
|
190 bool operator> (Row c) const {return id> c.id;} |
|
191 bool operator==(Row c) const {return id==c.id;} |
|
192 bool operator!=(Row c) const {return id!=c.id;} |
|
193 }; |
|
194 |
|
195 class RowIt : public Row { |
|
196 const LpSolverBase *_lp; |
|
197 public: |
|
198 RowIt() {} |
|
199 RowIt(const LpSolverBase &lp) : _lp(&lp) |
|
200 { |
|
201 _lp->rows.firstFix(id); |
|
202 } |
|
203 RowIt(const Invalid&) : Row(INVALID) {} |
|
204 RowIt &operator++() |
|
205 { |
|
206 _lp->rows.nextFix(id); |
|
207 return *this; |
|
208 } |
|
209 }; |
|
210 |
|
211 static int id(const Row& row) { return row.id; } |
|
212 |
|
213 protected: |
|
214 |
|
215 int _lpId(const Col& c) const { |
|
216 return cols.floatingId(id(c)); |
|
217 } |
|
218 |
|
219 int _lpId(const Row& r) const { |
|
220 return rows.floatingId(id(r)); |
|
221 } |
|
222 |
|
223 Col _item(int i, Col) const { |
|
224 return Col(cols.fixId(i)); |
|
225 } |
|
226 |
|
227 Row _item(int i, Row) const { |
|
228 return Row(rows.fixId(i)); |
|
229 } |
|
230 |
|
231 |
|
232 public: |
|
233 |
|
234 ///Linear expression of variables and a constant component |
|
235 |
|
236 ///This data structure stores a linear expression of the variables |
|
237 ///(\ref Col "Col"s) and also has a constant component. |
|
238 /// |
|
239 ///There are several ways to access and modify the contents of this |
|
240 ///container. |
|
241 ///- Its it fully compatible with \c std::map<Col,double>, so for expamle |
|
242 ///if \c e is an Expr and \c v and \c w are of type \ref Col, then you can |
|
243 ///read and modify the coefficients like |
|
244 ///these. |
|
245 ///\code |
|
246 ///e[v]=5; |
|
247 ///e[v]+=12; |
|
248 ///e.erase(v); |
|
249 ///\endcode |
|
250 ///or you can also iterate through its elements. |
|
251 ///\code |
|
252 ///double s=0; |
|
253 ///for(LpSolverBase::Expr::iterator i=e.begin();i!=e.end();++i) |
|
254 /// s+=i->second; |
|
255 ///\endcode |
|
256 ///(This code computes the sum of all coefficients). |
|
257 ///- Numbers (<tt>double</tt>'s) |
|
258 ///and variables (\ref Col "Col"s) directly convert to an |
|
259 ///\ref Expr and the usual linear operations are defined, so |
|
260 ///\code |
|
261 ///v+w |
|
262 ///2*v-3.12*(v-w/2)+2 |
|
263 ///v*2.1+(3*v+(v*12+w+6)*3)/2 |
|
264 ///\endcode |
|
265 ///are valid \ref Expr "Expr"essions. |
|
266 ///The usual assignment operations are also defined. |
|
267 ///\code |
|
268 ///e=v+w; |
|
269 ///e+=2*v-3.12*(v-w/2)+2; |
|
270 ///e*=3.4; |
|
271 ///e/=5; |
|
272 ///\endcode |
|
273 ///- The constant member can be set and read by \ref constComp() |
|
274 ///\code |
|
275 ///e.constComp()=12; |
|
276 ///double c=e.constComp(); |
|
277 ///\endcode |
|
278 /// |
|
279 ///\note \ref clear() not only sets all coefficients to 0 but also |
|
280 ///clears the constant components. |
|
281 /// |
|
282 ///\sa Constr |
|
283 /// |
|
284 class Expr : public std::map<Col,Value> |
|
285 { |
|
286 public: |
|
287 typedef LpSolverBase::Col Key; |
|
288 typedef LpSolverBase::Value Value; |
|
289 |
|
290 protected: |
|
291 typedef std::map<Col,Value> Base; |
|
292 |
|
293 Value const_comp; |
|
294 public: |
|
295 typedef True IsLinExpression; |
|
296 ///\e |
|
297 Expr() : Base(), const_comp(0) { } |
|
298 ///\e |
|
299 Expr(const Key &v) : const_comp(0) { |
|
300 Base::insert(std::make_pair(v, 1)); |
|
301 } |
|
302 ///\e |
|
303 Expr(const Value &v) : const_comp(v) {} |
|
304 ///\e |
|
305 void set(const Key &v,const Value &c) { |
|
306 Base::insert(std::make_pair(v, c)); |
|
307 } |
|
308 ///\e |
|
309 Value &constComp() { return const_comp; } |
|
310 ///\e |
|
311 const Value &constComp() const { return const_comp; } |
|
312 |
|
313 ///Removes the components with zero coefficient. |
|
314 void simplify() { |
|
315 for (Base::iterator i=Base::begin(); i!=Base::end();) { |
|
316 Base::iterator j=i; |
|
317 ++j; |
|
318 if ((*i).second==0) Base::erase(i); |
|
319 i=j; |
|
320 } |
|
321 } |
|
322 |
|
323 void simplify() const { |
|
324 const_cast<Expr*>(this)->simplify(); |
|
325 } |
|
326 |
|
327 ///Removes the coefficients closer to zero than \c tolerance. |
|
328 void simplify(double &tolerance) { |
|
329 for (Base::iterator i=Base::begin(); i!=Base::end();) { |
|
330 Base::iterator j=i; |
|
331 ++j; |
|
332 if (std::fabs((*i).second)<tolerance) Base::erase(i); |
|
333 i=j; |
|
334 } |
|
335 } |
|
336 |
|
337 ///Sets all coefficients and the constant component to 0. |
|
338 void clear() { |
|
339 Base::clear(); |
|
340 const_comp=0; |
|
341 } |
|
342 |
|
343 ///\e |
|
344 Expr &operator+=(const Expr &e) { |
|
345 for (Base::const_iterator j=e.begin(); j!=e.end(); ++j) |
|
346 (*this)[j->first]+=j->second; |
|
347 const_comp+=e.const_comp; |
|
348 return *this; |
|
349 } |
|
350 ///\e |
|
351 Expr &operator-=(const Expr &e) { |
|
352 for (Base::const_iterator j=e.begin(); j!=e.end(); ++j) |
|
353 (*this)[j->first]-=j->second; |
|
354 const_comp-=e.const_comp; |
|
355 return *this; |
|
356 } |
|
357 ///\e |
|
358 Expr &operator*=(const Value &c) { |
|
359 for (Base::iterator j=Base::begin(); j!=Base::end(); ++j) |
|
360 j->second*=c; |
|
361 const_comp*=c; |
|
362 return *this; |
|
363 } |
|
364 ///\e |
|
365 Expr &operator/=(const Value &c) { |
|
366 for (Base::iterator j=Base::begin(); j!=Base::end(); ++j) |
|
367 j->second/=c; |
|
368 const_comp/=c; |
|
369 return *this; |
|
370 } |
|
371 |
|
372 }; |
|
373 |
|
374 ///Linear constraint |
|
375 |
|
376 ///This data stucture represents a linear constraint in the LP. |
|
377 ///Basically it is a linear expression with a lower or an upper bound |
|
378 ///(or both). These parts of the constraint can be obtained by the member |
|
379 ///functions \ref expr(), \ref lowerBound() and \ref upperBound(), |
|
380 ///respectively. |
|
381 ///There are two ways to construct a constraint. |
|
382 ///- You can set the linear expression and the bounds directly |
|
383 /// by the functions above. |
|
384 ///- The operators <tt>\<=</tt>, <tt>==</tt> and <tt>\>=</tt> |
|
385 /// are defined between expressions, or even between constraints whenever |
|
386 /// it makes sense. Therefore if \c e and \c f are linear expressions and |
|
387 /// \c s and \c t are numbers, then the followings are valid expressions |
|
388 /// and thus they can be used directly e.g. in \ref addRow() whenever |
|
389 /// it makes sense. |
|
390 ///\code |
|
391 /// e<=s |
|
392 /// e<=f |
|
393 /// e==f |
|
394 /// s<=e<=t |
|
395 /// e>=t |
|
396 ///\endcode |
|
397 ///\warning The validity of a constraint is checked only at run time, so |
|
398 ///e.g. \ref addRow(<tt>x[1]\<=x[2]<=5</tt>) will compile, but will throw |
|
399 ///an assertion. |
|
400 class Constr |
|
401 { |
|
402 public: |
|
403 typedef LpSolverBase::Expr Expr; |
|
404 typedef Expr::Key Key; |
|
405 typedef Expr::Value Value; |
|
406 |
|
407 protected: |
|
408 Expr _expr; |
|
409 Value _lb,_ub; |
|
410 public: |
|
411 ///\e |
|
412 Constr() : _expr(), _lb(NaN), _ub(NaN) {} |
|
413 ///\e |
|
414 Constr(Value lb,const Expr &e,Value ub) : |
|
415 _expr(e), _lb(lb), _ub(ub) {} |
|
416 ///\e |
|
417 Constr(const Expr &e,Value ub) : |
|
418 _expr(e), _lb(NaN), _ub(ub) {} |
|
419 ///\e |
|
420 Constr(Value lb,const Expr &e) : |
|
421 _expr(e), _lb(lb), _ub(NaN) {} |
|
422 ///\e |
|
423 Constr(const Expr &e) : |
|
424 _expr(e), _lb(NaN), _ub(NaN) {} |
|
425 ///\e |
|
426 void clear() |
|
427 { |
|
428 _expr.clear(); |
|
429 _lb=_ub=NaN; |
|
430 } |
|
431 |
|
432 ///Reference to the linear expression |
|
433 Expr &expr() { return _expr; } |
|
434 ///Cont reference to the linear expression |
|
435 const Expr &expr() const { return _expr; } |
|
436 ///Reference to the lower bound. |
|
437 |
|
438 ///\return |
|
439 ///- \ref INF "INF": the constraint is lower unbounded. |
|
440 ///- \ref NaN "NaN": lower bound has not been set. |
|
441 ///- finite number: the lower bound |
|
442 Value &lowerBound() { return _lb; } |
|
443 ///The const version of \ref lowerBound() |
|
444 const Value &lowerBound() const { return _lb; } |
|
445 ///Reference to the upper bound. |
|
446 |
|
447 ///\return |
|
448 ///- \ref INF "INF": the constraint is upper unbounded. |
|
449 ///- \ref NaN "NaN": upper bound has not been set. |
|
450 ///- finite number: the upper bound |
|
451 Value &upperBound() { return _ub; } |
|
452 ///The const version of \ref upperBound() |
|
453 const Value &upperBound() const { return _ub; } |
|
454 ///Is the constraint lower bounded? |
|
455 bool lowerBounded() const { |
|
456 return isFinite(_lb); |
|
457 } |
|
458 ///Is the constraint upper bounded? |
|
459 bool upperBounded() const { |
|
460 return isFinite(_ub); |
|
461 } |
|
462 |
|
463 }; |
|
464 |
|
465 ///Linear expression of rows |
|
466 |
|
467 ///This data structure represents a column of the matrix, |
|
468 ///thas is it strores a linear expression of the dual variables |
|
469 ///(\ref Row "Row"s). |
|
470 /// |
|
471 ///There are several ways to access and modify the contents of this |
|
472 ///container. |
|
473 ///- Its it fully compatible with \c std::map<Row,double>, so for expamle |
|
474 ///if \c e is an DualExpr and \c v |
|
475 ///and \c w are of type \ref Row, then you can |
|
476 ///read and modify the coefficients like |
|
477 ///these. |
|
478 ///\code |
|
479 ///e[v]=5; |
|
480 ///e[v]+=12; |
|
481 ///e.erase(v); |
|
482 ///\endcode |
|
483 ///or you can also iterate through its elements. |
|
484 ///\code |
|
485 ///double s=0; |
|
486 ///for(LpSolverBase::DualExpr::iterator i=e.begin();i!=e.end();++i) |
|
487 /// s+=i->second; |
|
488 ///\endcode |
|
489 ///(This code computes the sum of all coefficients). |
|
490 ///- Numbers (<tt>double</tt>'s) |
|
491 ///and variables (\ref Row "Row"s) directly convert to an |
|
492 ///\ref DualExpr and the usual linear operations are defined, so |
|
493 ///\code |
|
494 ///v+w |
|
495 ///2*v-3.12*(v-w/2) |
|
496 ///v*2.1+(3*v+(v*12+w)*3)/2 |
|
497 ///\endcode |
|
498 ///are valid \ref DualExpr "DualExpr"essions. |
|
499 ///The usual assignment operations are also defined. |
|
500 ///\code |
|
501 ///e=v+w; |
|
502 ///e+=2*v-3.12*(v-w/2); |
|
503 ///e*=3.4; |
|
504 ///e/=5; |
|
505 ///\endcode |
|
506 /// |
|
507 ///\sa Expr |
|
508 /// |
|
509 class DualExpr : public std::map<Row,Value> |
|
510 { |
|
511 public: |
|
512 typedef LpSolverBase::Row Key; |
|
513 typedef LpSolverBase::Value Value; |
|
514 |
|
515 protected: |
|
516 typedef std::map<Row,Value> Base; |
|
517 |
|
518 public: |
|
519 typedef True IsLinExpression; |
|
520 ///\e |
|
521 DualExpr() : Base() { } |
|
522 ///\e |
|
523 DualExpr(const Key &v) { |
|
524 Base::insert(std::make_pair(v, 1)); |
|
525 } |
|
526 ///\e |
|
527 void set(const Key &v,const Value &c) { |
|
528 Base::insert(std::make_pair(v, c)); |
|
529 } |
|
530 |
|
531 ///Removes the components with zero coefficient. |
|
532 void simplify() { |
|
533 for (Base::iterator i=Base::begin(); i!=Base::end();) { |
|
534 Base::iterator j=i; |
|
535 ++j; |
|
536 if ((*i).second==0) Base::erase(i); |
|
537 i=j; |
|
538 } |
|
539 } |
|
540 |
|
541 void simplify() const { |
|
542 const_cast<DualExpr*>(this)->simplify(); |
|
543 } |
|
544 |
|
545 ///Removes the coefficients closer to zero than \c tolerance. |
|
546 void simplify(double &tolerance) { |
|
547 for (Base::iterator i=Base::begin(); i!=Base::end();) { |
|
548 Base::iterator j=i; |
|
549 ++j; |
|
550 if (std::fabs((*i).second)<tolerance) Base::erase(i); |
|
551 i=j; |
|
552 } |
|
553 } |
|
554 |
|
555 ///Sets all coefficients to 0. |
|
556 void clear() { |
|
557 Base::clear(); |
|
558 } |
|
559 |
|
560 ///\e |
|
561 DualExpr &operator+=(const DualExpr &e) { |
|
562 for (Base::const_iterator j=e.begin(); j!=e.end(); ++j) |
|
563 (*this)[j->first]+=j->second; |
|
564 return *this; |
|
565 } |
|
566 ///\e |
|
567 DualExpr &operator-=(const DualExpr &e) { |
|
568 for (Base::const_iterator j=e.begin(); j!=e.end(); ++j) |
|
569 (*this)[j->first]-=j->second; |
|
570 return *this; |
|
571 } |
|
572 ///\e |
|
573 DualExpr &operator*=(const Value &c) { |
|
574 for (Base::iterator j=Base::begin(); j!=Base::end(); ++j) |
|
575 j->second*=c; |
|
576 return *this; |
|
577 } |
|
578 ///\e |
|
579 DualExpr &operator/=(const Value &c) { |
|
580 for (Base::iterator j=Base::begin(); j!=Base::end(); ++j) |
|
581 j->second/=c; |
|
582 return *this; |
|
583 } |
|
584 }; |
|
585 |
|
586 |
|
587 private: |
|
588 |
|
589 template <typename _Expr> |
|
590 class MappedOutputIterator { |
|
591 public: |
|
592 |
|
593 typedef std::insert_iterator<_Expr> Base; |
|
594 |
|
595 typedef std::output_iterator_tag iterator_category; |
|
596 typedef void difference_type; |
|
597 typedef void value_type; |
|
598 typedef void reference; |
|
599 typedef void pointer; |
|
600 |
|
601 MappedOutputIterator(const Base& _base, const LpSolverBase& _lp) |
|
602 : base(_base), lp(_lp) {} |
|
603 |
|
604 MappedOutputIterator& operator*() { |
|
605 return *this; |
|
606 } |
|
607 |
|
608 MappedOutputIterator& operator=(const std::pair<int, Value>& value) { |
|
609 *base = std::make_pair(lp._item(value.first, typename _Expr::Key()), |
|
610 value.second); |
|
611 return *this; |
|
612 } |
|
613 |
|
614 MappedOutputIterator& operator++() { |
|
615 ++base; |
|
616 return *this; |
|
617 } |
|
618 |
|
619 MappedOutputIterator operator++(int) { |
|
620 MappedOutputIterator tmp(*this); |
|
621 ++base; |
|
622 return tmp; |
|
623 } |
|
624 |
|
625 bool operator==(const MappedOutputIterator& it) const { |
|
626 return base == it.base; |
|
627 } |
|
628 |
|
629 bool operator!=(const MappedOutputIterator& it) const { |
|
630 return base != it.base; |
|
631 } |
|
632 |
|
633 private: |
|
634 Base base; |
|
635 const LpSolverBase& lp; |
|
636 }; |
|
637 |
|
638 template <typename Expr> |
|
639 class MappedInputIterator { |
|
640 public: |
|
641 |
|
642 typedef typename Expr::const_iterator Base; |
|
643 |
|
644 typedef typename Base::iterator_category iterator_category; |
|
645 typedef typename Base::difference_type difference_type; |
|
646 typedef const std::pair<int, Value> value_type; |
|
647 typedef value_type reference; |
|
648 class pointer { |
|
649 public: |
|
650 pointer(value_type& _value) : value(_value) {} |
|
651 value_type* operator->() { return &value; } |
|
652 private: |
|
653 value_type value; |
|
654 }; |
|
655 |
|
656 MappedInputIterator(const Base& _base, const LpSolverBase& _lp) |
|
657 : base(_base), lp(_lp) {} |
|
658 |
|
659 reference operator*() { |
|
660 return std::make_pair(lp._lpId(base->first), base->second); |
|
661 } |
|
662 |
|
663 pointer operator->() { |
|
664 return pointer(operator*()); |
|
665 } |
|
666 |
|
667 MappedInputIterator& operator++() { |
|
668 ++base; |
|
669 return *this; |
|
670 } |
|
671 |
|
672 MappedInputIterator operator++(int) { |
|
673 MappedInputIterator tmp(*this); |
|
674 ++base; |
|
675 return tmp; |
|
676 } |
|
677 |
|
678 bool operator==(const MappedInputIterator& it) const { |
|
679 return base == it.base; |
|
680 } |
|
681 |
|
682 bool operator!=(const MappedInputIterator& it) const { |
|
683 return base != it.base; |
|
684 } |
|
685 |
|
686 private: |
|
687 Base base; |
|
688 const LpSolverBase& lp; |
|
689 }; |
|
690 |
|
691 protected: |
|
692 |
|
693 /// STL compatible iterator for lp col |
|
694 typedef MappedInputIterator<Expr> ConstRowIterator; |
|
695 /// STL compatible iterator for lp row |
|
696 typedef MappedInputIterator<DualExpr> ConstColIterator; |
|
697 |
|
698 /// STL compatible iterator for lp col |
|
699 typedef MappedOutputIterator<Expr> RowIterator; |
|
700 /// STL compatible iterator for lp row |
|
701 typedef MappedOutputIterator<DualExpr> ColIterator; |
|
702 |
|
703 //Abstract virtual functions |
|
704 virtual LpSolverBase* _newLp() = 0; |
|
705 virtual LpSolverBase* _copyLp(){ |
|
706 LpSolverBase* newlp = _newLp(); |
|
707 |
|
708 std::map<Col, Col> ref; |
|
709 for (LpSolverBase::ColIt it(*this); it != INVALID; ++it) { |
|
710 Col ccol = newlp->addCol(); |
|
711 ref[it] = ccol; |
|
712 newlp->colName(ccol, colName(it)); |
|
713 newlp->colLowerBound(ccol, colLowerBound(it)); |
|
714 newlp->colUpperBound(ccol, colUpperBound(it)); |
|
715 } |
|
716 |
|
717 for (LpSolverBase::RowIt it(*this); it != INVALID; ++it) { |
|
718 Expr e = row(it), ce; |
|
719 for (Expr::iterator jt = e.begin(); jt != e.end(); ++jt) { |
|
720 ce[ref[jt->first]] = jt->second; |
|
721 } |
|
722 ce += e.constComp(); |
|
723 Row r = newlp->addRow(ce); |
|
724 |
|
725 double lower, upper; |
|
726 getRowBounds(it, lower, upper); |
|
727 newlp->rowBounds(r, lower, upper); |
|
728 } |
|
729 |
|
730 return newlp; |
|
731 }; |
|
732 |
|
733 virtual int _addCol() = 0; |
|
734 virtual int _addRow() = 0; |
|
735 |
|
736 virtual void _eraseCol(int col) = 0; |
|
737 virtual void _eraseRow(int row) = 0; |
|
738 |
|
739 virtual void _getColName(int col, std::string & name) const = 0; |
|
740 virtual void _setColName(int col, const std::string & name) = 0; |
|
741 virtual int _colByName(const std::string& name) const = 0; |
|
742 |
|
743 virtual void _setRowCoeffs(int i, ConstRowIterator b, |
|
744 ConstRowIterator e) = 0; |
|
745 virtual void _getRowCoeffs(int i, RowIterator b) const = 0; |
|
746 virtual void _setColCoeffs(int i, ConstColIterator b, |
|
747 ConstColIterator e) = 0; |
|
748 virtual void _getColCoeffs(int i, ColIterator b) const = 0; |
|
749 virtual void _setCoeff(int row, int col, Value value) = 0; |
|
750 virtual Value _getCoeff(int row, int col) const = 0; |
|
751 virtual void _setColLowerBound(int i, Value value) = 0; |
|
752 virtual Value _getColLowerBound(int i) const = 0; |
|
753 virtual void _setColUpperBound(int i, Value value) = 0; |
|
754 virtual Value _getColUpperBound(int i) const = 0; |
|
755 virtual void _setRowBounds(int i, Value lower, Value upper) = 0; |
|
756 virtual void _getRowBounds(int i, Value &lower, Value &upper) const = 0; |
|
757 |
|
758 virtual void _setObjCoeff(int i, Value obj_coef) = 0; |
|
759 virtual Value _getObjCoeff(int i) const = 0; |
|
760 virtual void _clearObj()=0; |
|
761 |
|
762 virtual SolveExitStatus _solve() = 0; |
|
763 virtual Value _getPrimal(int i) const = 0; |
|
764 virtual Value _getDual(int i) const = 0; |
|
765 virtual Value _getPrimalValue() const = 0; |
|
766 virtual bool _isBasicCol(int i) const = 0; |
|
767 virtual SolutionStatus _getPrimalStatus() const = 0; |
|
768 virtual SolutionStatus _getDualStatus() const = 0; |
|
769 virtual ProblemTypes _getProblemType() const = 0; |
|
770 |
|
771 virtual void _setMax() = 0; |
|
772 virtual void _setMin() = 0; |
|
773 |
|
774 |
|
775 virtual bool _isMax() const = 0; |
|
776 |
|
777 //Own protected stuff |
|
778 |
|
779 //Constant component of the objective function |
|
780 Value obj_const_comp; |
|
781 |
|
782 public: |
|
783 |
|
784 ///\e |
|
785 LpSolverBase() : obj_const_comp(0) {} |
|
786 |
|
787 ///\e |
|
788 virtual ~LpSolverBase() {} |
|
789 |
|
790 ///Creates a new LP problem |
|
791 LpSolverBase* newLp() {return _newLp();} |
|
792 ///Makes a copy of the LP problem |
|
793 LpSolverBase* copyLp() {return _copyLp();} |
|
794 |
|
795 ///\name Build up and modify the LP |
|
796 |
|
797 ///@{ |
|
798 |
|
799 ///Add a new empty column (i.e a new variable) to the LP |
|
800 Col addCol() { Col c; _addCol(); c.id = cols.addId(); return c;} |
|
801 |
|
802 ///\brief Adds several new columns |
|
803 ///(i.e a variables) at once |
|
804 /// |
|
805 ///This magic function takes a container as its argument |
|
806 ///and fills its elements |
|
807 ///with new columns (i.e. variables) |
|
808 ///\param t can be |
|
809 ///- a standard STL compatible iterable container with |
|
810 ///\ref Col as its \c values_type |
|
811 ///like |
|
812 ///\code |
|
813 ///std::vector<LpSolverBase::Col> |
|
814 ///std::list<LpSolverBase::Col> |
|
815 ///\endcode |
|
816 ///- a standard STL compatible iterable container with |
|
817 ///\ref Col as its \c mapped_type |
|
818 ///like |
|
819 ///\code |
|
820 ///std::map<AnyType,LpSolverBase::Col> |
|
821 ///\endcode |
|
822 ///- an iterable lemon \ref concepts::WriteMap "write map" like |
|
823 ///\code |
|
824 ///ListGraph::NodeMap<LpSolverBase::Col> |
|
825 ///ListGraph::EdgeMap<LpSolverBase::Col> |
|
826 ///\endcode |
|
827 ///\return The number of the created column. |
|
828 #ifdef DOXYGEN |
|
829 template<class T> |
|
830 int addColSet(T &t) { return 0;} |
|
831 #else |
|
832 template<class T> |
|
833 typename enable_if<typename T::value_type::LpSolverCol,int>::type |
|
834 addColSet(T &t,dummy<0> = 0) { |
|
835 int s=0; |
|
836 for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addCol();s++;} |
|
837 return s; |
|
838 } |
|
839 template<class T> |
|
840 typename enable_if<typename T::value_type::second_type::LpSolverCol, |
|
841 int>::type |
|
842 addColSet(T &t,dummy<1> = 1) { |
|
843 int s=0; |
|
844 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
845 i->second=addCol(); |
|
846 s++; |
|
847 } |
|
848 return s; |
|
849 } |
|
850 template<class T> |
|
851 typename enable_if<typename T::MapIt::Value::LpSolverCol, |
|
852 int>::type |
|
853 addColSet(T &t,dummy<2> = 2) { |
|
854 int s=0; |
|
855 for(typename T::MapIt i(t); i!=INVALID; ++i) |
|
856 { |
|
857 i.set(addCol()); |
|
858 s++; |
|
859 } |
|
860 return s; |
|
861 } |
|
862 #endif |
|
863 |
|
864 ///Set a column (i.e a dual constraint) of the LP |
|
865 |
|
866 ///\param c is the column to be modified |
|
867 ///\param e is a dual linear expression (see \ref DualExpr) |
|
868 ///a better one. |
|
869 void col(Col c,const DualExpr &e) { |
|
870 e.simplify(); |
|
871 _setColCoeffs(_lpId(c), ConstColIterator(e.begin(), *this), |
|
872 ConstColIterator(e.end(), *this)); |
|
873 } |
|
874 |
|
875 ///Get a column (i.e a dual constraint) of the LP |
|
876 |
|
877 ///\param r is the column to get |
|
878 ///\return the dual expression associated to the column |
|
879 DualExpr col(Col c) const { |
|
880 DualExpr e; |
|
881 _getColCoeffs(_lpId(c), ColIterator(std::inserter(e, e.end()), *this)); |
|
882 return e; |
|
883 } |
|
884 |
|
885 ///Add a new column to the LP |
|
886 |
|
887 ///\param e is a dual linear expression (see \ref DualExpr) |
|
888 ///\param obj is the corresponding component of the objective |
|
889 ///function. It is 0 by default. |
|
890 ///\return The created column. |
|
891 Col addCol(const DualExpr &e, Value o = 0) { |
|
892 Col c=addCol(); |
|
893 col(c,e); |
|
894 objCoeff(c,o); |
|
895 return c; |
|
896 } |
|
897 |
|
898 ///Add a new empty row (i.e a new constraint) to the LP |
|
899 |
|
900 ///This function adds a new empty row (i.e a new constraint) to the LP. |
|
901 ///\return The created row |
|
902 Row addRow() { Row r; _addRow(); r.id = rows.addId(); return r;} |
|
903 |
|
904 ///\brief Add several new rows |
|
905 ///(i.e a constraints) at once |
|
906 /// |
|
907 ///This magic function takes a container as its argument |
|
908 ///and fills its elements |
|
909 ///with new row (i.e. variables) |
|
910 ///\param t can be |
|
911 ///- a standard STL compatible iterable container with |
|
912 ///\ref Row as its \c values_type |
|
913 ///like |
|
914 ///\code |
|
915 ///std::vector<LpSolverBase::Row> |
|
916 ///std::list<LpSolverBase::Row> |
|
917 ///\endcode |
|
918 ///- a standard STL compatible iterable container with |
|
919 ///\ref Row as its \c mapped_type |
|
920 ///like |
|
921 ///\code |
|
922 ///std::map<AnyType,LpSolverBase::Row> |
|
923 ///\endcode |
|
924 ///- an iterable lemon \ref concepts::WriteMap "write map" like |
|
925 ///\code |
|
926 ///ListGraph::NodeMap<LpSolverBase::Row> |
|
927 ///ListGraph::EdgeMap<LpSolverBase::Row> |
|
928 ///\endcode |
|
929 ///\return The number of rows created. |
|
930 #ifdef DOXYGEN |
|
931 template<class T> |
|
932 int addRowSet(T &t) { return 0;} |
|
933 #else |
|
934 template<class T> |
|
935 typename enable_if<typename T::value_type::LpSolverRow,int>::type |
|
936 addRowSet(T &t,dummy<0> = 0) { |
|
937 int s=0; |
|
938 for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addRow();s++;} |
|
939 return s; |
|
940 } |
|
941 template<class T> |
|
942 typename enable_if<typename T::value_type::second_type::LpSolverRow, |
|
943 int>::type |
|
944 addRowSet(T &t,dummy<1> = 1) { |
|
945 int s=0; |
|
946 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
947 i->second=addRow(); |
|
948 s++; |
|
949 } |
|
950 return s; |
|
951 } |
|
952 template<class T> |
|
953 typename enable_if<typename T::MapIt::Value::LpSolverRow, |
|
954 int>::type |
|
955 addRowSet(T &t,dummy<2> = 2) { |
|
956 int s=0; |
|
957 for(typename T::MapIt i(t); i!=INVALID; ++i) |
|
958 { |
|
959 i.set(addRow()); |
|
960 s++; |
|
961 } |
|
962 return s; |
|
963 } |
|
964 #endif |
|
965 |
|
966 ///Set a row (i.e a constraint) of the LP |
|
967 |
|
968 ///\param r is the row to be modified |
|
969 ///\param l is lower bound (-\ref INF means no bound) |
|
970 ///\param e is a linear expression (see \ref Expr) |
|
971 ///\param u is the upper bound (\ref INF means no bound) |
|
972 ///\bug This is a temporary function. The interface will change to |
|
973 ///a better one. |
|
974 ///\todo Option to control whether a constraint with a single variable is |
|
975 ///added or not. |
|
976 void row(Row r, Value l, const Expr &e, Value u) { |
|
977 e.simplify(); |
|
978 _setRowCoeffs(_lpId(r), ConstRowIterator(e.begin(), *this), |
|
979 ConstRowIterator(e.end(), *this)); |
|
980 _setRowBounds(_lpId(r),l-e.constComp(),u-e.constComp()); |
|
981 } |
|
982 |
|
983 ///Set a row (i.e a constraint) of the LP |
|
984 |
|
985 ///\param r is the row to be modified |
|
986 ///\param c is a linear expression (see \ref Constr) |
|
987 void row(Row r, const Constr &c) { |
|
988 row(r, c.lowerBounded()?c.lowerBound():-INF, |
|
989 c.expr(), c.upperBounded()?c.upperBound():INF); |
|
990 } |
|
991 |
|
992 |
|
993 ///Get a row (i.e a constraint) of the LP |
|
994 |
|
995 ///\param r is the row to get |
|
996 ///\return the expression associated to the row |
|
997 Expr row(Row r) const { |
|
998 Expr e; |
|
999 _getRowCoeffs(_lpId(r), RowIterator(std::inserter(e, e.end()), *this)); |
|
1000 return e; |
|
1001 } |
|
1002 |
|
1003 ///Add a new row (i.e a new constraint) to the LP |
|
1004 |
|
1005 ///\param l is the lower bound (-\ref INF means no bound) |
|
1006 ///\param e is a linear expression (see \ref Expr) |
|
1007 ///\param u is the upper bound (\ref INF means no bound) |
|
1008 ///\return The created row. |
|
1009 ///\bug This is a temporary function. The interface will change to |
|
1010 ///a better one. |
|
1011 Row addRow(Value l,const Expr &e, Value u) { |
|
1012 Row r=addRow(); |
|
1013 row(r,l,e,u); |
|
1014 return r; |
|
1015 } |
|
1016 |
|
1017 ///Add a new row (i.e a new constraint) to the LP |
|
1018 |
|
1019 ///\param c is a linear expression (see \ref Constr) |
|
1020 ///\return The created row. |
|
1021 Row addRow(const Constr &c) { |
|
1022 Row r=addRow(); |
|
1023 row(r,c); |
|
1024 return r; |
|
1025 } |
|
1026 ///Erase a coloumn (i.e a variable) from the LP |
|
1027 |
|
1028 ///\param c is the coloumn to be deleted |
|
1029 ///\todo Please check this |
|
1030 void eraseCol(Col c) { |
|
1031 _eraseCol(_lpId(c)); |
|
1032 cols.eraseId(c.id); |
|
1033 } |
|
1034 ///Erase a row (i.e a constraint) from the LP |
|
1035 |
|
1036 ///\param r is the row to be deleted |
|
1037 ///\todo Please check this |
|
1038 void eraseRow(Row r) { |
|
1039 _eraseRow(_lpId(r)); |
|
1040 rows.eraseId(r.id); |
|
1041 } |
|
1042 |
|
1043 /// Get the name of a column |
|
1044 |
|
1045 ///\param c is the coresponding coloumn |
|
1046 ///\return The name of the colunm |
|
1047 std::string colName(Col c) const { |
|
1048 std::string name; |
|
1049 _getColName(_lpId(c), name); |
|
1050 return name; |
|
1051 } |
|
1052 |
|
1053 /// Set the name of a column |
|
1054 |
|
1055 ///\param c is the coresponding coloumn |
|
1056 ///\param name The name to be given |
|
1057 void colName(Col c, const std::string& name) { |
|
1058 _setColName(_lpId(c), name); |
|
1059 } |
|
1060 |
|
1061 /// Get the column by its name |
|
1062 |
|
1063 ///\param name The name of the column |
|
1064 ///\return the proper column or \c INVALID |
|
1065 Col colByName(const std::string& name) const { |
|
1066 int k = _colByName(name); |
|
1067 return k != -1 ? Col(cols.fixId(k)) : Col(INVALID); |
|
1068 } |
|
1069 |
|
1070 /// Set an element of the coefficient matrix of the LP |
|
1071 |
|
1072 ///\param r is the row of the element to be modified |
|
1073 ///\param c is the coloumn of the element to be modified |
|
1074 ///\param val is the new value of the coefficient |
|
1075 |
|
1076 void coeff(Row r, Col c, Value val) { |
|
1077 _setCoeff(_lpId(r),_lpId(c), val); |
|
1078 } |
|
1079 |
|
1080 /// Get an element of the coefficient matrix of the LP |
|
1081 |
|
1082 ///\param r is the row of the element in question |
|
1083 ///\param c is the coloumn of the element in question |
|
1084 ///\return the corresponding coefficient |
|
1085 |
|
1086 Value coeff(Row r, Col c) const { |
|
1087 return _getCoeff(_lpId(r),_lpId(c)); |
|
1088 } |
|
1089 |
|
1090 /// Set the lower bound of a column (i.e a variable) |
|
1091 |
|
1092 /// The lower bound of a variable (column) has to be given by an |
|
1093 /// extended number of type Value, i.e. a finite number of type |
|
1094 /// Value or -\ref INF. |
|
1095 void colLowerBound(Col c, Value value) { |
|
1096 _setColLowerBound(_lpId(c),value); |
|
1097 } |
|
1098 |
|
1099 /// Get the lower bound of a column (i.e a variable) |
|
1100 |
|
1101 /// This function returns the lower bound for column (variable) \t c |
|
1102 /// (this might be -\ref INF as well). |
|
1103 ///\return The lower bound for coloumn \t c |
|
1104 Value colLowerBound(Col c) const { |
|
1105 return _getColLowerBound(_lpId(c)); |
|
1106 } |
|
1107 |
|
1108 ///\brief Set the lower bound of several columns |
|
1109 ///(i.e a variables) at once |
|
1110 /// |
|
1111 ///This magic function takes a container as its argument |
|
1112 ///and applies the function on all of its elements. |
|
1113 /// The lower bound of a variable (column) has to be given by an |
|
1114 /// extended number of type Value, i.e. a finite number of type |
|
1115 /// Value or -\ref INF. |
|
1116 #ifdef DOXYGEN |
|
1117 template<class T> |
|
1118 void colLowerBound(T &t, Value value) { return 0;} |
|
1119 #else |
|
1120 template<class T> |
|
1121 typename enable_if<typename T::value_type::LpSolverCol,void>::type |
|
1122 colLowerBound(T &t, Value value,dummy<0> = 0) { |
|
1123 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
1124 colLowerBound(*i, value); |
|
1125 } |
|
1126 } |
|
1127 template<class T> |
|
1128 typename enable_if<typename T::value_type::second_type::LpSolverCol, |
|
1129 void>::type |
|
1130 colLowerBound(T &t, Value value,dummy<1> = 1) { |
|
1131 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
1132 colLowerBound(i->second, value); |
|
1133 } |
|
1134 } |
|
1135 template<class T> |
|
1136 typename enable_if<typename T::MapIt::Value::LpSolverCol, |
|
1137 void>::type |
|
1138 colLowerBound(T &t, Value value,dummy<2> = 2) { |
|
1139 for(typename T::MapIt i(t); i!=INVALID; ++i){ |
|
1140 colLowerBound(*i, value); |
|
1141 } |
|
1142 } |
|
1143 #endif |
|
1144 |
|
1145 /// Set the upper bound of a column (i.e a variable) |
|
1146 |
|
1147 /// The upper bound of a variable (column) has to be given by an |
|
1148 /// extended number of type Value, i.e. a finite number of type |
|
1149 /// Value or \ref INF. |
|
1150 void colUpperBound(Col c, Value value) { |
|
1151 _setColUpperBound(_lpId(c),value); |
|
1152 }; |
|
1153 |
|
1154 /// Get the upper bound of a column (i.e a variable) |
|
1155 |
|
1156 /// This function returns the upper bound for column (variable) \t c |
|
1157 /// (this might be \ref INF as well). |
|
1158 ///\return The upper bound for coloumn \t c |
|
1159 Value colUpperBound(Col c) const { |
|
1160 return _getColUpperBound(_lpId(c)); |
|
1161 } |
|
1162 |
|
1163 ///\brief Set the upper bound of several columns |
|
1164 ///(i.e a variables) at once |
|
1165 /// |
|
1166 ///This magic function takes a container as its argument |
|
1167 ///and applies the function on all of its elements. |
|
1168 /// The upper bound of a variable (column) has to be given by an |
|
1169 /// extended number of type Value, i.e. a finite number of type |
|
1170 /// Value or \ref INF. |
|
1171 #ifdef DOXYGEN |
|
1172 template<class T> |
|
1173 void colUpperBound(T &t, Value value) { return 0;} |
|
1174 #else |
|
1175 template<class T> |
|
1176 typename enable_if<typename T::value_type::LpSolverCol,void>::type |
|
1177 colUpperBound(T &t, Value value,dummy<0> = 0) { |
|
1178 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
1179 colUpperBound(*i, value); |
|
1180 } |
|
1181 } |
|
1182 template<class T> |
|
1183 typename enable_if<typename T::value_type::second_type::LpSolverCol, |
|
1184 void>::type |
|
1185 colUpperBound(T &t, Value value,dummy<1> = 1) { |
|
1186 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
1187 colUpperBound(i->second, value); |
|
1188 } |
|
1189 } |
|
1190 template<class T> |
|
1191 typename enable_if<typename T::MapIt::Value::LpSolverCol, |
|
1192 void>::type |
|
1193 colUpperBound(T &t, Value value,dummy<2> = 2) { |
|
1194 for(typename T::MapIt i(t); i!=INVALID; ++i){ |
|
1195 colUpperBound(*i, value); |
|
1196 } |
|
1197 } |
|
1198 #endif |
|
1199 |
|
1200 /// Set the lower and the upper bounds of a column (i.e a variable) |
|
1201 |
|
1202 /// The lower and the upper bounds of |
|
1203 /// a variable (column) have to be given by an |
|
1204 /// extended number of type Value, i.e. a finite number of type |
|
1205 /// Value, -\ref INF or \ref INF. |
|
1206 void colBounds(Col c, Value lower, Value upper) { |
|
1207 _setColLowerBound(_lpId(c),lower); |
|
1208 _setColUpperBound(_lpId(c),upper); |
|
1209 } |
|
1210 |
|
1211 ///\brief Set the lower and the upper bound of several columns |
|
1212 ///(i.e a variables) at once |
|
1213 /// |
|
1214 ///This magic function takes a container as its argument |
|
1215 ///and applies the function on all of its elements. |
|
1216 /// The lower and the upper bounds of |
|
1217 /// a variable (column) have to be given by an |
|
1218 /// extended number of type Value, i.e. a finite number of type |
|
1219 /// Value, -\ref INF or \ref INF. |
|
1220 #ifdef DOXYGEN |
|
1221 template<class T> |
|
1222 void colBounds(T &t, Value lower, Value upper) { return 0;} |
|
1223 #else |
|
1224 template<class T> |
|
1225 typename enable_if<typename T::value_type::LpSolverCol,void>::type |
|
1226 colBounds(T &t, Value lower, Value upper,dummy<0> = 0) { |
|
1227 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
1228 colBounds(*i, lower, upper); |
|
1229 } |
|
1230 } |
|
1231 template<class T> |
|
1232 typename enable_if<typename T::value_type::second_type::LpSolverCol, |
|
1233 void>::type |
|
1234 colBounds(T &t, Value lower, Value upper,dummy<1> = 1) { |
|
1235 for(typename T::iterator i=t.begin();i!=t.end();++i) { |
|
1236 colBounds(i->second, lower, upper); |
|
1237 } |
|
1238 } |
|
1239 template<class T> |
|
1240 typename enable_if<typename T::MapIt::Value::LpSolverCol, |
|
1241 void>::type |
|
1242 colBounds(T &t, Value lower, Value upper,dummy<2> = 2) { |
|
1243 for(typename T::MapIt i(t); i!=INVALID; ++i){ |
|
1244 colBounds(*i, lower, upper); |
|
1245 } |
|
1246 } |
|
1247 #endif |
|
1248 |
|
1249 |
|
1250 /// Set the lower and the upper bounds of a row (i.e a constraint) |
|
1251 |
|
1252 /// The lower and the upper bound of a constraint (row) have to be |
|
1253 /// given by an extended number of type Value, i.e. a finite |
|
1254 /// number of type Value, -\ref INF or \ref INF. There is no |
|
1255 /// separate function for the lower and the upper bound because |
|
1256 /// that would have been hard to implement for CPLEX. |
|
1257 void rowBounds(Row c, Value lower, Value upper) { |
|
1258 _setRowBounds(_lpId(c),lower, upper); |
|
1259 } |
|
1260 |
|
1261 /// Get the lower and the upper bounds of a row (i.e a constraint) |
|
1262 |
|
1263 /// The lower and the upper bound of |
|
1264 /// a constraint (row) are |
|
1265 /// extended numbers of type Value, i.e. finite numbers of type |
|
1266 /// Value, -\ref INF or \ref INF. |
|
1267 /// \todo There is no separate function for the |
|
1268 /// lower and the upper bound because we had problems with the |
|
1269 /// implementation of the setting functions for CPLEX: |
|
1270 /// check out whether this can be done for these functions. |
|
1271 void getRowBounds(Row c, Value &lower, Value &upper) const { |
|
1272 _getRowBounds(_lpId(c),lower, upper); |
|
1273 } |
|
1274 |
|
1275 ///Set an element of the objective function |
|
1276 void objCoeff(Col c, Value v) {_setObjCoeff(_lpId(c),v); }; |
|
1277 |
|
1278 ///Get an element of the objective function |
|
1279 Value objCoeff(Col c) const { return _getObjCoeff(_lpId(c)); }; |
|
1280 |
|
1281 ///Set the objective function |
|
1282 |
|
1283 ///\param e is a linear expression of type \ref Expr. |
|
1284 void obj(Expr e) { |
|
1285 _clearObj(); |
|
1286 for (Expr::iterator i=e.begin(); i!=e.end(); ++i) |
|
1287 objCoeff((*i).first,(*i).second); |
|
1288 obj_const_comp=e.constComp(); |
|
1289 } |
|
1290 |
|
1291 ///Get the objective function |
|
1292 |
|
1293 ///\return the objective function as a linear expression of type \ref Expr. |
|
1294 Expr obj() const { |
|
1295 Expr e; |
|
1296 for (ColIt it(*this); it != INVALID; ++it) { |
|
1297 double c = objCoeff(it); |
|
1298 if (c != 0.0) { |
|
1299 e.insert(std::make_pair(it, c)); |
|
1300 } |
|
1301 } |
|
1302 return e; |
|
1303 } |
|
1304 |
|
1305 |
|
1306 ///Maximize |
|
1307 void max() { _setMax(); } |
|
1308 ///Minimize |
|
1309 void min() { _setMin(); } |
|
1310 |
|
1311 ///Query function: is this a maximization problem? |
|
1312 bool isMax() const {return _isMax(); } |
|
1313 |
|
1314 ///Query function: is this a minimization problem? |
|
1315 bool isMin() const {return !isMax(); } |
|
1316 |
|
1317 ///@} |
|
1318 |
|
1319 |
|
1320 ///\name Solve the LP |
|
1321 |
|
1322 ///@{ |
|
1323 |
|
1324 ///\e Solve the LP problem at hand |
|
1325 /// |
|
1326 ///\return The result of the optimization procedure. Possible |
|
1327 ///values and their meanings can be found in the documentation of |
|
1328 ///\ref SolveExitStatus. |
|
1329 /// |
|
1330 ///\todo Which method is used to solve the problem |
|
1331 SolveExitStatus solve() { return _solve(); } |
|
1332 |
|
1333 ///@} |
|
1334 |
|
1335 ///\name Obtain the solution |
|
1336 |
|
1337 ///@{ |
|
1338 |
|
1339 /// The status of the primal problem (the original LP problem) |
|
1340 SolutionStatus primalStatus() const { |
|
1341 return _getPrimalStatus(); |
|
1342 } |
|
1343 |
|
1344 /// The status of the dual (of the original LP) problem |
|
1345 SolutionStatus dualStatus() const { |
|
1346 return _getDualStatus(); |
|
1347 } |
|
1348 |
|
1349 ///The type of the original LP problem |
|
1350 ProblemTypes problemType() const { |
|
1351 return _getProblemType(); |
|
1352 } |
|
1353 |
|
1354 ///\e |
|
1355 Value primal(Col c) const { return _getPrimal(_lpId(c)); } |
|
1356 ///\e |
|
1357 Value primal(const Expr& e) const { |
|
1358 double res = e.constComp(); |
|
1359 for (std::map<Col, double>::const_iterator it = e.begin(); |
|
1360 it != e.end(); ++it) { |
|
1361 res += _getPrimal(_lpId(it->first)) * it->second; |
|
1362 } |
|
1363 return res; |
|
1364 } |
|
1365 |
|
1366 ///\e |
|
1367 Value dual(Row r) const { return _getDual(_lpId(r)); } |
|
1368 ///\e |
|
1369 Value dual(const DualExpr& e) const { |
|
1370 double res = 0.0; |
|
1371 for (std::map<Row, double>::const_iterator it = e.begin(); |
|
1372 it != e.end(); ++it) { |
|
1373 res += _getPrimal(_lpId(it->first)) * it->second; |
|
1374 } |
|
1375 return res; |
|
1376 } |
|
1377 |
|
1378 ///\e |
|
1379 bool isBasicCol(Col c) const { return _isBasicCol(_lpId(c)); } |
|
1380 |
|
1381 ///\e |
|
1382 |
|
1383 ///\return |
|
1384 ///- \ref INF or -\ref INF means either infeasibility or unboundedness |
|
1385 /// of the primal problem, depending on whether we minimize or maximize. |
|
1386 ///- \ref NaN if no primal solution is found. |
|
1387 ///- The (finite) objective value if an optimal solution is found. |
|
1388 Value primalValue() const { return _getPrimalValue()+obj_const_comp;} |
|
1389 ///@} |
|
1390 |
|
1391 }; |
|
1392 |
|
1393 |
|
1394 /// \ingroup lp_group |
|
1395 /// |
|
1396 /// \brief Common base class for MIP solvers |
|
1397 /// \todo Much more docs |
|
1398 class MipSolverBase : virtual public LpSolverBase{ |
|
1399 public: |
|
1400 |
|
1401 ///Possible variable (coloumn) types (e.g. real, integer, binary etc.) |
|
1402 enum ColTypes { |
|
1403 ///Continuous variable |
|
1404 REAL = 0, |
|
1405 ///Integer variable |
|
1406 |
|
1407 ///Unfortunately, cplex 7.5 somewhere writes something like |
|
1408 ///#define INTEGER 'I' |
|
1409 INT = 1 |
|
1410 ///\todo No support for other types yet. |
|
1411 }; |
|
1412 |
|
1413 ///Sets the type of the given coloumn to the given type |
|
1414 /// |
|
1415 ///Sets the type of the given coloumn to the given type. |
|
1416 void colType(Col c, ColTypes col_type) { |
|
1417 _colType(_lpId(c),col_type); |
|
1418 } |
|
1419 |
|
1420 ///Gives back the type of the column. |
|
1421 /// |
|
1422 ///Gives back the type of the column. |
|
1423 ColTypes colType(Col c) const { |
|
1424 return _colType(_lpId(c)); |
|
1425 } |
|
1426 |
|
1427 ///Sets the type of the given Col to integer or remove that property. |
|
1428 /// |
|
1429 ///Sets the type of the given Col to integer or remove that property. |
|
1430 void integer(Col c, bool enable) { |
|
1431 if (enable) |
|
1432 colType(c,INT); |
|
1433 else |
|
1434 colType(c,REAL); |
|
1435 } |
|
1436 |
|
1437 ///Gives back whether the type of the column is integer or not. |
|
1438 /// |
|
1439 ///Gives back the type of the column. |
|
1440 ///\return true if the column has integer type and false if not. |
|
1441 bool integer(Col c) const { |
|
1442 return (colType(c)==INT); |
|
1443 } |
|
1444 |
|
1445 /// The status of the MIP problem |
|
1446 SolutionStatus mipStatus() const { |
|
1447 return _getMipStatus(); |
|
1448 } |
|
1449 |
|
1450 protected: |
|
1451 |
|
1452 virtual ColTypes _colType(int col) const = 0; |
|
1453 virtual void _colType(int col, ColTypes col_type) = 0; |
|
1454 virtual SolutionStatus _getMipStatus() const = 0; |
|
1455 |
|
1456 }; |
|
1457 |
|
1458 ///\relates LpSolverBase::Expr |
|
1459 /// |
|
1460 inline LpSolverBase::Expr operator+(const LpSolverBase::Expr &a, |
|
1461 const LpSolverBase::Expr &b) |
|
1462 { |
|
1463 LpSolverBase::Expr tmp(a); |
|
1464 tmp+=b; |
|
1465 return tmp; |
|
1466 } |
|
1467 ///\e |
|
1468 |
|
1469 ///\relates LpSolverBase::Expr |
|
1470 /// |
|
1471 inline LpSolverBase::Expr operator-(const LpSolverBase::Expr &a, |
|
1472 const LpSolverBase::Expr &b) |
|
1473 { |
|
1474 LpSolverBase::Expr tmp(a); |
|
1475 tmp-=b; |
|
1476 return tmp; |
|
1477 } |
|
1478 ///\e |
|
1479 |
|
1480 ///\relates LpSolverBase::Expr |
|
1481 /// |
|
1482 inline LpSolverBase::Expr operator*(const LpSolverBase::Expr &a, |
|
1483 const LpSolverBase::Value &b) |
|
1484 { |
|
1485 LpSolverBase::Expr tmp(a); |
|
1486 tmp*=b; |
|
1487 return tmp; |
|
1488 } |
|
1489 |
|
1490 ///\e |
|
1491 |
|
1492 ///\relates LpSolverBase::Expr |
|
1493 /// |
|
1494 inline LpSolverBase::Expr operator*(const LpSolverBase::Value &a, |
|
1495 const LpSolverBase::Expr &b) |
|
1496 { |
|
1497 LpSolverBase::Expr tmp(b); |
|
1498 tmp*=a; |
|
1499 return tmp; |
|
1500 } |
|
1501 ///\e |
|
1502 |
|
1503 ///\relates LpSolverBase::Expr |
|
1504 /// |
|
1505 inline LpSolverBase::Expr operator/(const LpSolverBase::Expr &a, |
|
1506 const LpSolverBase::Value &b) |
|
1507 { |
|
1508 LpSolverBase::Expr tmp(a); |
|
1509 tmp/=b; |
|
1510 return tmp; |
|
1511 } |
|
1512 |
|
1513 ///\e |
|
1514 |
|
1515 ///\relates LpSolverBase::Constr |
|
1516 /// |
|
1517 inline LpSolverBase::Constr operator<=(const LpSolverBase::Expr &e, |
|
1518 const LpSolverBase::Expr &f) |
|
1519 { |
|
1520 return LpSolverBase::Constr(-LpSolverBase::INF,e-f,0); |
|
1521 } |
|
1522 |
|
1523 ///\e |
|
1524 |
|
1525 ///\relates LpSolverBase::Constr |
|
1526 /// |
|
1527 inline LpSolverBase::Constr operator<=(const LpSolverBase::Value &e, |
|
1528 const LpSolverBase::Expr &f) |
|
1529 { |
|
1530 return LpSolverBase::Constr(e,f); |
|
1531 } |
|
1532 |
|
1533 ///\e |
|
1534 |
|
1535 ///\relates LpSolverBase::Constr |
|
1536 /// |
|
1537 inline LpSolverBase::Constr operator<=(const LpSolverBase::Expr &e, |
|
1538 const LpSolverBase::Value &f) |
|
1539 { |
|
1540 return LpSolverBase::Constr(-LpSolverBase::INF,e,f); |
|
1541 } |
|
1542 |
|
1543 ///\e |
|
1544 |
|
1545 ///\relates LpSolverBase::Constr |
|
1546 /// |
|
1547 inline LpSolverBase::Constr operator>=(const LpSolverBase::Expr &e, |
|
1548 const LpSolverBase::Expr &f) |
|
1549 { |
|
1550 return LpSolverBase::Constr(-LpSolverBase::INF,f-e,0); |
|
1551 } |
|
1552 |
|
1553 |
|
1554 ///\e |
|
1555 |
|
1556 ///\relates LpSolverBase::Constr |
|
1557 /// |
|
1558 inline LpSolverBase::Constr operator>=(const LpSolverBase::Value &e, |
|
1559 const LpSolverBase::Expr &f) |
|
1560 { |
|
1561 return LpSolverBase::Constr(f,e); |
|
1562 } |
|
1563 |
|
1564 |
|
1565 ///\e |
|
1566 |
|
1567 ///\relates LpSolverBase::Constr |
|
1568 /// |
|
1569 inline LpSolverBase::Constr operator>=(const LpSolverBase::Expr &e, |
|
1570 const LpSolverBase::Value &f) |
|
1571 { |
|
1572 return LpSolverBase::Constr(f,e,LpSolverBase::INF); |
|
1573 } |
|
1574 |
|
1575 ///\e |
|
1576 |
|
1577 ///\relates LpSolverBase::Constr |
|
1578 /// |
|
1579 inline LpSolverBase::Constr operator==(const LpSolverBase::Expr &e, |
|
1580 const LpSolverBase::Value &f) |
|
1581 { |
|
1582 return LpSolverBase::Constr(f,e,f); |
|
1583 } |
|
1584 |
|
1585 ///\e |
|
1586 |
|
1587 ///\relates LpSolverBase::Constr |
|
1588 /// |
|
1589 inline LpSolverBase::Constr operator==(const LpSolverBase::Expr &e, |
|
1590 const LpSolverBase::Expr &f) |
|
1591 { |
|
1592 return LpSolverBase::Constr(0,e-f,0); |
|
1593 } |
|
1594 |
|
1595 ///\e |
|
1596 |
|
1597 ///\relates LpSolverBase::Constr |
|
1598 /// |
|
1599 inline LpSolverBase::Constr operator<=(const LpSolverBase::Value &n, |
|
1600 const LpSolverBase::Constr&c) |
|
1601 { |
|
1602 LpSolverBase::Constr tmp(c); |
|
1603 LEMON_ASSERT(LpSolverBase::isNaN(tmp.lowerBound()), "Wrong LP constraint"); |
|
1604 tmp.lowerBound()=n; |
|
1605 return tmp; |
|
1606 } |
|
1607 ///\e |
|
1608 |
|
1609 ///\relates LpSolverBase::Constr |
|
1610 /// |
|
1611 inline LpSolverBase::Constr operator<=(const LpSolverBase::Constr& c, |
|
1612 const LpSolverBase::Value &n) |
|
1613 { |
|
1614 LpSolverBase::Constr tmp(c); |
|
1615 LEMON_ASSERT(LpSolverBase::isNaN(tmp.upperBound()), "Wrong LP constraint"); |
|
1616 tmp.upperBound()=n; |
|
1617 return tmp; |
|
1618 } |
|
1619 |
|
1620 ///\e |
|
1621 |
|
1622 ///\relates LpSolverBase::Constr |
|
1623 /// |
|
1624 inline LpSolverBase::Constr operator>=(const LpSolverBase::Value &n, |
|
1625 const LpSolverBase::Constr&c) |
|
1626 { |
|
1627 LpSolverBase::Constr tmp(c); |
|
1628 LEMON_ASSERT(LpSolverBase::isNaN(tmp.upperBound()), "Wrong LP constraint"); |
|
1629 tmp.upperBound()=n; |
|
1630 return tmp; |
|
1631 } |
|
1632 ///\e |
|
1633 |
|
1634 ///\relates LpSolverBase::Constr |
|
1635 /// |
|
1636 inline LpSolverBase::Constr operator>=(const LpSolverBase::Constr& c, |
|
1637 const LpSolverBase::Value &n) |
|
1638 { |
|
1639 LpSolverBase::Constr tmp(c); |
|
1640 LEMON_ASSERT(LpSolverBase::isNaN(tmp.lowerBound()), "Wrong LP constraint"); |
|
1641 tmp.lowerBound()=n; |
|
1642 return tmp; |
|
1643 } |
|
1644 |
|
1645 ///\e |
|
1646 |
|
1647 ///\relates LpSolverBase::DualExpr |
|
1648 /// |
|
1649 inline LpSolverBase::DualExpr operator+(const LpSolverBase::DualExpr &a, |
|
1650 const LpSolverBase::DualExpr &b) |
|
1651 { |
|
1652 LpSolverBase::DualExpr tmp(a); |
|
1653 tmp+=b; |
|
1654 return tmp; |
|
1655 } |
|
1656 ///\e |
|
1657 |
|
1658 ///\relates LpSolverBase::DualExpr |
|
1659 /// |
|
1660 inline LpSolverBase::DualExpr operator-(const LpSolverBase::DualExpr &a, |
|
1661 const LpSolverBase::DualExpr &b) |
|
1662 { |
|
1663 LpSolverBase::DualExpr tmp(a); |
|
1664 tmp-=b; |
|
1665 return tmp; |
|
1666 } |
|
1667 ///\e |
|
1668 |
|
1669 ///\relates LpSolverBase::DualExpr |
|
1670 /// |
|
1671 inline LpSolverBase::DualExpr operator*(const LpSolverBase::DualExpr &a, |
|
1672 const LpSolverBase::Value &b) |
|
1673 { |
|
1674 LpSolverBase::DualExpr tmp(a); |
|
1675 tmp*=b; |
|
1676 return tmp; |
|
1677 } |
|
1678 |
|
1679 ///\e |
|
1680 |
|
1681 ///\relates LpSolverBase::DualExpr |
|
1682 /// |
|
1683 inline LpSolverBase::DualExpr operator*(const LpSolverBase::Value &a, |
|
1684 const LpSolverBase::DualExpr &b) |
|
1685 { |
|
1686 LpSolverBase::DualExpr tmp(b); |
|
1687 tmp*=a; |
|
1688 return tmp; |
|
1689 } |
|
1690 ///\e |
|
1691 |
|
1692 ///\relates LpSolverBase::DualExpr |
|
1693 /// |
|
1694 inline LpSolverBase::DualExpr operator/(const LpSolverBase::DualExpr &a, |
|
1695 const LpSolverBase::Value &b) |
|
1696 { |
|
1697 LpSolverBase::DualExpr tmp(a); |
|
1698 tmp/=b; |
|
1699 return tmp; |
|
1700 } |
|
1701 |
|
1702 |
|
1703 } //namespace lemon |
|
1704 |
|
1705 #endif //LEMON_LP_BASE_H |