lemon/lp_base.h
author Alpar Juttner <alpar@cs.elte.hu>
Wed, 17 Oct 2018 22:55:02 +0200
changeset 1195 73e29215aaa4
parent 1130 0759d974de81
permissions -rw-r--r--
Add citation for Vf2pp (#597)
deba@458
     1
/* -*- mode: C++; indent-tabs-mode: nil; -*-
deba@458
     2
 *
deba@458
     3
 * This file is a part of LEMON, a generic C++ optimization library.
deba@458
     4
 *
alpar@1092
     5
 * Copyright (C) 2003-2013
deba@458
     6
 * Egervary Jeno Kombinatorikus Optimalizalasi Kutatocsoport
deba@458
     7
 * (Egervary Research Group on Combinatorial Optimization, EGRES).
deba@458
     8
 *
deba@458
     9
 * Permission to use, modify and distribute this software is granted
deba@458
    10
 * provided that this copyright notice appears in all copies. For
deba@458
    11
 * precise terms see the accompanying LICENSE file.
deba@458
    12
 *
deba@458
    13
 * This software is provided "AS IS" with no warranty of any kind,
deba@458
    14
 * express or implied, and with no claim as to its suitability for any
deba@458
    15
 * purpose.
deba@458
    16
 *
deba@458
    17
 */
deba@458
    18
deba@458
    19
#ifndef LEMON_LP_BASE_H
deba@458
    20
#define LEMON_LP_BASE_H
deba@458
    21
deba@458
    22
#include<iostream>
deba@458
    23
#include<vector>
deba@458
    24
#include<map>
deba@458
    25
#include<limits>
deba@458
    26
#include<lemon/math.h>
deba@458
    27
deba@459
    28
#include<lemon/error.h>
deba@459
    29
#include<lemon/assert.h>
deba@459
    30
deba@458
    31
#include<lemon/core.h>
deba@459
    32
#include<lemon/bits/solver_bits.h>
deba@458
    33
ggab90@1130
    34
#include<lemon/bits/stl_iterators.h>
ggab90@1130
    35
deba@458
    36
///\file
deba@458
    37
///\brief The interface of the LP solver interface.
deba@458
    38
///\ingroup lp_group
deba@458
    39
namespace lemon {
deba@458
    40
deba@459
    41
  ///Common base class for LP and MIP solvers
deba@458
    42
deba@459
    43
  ///Usually this class is not used directly, please use one of the concrete
deba@459
    44
  ///implementations of the solver interface.
deba@458
    45
  ///\ingroup lp_group
deba@459
    46
  class LpBase {
deba@458
    47
deba@458
    48
  protected:
deba@458
    49
ggab90@1130
    50
    _solver_bits::VarIndex _rows;
ggab90@1130
    51
    _solver_bits::VarIndex _cols;
deba@458
    52
deba@458
    53
  public:
deba@458
    54
deba@458
    55
    ///Possible outcomes of an LP solving procedure
deba@458
    56
    enum SolveExitStatus {
kpeter@584
    57
      /// = 0. It means that the problem has been successfully solved: either
deba@458
    58
      ///an optimal solution has been found or infeasibility/unboundedness
deba@458
    59
      ///has been proved.
deba@458
    60
      SOLVED = 0,
kpeter@584
    61
      /// = 1. Any other case (including the case when some user specified
kpeter@584
    62
      ///limit has been exceeded).
deba@458
    63
      UNSOLVED = 1
deba@458
    64
    };
deba@458
    65
deba@459
    66
    ///Direction of the optimization
deba@459
    67
    enum Sense {
deba@459
    68
      /// Minimization
deba@459
    69
      MIN,
deba@459
    70
      /// Maximization
deba@459
    71
      MAX
deba@458
    72
    };
deba@458
    73
deba@576
    74
    ///Enum for \c messageLevel() parameter
deba@576
    75
    enum MessageLevel {
kpeter@584
    76
      /// No output (default value).
deba@576
    77
      MESSAGE_NOTHING,
kpeter@584
    78
      /// Error messages only.
deba@576
    79
      MESSAGE_ERROR,
kpeter@584
    80
      /// Warnings.
deba@576
    81
      MESSAGE_WARNING,
kpeter@584
    82
      /// Normal output.
deba@576
    83
      MESSAGE_NORMAL,
kpeter@584
    84
      /// Verbose output.
deba@576
    85
      MESSAGE_VERBOSE
deba@576
    86
    };
alpar@877
    87
deba@576
    88
deba@458
    89
    ///The floating point type used by the solver
deba@458
    90
    typedef double Value;
deba@458
    91
    ///The infinity constant
deba@458
    92
    static const Value INF;
deba@458
    93
    ///The not a number constant
deba@458
    94
    static const Value NaN;
deba@458
    95
deba@458
    96
    friend class Col;
deba@458
    97
    friend class ColIt;
deba@458
    98
    friend class Row;
deba@459
    99
    friend class RowIt;
deba@458
   100
deba@458
   101
    ///Refer to a column of the LP.
deba@458
   102
deba@458
   103
    ///This type is used to refer to a column of the LP.
deba@458
   104
    ///
deba@458
   105
    ///Its value remains valid and correct even after the addition or erase of
deba@458
   106
    ///other columns.
deba@458
   107
    ///
deba@459
   108
    ///\note This class is similar to other Item types in LEMON, like
deba@459
   109
    ///Node and Arc types in digraph.
deba@458
   110
    class Col {
deba@459
   111
      friend class LpBase;
deba@458
   112
    protected:
deba@459
   113
      int _id;
deba@459
   114
      explicit Col(int id) : _id(id) {}
deba@458
   115
    public:
deba@458
   116
      typedef Value ExprValue;
deba@459
   117
      typedef True LpCol;
deba@459
   118
      /// Default constructor
alpar@877
   119
deba@459
   120
      /// \warning The default constructor sets the Col to an
deba@459
   121
      /// undefined value.
deba@458
   122
      Col() {}
deba@459
   123
      /// Invalid constructor \& conversion.
alpar@877
   124
deba@459
   125
      /// This constructor initializes the Col to be invalid.
alpar@877
   126
      /// \sa Invalid for more details.
deba@459
   127
      Col(const Invalid&) : _id(-1) {}
deba@459
   128
      /// Equality operator
deba@459
   129
deba@459
   130
      /// Two \ref Col "Col"s are equal if and only if they point to
deba@459
   131
      /// the same LP column or both are invalid.
deba@459
   132
      bool operator==(Col c) const  {return _id == c._id;}
deba@459
   133
      /// Inequality operator
deba@459
   134
deba@459
   135
      /// \sa operator==(Col c)
deba@459
   136
      ///
deba@459
   137
      bool operator!=(Col c) const  {return _id != c._id;}
deba@459
   138
      /// Artificial ordering operator.
deba@459
   139
deba@459
   140
      /// To allow the use of this object in std::map or similar
deba@459
   141
      /// associative container we require this.
deba@459
   142
      ///
deba@459
   143
      /// \note This operator only have to define some strict ordering of
deba@459
   144
      /// the items; this order has nothing to do with the iteration
deba@459
   145
      /// ordering of the items.
deba@459
   146
      bool operator<(Col c) const  {return _id < c._id;}
deba@458
   147
    };
deba@458
   148
deba@459
   149
    ///Iterator for iterate over the columns of an LP problem
deba@459
   150
kpeter@786
   151
    /// Its usage is quite simple, for example, you can count the number
deba@459
   152
    /// of columns in an LP \c lp:
deba@459
   153
    ///\code
deba@459
   154
    /// int count=0;
deba@459
   155
    /// for (LpBase::ColIt c(lp); c!=INVALID; ++c) ++count;
deba@459
   156
    ///\endcode
deba@458
   157
    class ColIt : public Col {
deba@459
   158
      const LpBase *_solver;
deba@458
   159
    public:
deba@459
   160
      /// Default constructor
alpar@877
   161
deba@459
   162
      /// \warning The default constructor sets the iterator
deba@459
   163
      /// to an undefined value.
deba@458
   164
      ColIt() {}
deba@459
   165
      /// Sets the iterator to the first Col
alpar@877
   166
deba@459
   167
      /// Sets the iterator to the first Col.
deba@459
   168
      ///
deba@459
   169
      ColIt(const LpBase &solver) : _solver(&solver)
deba@458
   170
      {
ggab90@1130
   171
        _solver->_cols.firstItem(_id);
deba@458
   172
      }
deba@459
   173
      /// Invalid constructor \& conversion
alpar@877
   174
deba@459
   175
      /// Initialize the iterator to be invalid.
deba@459
   176
      /// \sa Invalid for more details.
deba@458
   177
      ColIt(const Invalid&) : Col(INVALID) {}
deba@459
   178
      /// Next column
alpar@877
   179
deba@459
   180
      /// Assign the iterator to the next column.
deba@459
   181
      ///
deba@458
   182
      ColIt &operator++()
deba@458
   183
      {
ggab90@1130
   184
        _solver->_cols.nextItem(_id);
deba@458
   185
        return *this;
deba@458
   186
      }
deba@458
   187
    };
deba@458
   188
ggab90@1130
   189
    /// \brief Gets the collection of the columns of the LP problem.
ggab90@1130
   190
    ///
ggab90@1130
   191
    /// This function can be used for iterating on
ggab90@1130
   192
    /// the columns of the LP problem. It returns a wrapped ColIt, which looks
ggab90@1130
   193
    /// like an STL container (by having begin() and end())
ggab90@1130
   194
    /// which you can use in range-based for loops, STL algorithms, etc.
ggab90@1130
   195
    /// For example you can write:
ggab90@1130
   196
    ///\code
ggab90@1130
   197
    /// for(auto c: lp.cols())
ggab90@1130
   198
    ///   doSomething(c);
alpar@1144
   199
    ///\endcode
ggab90@1130
   200
    LemonRangeWrapper1<ColIt, LpBase> cols() {
ggab90@1130
   201
      return LemonRangeWrapper1<ColIt, LpBase>(*this);
ggab90@1130
   202
    }
ggab90@1130
   203
ggab90@1130
   204
deba@459
   205
    /// \brief Returns the ID of the column.
deba@459
   206
    static int id(const Col& col) { return col._id; }
deba@459
   207
    /// \brief Returns the column with the given ID.
deba@459
   208
    ///
deba@459
   209
    /// \pre The argument should be a valid column ID in the LP problem.
deba@459
   210
    static Col colFromId(int id) { return Col(id); }
deba@458
   211
deba@458
   212
    ///Refer to a row of the LP.
deba@458
   213
deba@458
   214
    ///This type is used to refer to a row of the LP.
deba@458
   215
    ///
deba@458
   216
    ///Its value remains valid and correct even after the addition or erase of
deba@458
   217
    ///other rows.
deba@458
   218
    ///
deba@459
   219
    ///\note This class is similar to other Item types in LEMON, like
deba@459
   220
    ///Node and Arc types in digraph.
deba@458
   221
    class Row {
deba@459
   222
      friend class LpBase;
deba@458
   223
    protected:
deba@459
   224
      int _id;
deba@459
   225
      explicit Row(int id) : _id(id) {}
deba@458
   226
    public:
deba@458
   227
      typedef Value ExprValue;
deba@459
   228
      typedef True LpRow;
deba@459
   229
      /// Default constructor
alpar@877
   230
deba@459
   231
      /// \warning The default constructor sets the Row to an
deba@459
   232
      /// undefined value.
deba@458
   233
      Row() {}
deba@459
   234
      /// Invalid constructor \& conversion.
alpar@877
   235
deba@459
   236
      /// This constructor initializes the Row to be invalid.
alpar@877
   237
      /// \sa Invalid for more details.
deba@459
   238
      Row(const Invalid&) : _id(-1) {}
deba@459
   239
      /// Equality operator
deba@458
   240
deba@459
   241
      /// Two \ref Row "Row"s are equal if and only if they point to
deba@459
   242
      /// the same LP row or both are invalid.
deba@459
   243
      bool operator==(Row r) const  {return _id == r._id;}
deba@459
   244
      /// Inequality operator
alpar@877
   245
deba@459
   246
      /// \sa operator==(Row r)
deba@459
   247
      ///
deba@459
   248
      bool operator!=(Row r) const  {return _id != r._id;}
deba@459
   249
      /// Artificial ordering operator.
deba@459
   250
deba@459
   251
      /// To allow the use of this object in std::map or similar
deba@459
   252
      /// associative container we require this.
deba@459
   253
      ///
deba@459
   254
      /// \note This operator only have to define some strict ordering of
deba@459
   255
      /// the items; this order has nothing to do with the iteration
deba@459
   256
      /// ordering of the items.
deba@459
   257
      bool operator<(Row r) const  {return _id < r._id;}
deba@458
   258
    };
deba@458
   259
deba@459
   260
    ///Iterator for iterate over the rows of an LP problem
deba@459
   261
kpeter@786
   262
    /// Its usage is quite simple, for example, you can count the number
deba@459
   263
    /// of rows in an LP \c lp:
deba@459
   264
    ///\code
deba@459
   265
    /// int count=0;
deba@459
   266
    /// for (LpBase::RowIt c(lp); c!=INVALID; ++c) ++count;
deba@459
   267
    ///\endcode
deba@458
   268
    class RowIt : public Row {
deba@459
   269
      const LpBase *_solver;
deba@458
   270
    public:
deba@459
   271
      /// Default constructor
alpar@877
   272
deba@459
   273
      /// \warning The default constructor sets the iterator
deba@459
   274
      /// to an undefined value.
deba@458
   275
      RowIt() {}
deba@459
   276
      /// Sets the iterator to the first Row
alpar@877
   277
deba@459
   278
      /// Sets the iterator to the first Row.
deba@459
   279
      ///
deba@459
   280
      RowIt(const LpBase &solver) : _solver(&solver)
deba@458
   281
      {
ggab90@1130
   282
        _solver->_rows.firstItem(_id);
deba@458
   283
      }
deba@459
   284
      /// Invalid constructor \& conversion
alpar@877
   285
deba@459
   286
      /// Initialize the iterator to be invalid.
deba@459
   287
      /// \sa Invalid for more details.
deba@458
   288
      RowIt(const Invalid&) : Row(INVALID) {}
deba@459
   289
      /// Next row
alpar@877
   290
deba@459
   291
      /// Assign the iterator to the next row.
deba@459
   292
      ///
deba@458
   293
      RowIt &operator++()
deba@458
   294
      {
ggab90@1130
   295
        _solver->_rows.nextItem(_id);
deba@458
   296
        return *this;
deba@458
   297
      }
deba@458
   298
    };
ggab90@1130
   299
    
ggab90@1130
   300
    /// \brief Gets the collection of the rows of the LP problem.
ggab90@1130
   301
    ///
ggab90@1130
   302
    /// This function can be used for iterating on
ggab90@1130
   303
    /// the rows of the LP problem. It returns a wrapped RowIt, which looks
ggab90@1130
   304
    /// like an STL container (by having begin() and end())
ggab90@1130
   305
    /// which you can use in range-based for loops, STL algorithms, etc.
ggab90@1130
   306
    /// For example you can write:
ggab90@1130
   307
    ///\code
ggab90@1130
   308
    /// for(auto c: lp.rows())
ggab90@1130
   309
    ///   doSomething(c);
alpar@1144
   310
    ///\endcode
ggab90@1130
   311
    LemonRangeWrapper1<RowIt, LpBase> rows() {
ggab90@1130
   312
      return LemonRangeWrapper1<RowIt, LpBase>(*this);
ggab90@1130
   313
    }
ggab90@1130
   314
    
deba@458
   315
deba@459
   316
    /// \brief Returns the ID of the row.
deba@459
   317
    static int id(const Row& row) { return row._id; }
deba@459
   318
    /// \brief Returns the row with the given ID.
deba@459
   319
    ///
deba@459
   320
    /// \pre The argument should be a valid row ID in the LP problem.
deba@459
   321
    static Row rowFromId(int id) { return Row(id); }
deba@458
   322
deba@458
   323
  public:
deba@458
   324
deba@458
   325
    ///Linear expression of variables and a constant component
deba@458
   326
deba@458
   327
    ///This data structure stores a linear expression of the variables
deba@458
   328
    ///(\ref Col "Col"s) and also has a constant component.
deba@458
   329
    ///
deba@458
   330
    ///There are several ways to access and modify the contents of this
deba@458
   331
    ///container.
deba@458
   332
    ///\code
deba@458
   333
    ///e[v]=5;
deba@458
   334
    ///e[v]+=12;
deba@458
   335
    ///e.erase(v);
deba@458
   336
    ///\endcode
deba@458
   337
    ///or you can also iterate through its elements.
deba@458
   338
    ///\code
deba@458
   339
    ///double s=0;
deba@459
   340
    ///for(LpBase::Expr::ConstCoeffIt i(e);i!=INVALID;++i)
deba@459
   341
    ///  s+=*i * primal(i);
deba@458
   342
    ///\endcode
deba@459
   343
    ///(This code computes the primal value of the expression).
deba@458
   344
    ///- Numbers (<tt>double</tt>'s)
deba@458
   345
    ///and variables (\ref Col "Col"s) directly convert to an
deba@458
   346
    ///\ref Expr and the usual linear operations are defined, so
deba@458
   347
    ///\code
deba@458
   348
    ///v+w
deba@458
   349
    ///2*v-3.12*(v-w/2)+2
deba@458
   350
    ///v*2.1+(3*v+(v*12+w+6)*3)/2
deba@458
   351
    ///\endcode
deba@459
   352
    ///are valid expressions.
deba@458
   353
    ///The usual assignment operations are also defined.
deba@458
   354
    ///\code
deba@458
   355
    ///e=v+w;
deba@458
   356
    ///e+=2*v-3.12*(v-w/2)+2;
deba@458
   357
    ///e*=3.4;
deba@458
   358
    ///e/=5;
deba@458
   359
    ///\endcode
deba@459
   360
    ///- The constant member can be set and read by dereference
deba@459
   361
    ///  operator (unary *)
deba@459
   362
    ///
deba@458
   363
    ///\code
deba@459
   364
    ///*e=12;
deba@459
   365
    ///double c=*e;
deba@458
   366
    ///\endcode
deba@458
   367
    ///
deba@458
   368
    ///\sa Constr
deba@459
   369
    class Expr {
deba@459
   370
      friend class LpBase;
deba@458
   371
    public:
deba@459
   372
      /// The key type of the expression
deba@459
   373
      typedef LpBase::Col Key;
deba@459
   374
      /// The value type of the expression
deba@459
   375
      typedef LpBase::Value Value;
deba@458
   376
deba@458
   377
    protected:
deba@459
   378
      Value const_comp;
deba@459
   379
      std::map<int, Value> comps;
deba@458
   380
deba@458
   381
    public:
deba@459
   382
      typedef True SolverExpr;
deba@459
   383
      /// Default constructor
alpar@877
   384
deba@459
   385
      /// Construct an empty expression, the coefficients and
deba@459
   386
      /// the constant component are initialized to zero.
deba@459
   387
      Expr() : const_comp(0) {}
deba@459
   388
      /// Construct an expression from a column
deba@459
   389
deba@459
   390
      /// Construct an expression, which has a term with \c c variable
deba@459
   391
      /// and 1.0 coefficient.
deba@459
   392
      Expr(const Col &c) : const_comp(0) {
deba@459
   393
        typedef std::map<int, Value>::value_type pair_type;
deba@459
   394
        comps.insert(pair_type(id(c), 1));
deba@458
   395
      }
deba@459
   396
      /// Construct an expression from a constant
deba@459
   397
deba@459
   398
      /// Construct an expression, which's constant component is \c v.
deba@459
   399
      ///
deba@458
   400
      Expr(const Value &v) : const_comp(v) {}
deba@459
   401
      /// Returns the coefficient of the column
deba@459
   402
      Value operator[](const Col& c) const {
deba@459
   403
        std::map<int, Value>::const_iterator it=comps.find(id(c));
deba@459
   404
        if (it != comps.end()) {
deba@459
   405
          return it->second;
deba@459
   406
        } else {
deba@459
   407
          return 0;
deba@458
   408
        }
deba@458
   409
      }
deba@459
   410
      /// Returns the coefficient of the column
deba@459
   411
      Value& operator[](const Col& c) {
deba@459
   412
        return comps[id(c)];
deba@459
   413
      }
deba@459
   414
      /// Sets the coefficient of the column
deba@459
   415
      void set(const Col &c, const Value &v) {
deba@459
   416
        if (v != 0.0) {
deba@459
   417
          typedef std::map<int, Value>::value_type pair_type;
deba@459
   418
          comps.insert(pair_type(id(c), v));
deba@459
   419
        } else {
deba@459
   420
          comps.erase(id(c));
deba@459
   421
        }
deba@459
   422
      }
deba@459
   423
      /// Returns the constant component of the expression
deba@459
   424
      Value& operator*() { return const_comp; }
deba@459
   425
      /// Returns the constant component of the expression
deba@459
   426
      const Value& operator*() const { return const_comp; }
deba@459
   427
      /// \brief Removes the coefficients which's absolute value does
deba@459
   428
      /// not exceed \c epsilon. It also sets to zero the constant
deba@459
   429
      /// component, if it does not exceed epsilon in absolute value.
deba@459
   430
      void simplify(Value epsilon = 0.0) {
deba@459
   431
        std::map<int, Value>::iterator it=comps.begin();
deba@459
   432
        while (it != comps.end()) {
deba@459
   433
          std::map<int, Value>::iterator jt=it;
deba@459
   434
          ++jt;
deba@459
   435
          if (std::fabs((*it).second) <= epsilon) comps.erase(it);
deba@459
   436
          it=jt;
deba@459
   437
        }
deba@459
   438
        if (std::fabs(const_comp) <= epsilon) const_comp = 0;
deba@458
   439
      }
deba@458
   440
deba@459
   441
      void simplify(Value epsilon = 0.0) const {
deba@459
   442
        const_cast<Expr*>(this)->simplify(epsilon);
deba@458
   443
      }
deba@458
   444
deba@458
   445
      ///Sets all coefficients and the constant component to 0.
deba@458
   446
      void clear() {
deba@459
   447
        comps.clear();
deba@458
   448
        const_comp=0;
deba@458
   449
      }
deba@458
   450
deba@459
   451
      ///Compound assignment
deba@458
   452
      Expr &operator+=(const Expr &e) {
deba@459
   453
        for (std::map<int, Value>::const_iterator it=e.comps.begin();
deba@459
   454
             it!=e.comps.end(); ++it)
deba@459
   455
          comps[it->first]+=it->second;
deba@458
   456
        const_comp+=e.const_comp;
deba@458
   457
        return *this;
deba@458
   458
      }
deba@459
   459
      ///Compound assignment
deba@458
   460
      Expr &operator-=(const Expr &e) {
deba@459
   461
        for (std::map<int, Value>::const_iterator it=e.comps.begin();
deba@459
   462
             it!=e.comps.end(); ++it)
deba@459
   463
          comps[it->first]-=it->second;
deba@458
   464
        const_comp-=e.const_comp;
deba@458
   465
        return *this;
deba@458
   466
      }
deba@459
   467
      ///Multiply with a constant
deba@459
   468
      Expr &operator*=(const Value &v) {
deba@459
   469
        for (std::map<int, Value>::iterator it=comps.begin();
deba@459
   470
             it!=comps.end(); ++it)
deba@459
   471
          it->second*=v;
deba@459
   472
        const_comp*=v;
deba@458
   473
        return *this;
deba@458
   474
      }
deba@459
   475
      ///Division with a constant
deba@458
   476
      Expr &operator/=(const Value &c) {
deba@459
   477
        for (std::map<int, Value>::iterator it=comps.begin();
deba@459
   478
             it!=comps.end(); ++it)
deba@459
   479
          it->second/=c;
deba@458
   480
        const_comp/=c;
deba@458
   481
        return *this;
deba@458
   482
      }
deba@458
   483
deba@459
   484
      ///Iterator over the expression
alpar@877
   485
alpar@877
   486
      ///The iterator iterates over the terms of the expression.
alpar@877
   487
      ///
deba@459
   488
      ///\code
deba@459
   489
      ///double s=0;
deba@459
   490
      ///for(LpBase::Expr::CoeffIt i(e);i!=INVALID;++i)
deba@459
   491
      ///  s+= *i * primal(i);
deba@459
   492
      ///\endcode
deba@459
   493
      class CoeffIt {
deba@459
   494
      private:
deba@459
   495
deba@459
   496
        std::map<int, Value>::iterator _it, _end;
deba@459
   497
deba@459
   498
      public:
deba@459
   499
deba@459
   500
        /// Sets the iterator to the first term
alpar@877
   501
deba@459
   502
        /// Sets the iterator to the first term of the expression.
deba@459
   503
        ///
deba@459
   504
        CoeffIt(Expr& e)
deba@459
   505
          : _it(e.comps.begin()), _end(e.comps.end()){}
deba@459
   506
deba@459
   507
        /// Convert the iterator to the column of the term
deba@459
   508
        operator Col() const {
deba@459
   509
          return colFromId(_it->first);
deba@459
   510
        }
deba@459
   511
deba@459
   512
        /// Returns the coefficient of the term
deba@459
   513
        Value& operator*() { return _it->second; }
deba@459
   514
deba@459
   515
        /// Returns the coefficient of the term
deba@459
   516
        const Value& operator*() const { return _it->second; }
deba@459
   517
        /// Next term
alpar@877
   518
deba@459
   519
        /// Assign the iterator to the next term.
deba@459
   520
        ///
deba@459
   521
        CoeffIt& operator++() { ++_it; return *this; }
deba@459
   522
deba@459
   523
        /// Equality operator
deba@459
   524
        bool operator==(Invalid) const { return _it == _end; }
deba@459
   525
        /// Inequality operator
deba@459
   526
        bool operator!=(Invalid) const { return _it != _end; }
deba@459
   527
      };
deba@459
   528
deba@459
   529
      /// Const iterator over the expression
alpar@877
   530
alpar@877
   531
      ///The iterator iterates over the terms of the expression.
alpar@877
   532
      ///
deba@459
   533
      ///\code
deba@459
   534
      ///double s=0;
deba@459
   535
      ///for(LpBase::Expr::ConstCoeffIt i(e);i!=INVALID;++i)
deba@459
   536
      ///  s+=*i * primal(i);
deba@459
   537
      ///\endcode
deba@459
   538
      class ConstCoeffIt {
deba@459
   539
      private:
deba@459
   540
deba@459
   541
        std::map<int, Value>::const_iterator _it, _end;
deba@459
   542
deba@459
   543
      public:
deba@459
   544
deba@459
   545
        /// Sets the iterator to the first term
alpar@877
   546
deba@459
   547
        /// Sets the iterator to the first term of the expression.
deba@459
   548
        ///
deba@459
   549
        ConstCoeffIt(const Expr& e)
deba@459
   550
          : _it(e.comps.begin()), _end(e.comps.end()){}
deba@459
   551
deba@459
   552
        /// Convert the iterator to the column of the term
deba@459
   553
        operator Col() const {
deba@459
   554
          return colFromId(_it->first);
deba@459
   555
        }
deba@459
   556
deba@459
   557
        /// Returns the coefficient of the term
deba@459
   558
        const Value& operator*() const { return _it->second; }
deba@459
   559
deba@459
   560
        /// Next term
alpar@877
   561
deba@459
   562
        /// Assign the iterator to the next term.
deba@459
   563
        ///
deba@459
   564
        ConstCoeffIt& operator++() { ++_it; return *this; }
deba@459
   565
deba@459
   566
        /// Equality operator
deba@459
   567
        bool operator==(Invalid) const { return _it == _end; }
deba@459
   568
        /// Inequality operator
deba@459
   569
        bool operator!=(Invalid) const { return _it != _end; }
deba@459
   570
      };
deba@459
   571
deba@458
   572
    };
deba@458
   573
deba@458
   574
    ///Linear constraint
deba@458
   575
deba@458
   576
    ///This data stucture represents a linear constraint in the LP.
deba@458
   577
    ///Basically it is a linear expression with a lower or an upper bound
deba@458
   578
    ///(or both). These parts of the constraint can be obtained by the member
deba@458
   579
    ///functions \ref expr(), \ref lowerBound() and \ref upperBound(),
deba@458
   580
    ///respectively.
deba@458
   581
    ///There are two ways to construct a constraint.
deba@458
   582
    ///- You can set the linear expression and the bounds directly
deba@458
   583
    ///  by the functions above.
deba@458
   584
    ///- The operators <tt>\<=</tt>, <tt>==</tt> and  <tt>\>=</tt>
deba@458
   585
    ///  are defined between expressions, or even between constraints whenever
deba@458
   586
    ///  it makes sense. Therefore if \c e and \c f are linear expressions and
deba@458
   587
    ///  \c s and \c t are numbers, then the followings are valid expressions
deba@458
   588
    ///  and thus they can be used directly e.g. in \ref addRow() whenever
deba@458
   589
    ///  it makes sense.
deba@458
   590
    ///\code
deba@458
   591
    ///  e<=s
deba@458
   592
    ///  e<=f
deba@458
   593
    ///  e==f
deba@458
   594
    ///  s<=e<=t
deba@458
   595
    ///  e>=t
deba@458
   596
    ///\endcode
deba@459
   597
    ///\warning The validity of a constraint is checked only at run
deba@459
   598
    ///time, so e.g. \ref addRow(<tt>x[1]\<=x[2]<=5</tt>) will
deba@459
   599
    ///compile, but will fail an assertion.
deba@458
   600
    class Constr
deba@458
   601
    {
deba@458
   602
    public:
deba@459
   603
      typedef LpBase::Expr Expr;
deba@458
   604
      typedef Expr::Key Key;
deba@458
   605
      typedef Expr::Value Value;
deba@458
   606
deba@458
   607
    protected:
deba@458
   608
      Expr _expr;
deba@458
   609
      Value _lb,_ub;
deba@458
   610
    public:
deba@458
   611
      ///\e
deba@458
   612
      Constr() : _expr(), _lb(NaN), _ub(NaN) {}
deba@458
   613
      ///\e
deba@459
   614
      Constr(Value lb, const Expr &e, Value ub) :
deba@458
   615
        _expr(e), _lb(lb), _ub(ub) {}
deba@458
   616
      Constr(const Expr &e) :
deba@458
   617
        _expr(e), _lb(NaN), _ub(NaN) {}
deba@458
   618
      ///\e
deba@458
   619
      void clear()
deba@458
   620
      {
deba@458
   621
        _expr.clear();
deba@458
   622
        _lb=_ub=NaN;
deba@458
   623
      }
deba@458
   624
deba@458
   625
      ///Reference to the linear expression
deba@458
   626
      Expr &expr() { return _expr; }
deba@458
   627
      ///Cont reference to the linear expression
deba@458
   628
      const Expr &expr() const { return _expr; }
deba@458
   629
      ///Reference to the lower bound.
deba@458
   630
deba@458
   631
      ///\return
deba@458
   632
      ///- \ref INF "INF": the constraint is lower unbounded.
deba@458
   633
      ///- \ref NaN "NaN": lower bound has not been set.
deba@458
   634
      ///- finite number: the lower bound
deba@458
   635
      Value &lowerBound() { return _lb; }
deba@458
   636
      ///The const version of \ref lowerBound()
deba@458
   637
      const Value &lowerBound() const { return _lb; }
deba@458
   638
      ///Reference to the upper bound.
deba@458
   639
deba@458
   640
      ///\return
deba@458
   641
      ///- \ref INF "INF": the constraint is upper unbounded.
deba@458
   642
      ///- \ref NaN "NaN": upper bound has not been set.
deba@458
   643
      ///- finite number: the upper bound
deba@458
   644
      Value &upperBound() { return _ub; }
deba@458
   645
      ///The const version of \ref upperBound()
deba@458
   646
      const Value &upperBound() const { return _ub; }
deba@458
   647
      ///Is the constraint lower bounded?
deba@458
   648
      bool lowerBounded() const {
alpar@487
   649
        return _lb != -INF && !isNaN(_lb);
deba@458
   650
      }
deba@458
   651
      ///Is the constraint upper bounded?
deba@458
   652
      bool upperBounded() const {
alpar@487
   653
        return _ub != INF && !isNaN(_ub);
deba@458
   654
      }
deba@458
   655
deba@458
   656
    };
deba@458
   657
deba@458
   658
    ///Linear expression of rows
deba@458
   659
deba@458
   660
    ///This data structure represents a column of the matrix,
deba@458
   661
    ///thas is it strores a linear expression of the dual variables
alpar@1144
   662
    ///(\ref LpBase::Row "Row"s).
deba@458
   663
    ///
deba@458
   664
    ///There are several ways to access and modify the contents of this
deba@458
   665
    ///container.
deba@458
   666
    ///\code
deba@458
   667
    ///e[v]=5;
deba@458
   668
    ///e[v]+=12;
deba@458
   669
    ///e.erase(v);
deba@458
   670
    ///\endcode
deba@458
   671
    ///or you can also iterate through its elements.
deba@458
   672
    ///\code
deba@458
   673
    ///double s=0;
deba@459
   674
    ///for(LpBase::DualExpr::ConstCoeffIt i(e);i!=INVALID;++i)
deba@459
   675
    ///  s+=*i;
deba@458
   676
    ///\endcode
deba@458
   677
    ///(This code computes the sum of all coefficients).
deba@458
   678
    ///- Numbers (<tt>double</tt>'s)
alpar@1144
   679
    ///and variables (\ref LpBase::Row "Row"s) directly convert to an
deba@458
   680
    ///\ref DualExpr and the usual linear operations are defined, so
deba@458
   681
    ///\code
deba@458
   682
    ///v+w
deba@458
   683
    ///2*v-3.12*(v-w/2)
deba@458
   684
    ///v*2.1+(3*v+(v*12+w)*3)/2
deba@458
   685
    ///\endcode
deba@459
   686
    ///are valid \ref DualExpr dual expressions.
deba@458
   687
    ///The usual assignment operations are also defined.
deba@458
   688
    ///\code
deba@458
   689
    ///e=v+w;
deba@458
   690
    ///e+=2*v-3.12*(v-w/2);
deba@458
   691
    ///e*=3.4;
deba@458
   692
    ///e/=5;
deba@458
   693
    ///\endcode
deba@458
   694
    ///
deba@458
   695
    ///\sa Expr
deba@459
   696
    class DualExpr {
deba@459
   697
      friend class LpBase;
deba@458
   698
    public:
deba@459
   699
      /// The key type of the expression
deba@459
   700
      typedef LpBase::Row Key;
deba@459
   701
      /// The value type of the expression
deba@459
   702
      typedef LpBase::Value Value;
deba@458
   703
deba@458
   704
    protected:
deba@459
   705
      std::map<int, Value> comps;
deba@458
   706
deba@458
   707
    public:
deba@459
   708
      typedef True SolverExpr;
deba@459
   709
      /// Default constructor
alpar@877
   710
deba@459
   711
      /// Construct an empty expression, the coefficients are
deba@459
   712
      /// initialized to zero.
deba@459
   713
      DualExpr() {}
deba@459
   714
      /// Construct an expression from a row
deba@459
   715
deba@459
   716
      /// Construct an expression, which has a term with \c r dual
deba@459
   717
      /// variable and 1.0 coefficient.
deba@459
   718
      DualExpr(const Row &r) {
deba@459
   719
        typedef std::map<int, Value>::value_type pair_type;
deba@459
   720
        comps.insert(pair_type(id(r), 1));
deba@458
   721
      }
deba@459
   722
      /// Returns the coefficient of the row
deba@459
   723
      Value operator[](const Row& r) const {
deba@459
   724
        std::map<int, Value>::const_iterator it = comps.find(id(r));
deba@459
   725
        if (it != comps.end()) {
deba@459
   726
          return it->second;
deba@459
   727
        } else {
deba@459
   728
          return 0;
deba@459
   729
        }
deba@458
   730
      }
deba@459
   731
      /// Returns the coefficient of the row
deba@459
   732
      Value& operator[](const Row& r) {
deba@459
   733
        return comps[id(r)];
deba@459
   734
      }
deba@459
   735
      /// Sets the coefficient of the row
deba@459
   736
      void set(const Row &r, const Value &v) {
deba@459
   737
        if (v != 0.0) {
deba@459
   738
          typedef std::map<int, Value>::value_type pair_type;
deba@459
   739
          comps.insert(pair_type(id(r), v));
deba@459
   740
        } else {
deba@459
   741
          comps.erase(id(r));
deba@459
   742
        }
deba@459
   743
      }
deba@459
   744
      /// \brief Removes the coefficients which's absolute value does
alpar@877
   745
      /// not exceed \c epsilon.
deba@459
   746
      void simplify(Value epsilon = 0.0) {
deba@459
   747
        std::map<int, Value>::iterator it=comps.begin();
deba@459
   748
        while (it != comps.end()) {
deba@459
   749
          std::map<int, Value>::iterator jt=it;
deba@459
   750
          ++jt;
deba@459
   751
          if (std::fabs((*it).second) <= epsilon) comps.erase(it);
deba@459
   752
          it=jt;
deba@458
   753
        }
deba@458
   754
      }
deba@458
   755
deba@459
   756
      void simplify(Value epsilon = 0.0) const {
deba@459
   757
        const_cast<DualExpr*>(this)->simplify(epsilon);
deba@458
   758
      }
deba@458
   759
deba@458
   760
      ///Sets all coefficients to 0.
deba@458
   761
      void clear() {
deba@459
   762
        comps.clear();
deba@459
   763
      }
deba@459
   764
      ///Compound assignment
deba@459
   765
      DualExpr &operator+=(const DualExpr &e) {
deba@459
   766
        for (std::map<int, Value>::const_iterator it=e.comps.begin();
deba@459
   767
             it!=e.comps.end(); ++it)
deba@459
   768
          comps[it->first]+=it->second;
deba@459
   769
        return *this;
deba@459
   770
      }
deba@459
   771
      ///Compound assignment
deba@459
   772
      DualExpr &operator-=(const DualExpr &e) {
deba@459
   773
        for (std::map<int, Value>::const_iterator it=e.comps.begin();
deba@459
   774
             it!=e.comps.end(); ++it)
deba@459
   775
          comps[it->first]-=it->second;
deba@459
   776
        return *this;
deba@459
   777
      }
deba@459
   778
      ///Multiply with a constant
deba@459
   779
      DualExpr &operator*=(const Value &v) {
deba@459
   780
        for (std::map<int, Value>::iterator it=comps.begin();
deba@459
   781
             it!=comps.end(); ++it)
deba@459
   782
          it->second*=v;
deba@459
   783
        return *this;
deba@459
   784
      }
deba@459
   785
      ///Division with a constant
deba@459
   786
      DualExpr &operator/=(const Value &v) {
deba@459
   787
        for (std::map<int, Value>::iterator it=comps.begin();
deba@459
   788
             it!=comps.end(); ++it)
deba@459
   789
          it->second/=v;
deba@459
   790
        return *this;
deba@458
   791
      }
deba@458
   792
deba@459
   793
      ///Iterator over the expression
alpar@877
   794
alpar@877
   795
      ///The iterator iterates over the terms of the expression.
alpar@877
   796
      ///
deba@459
   797
      ///\code
deba@459
   798
      ///double s=0;
deba@459
   799
      ///for(LpBase::DualExpr::CoeffIt i(e);i!=INVALID;++i)
deba@459
   800
      ///  s+= *i * dual(i);
deba@459
   801
      ///\endcode
deba@459
   802
      class CoeffIt {
deba@459
   803
      private:
deba@459
   804
deba@459
   805
        std::map<int, Value>::iterator _it, _end;
deba@459
   806
deba@459
   807
      public:
deba@459
   808
deba@459
   809
        /// Sets the iterator to the first term
alpar@877
   810
deba@459
   811
        /// Sets the iterator to the first term of the expression.
deba@459
   812
        ///
deba@459
   813
        CoeffIt(DualExpr& e)
deba@459
   814
          : _it(e.comps.begin()), _end(e.comps.end()){}
deba@459
   815
deba@459
   816
        /// Convert the iterator to the row of the term
deba@459
   817
        operator Row() const {
deba@459
   818
          return rowFromId(_it->first);
deba@459
   819
        }
deba@459
   820
deba@459
   821
        /// Returns the coefficient of the term
deba@459
   822
        Value& operator*() { return _it->second; }
deba@459
   823
deba@459
   824
        /// Returns the coefficient of the term
deba@459
   825
        const Value& operator*() const { return _it->second; }
deba@459
   826
deba@459
   827
        /// Next term
alpar@877
   828
deba@459
   829
        /// Assign the iterator to the next term.
deba@459
   830
        ///
deba@459
   831
        CoeffIt& operator++() { ++_it; return *this; }
deba@459
   832
deba@459
   833
        /// Equality operator
deba@459
   834
        bool operator==(Invalid) const { return _it == _end; }
deba@459
   835
        /// Inequality operator
deba@459
   836
        bool operator!=(Invalid) const { return _it != _end; }
deba@459
   837
      };
deba@459
   838
deba@459
   839
      ///Iterator over the expression
alpar@877
   840
alpar@877
   841
      ///The iterator iterates over the terms of the expression.
alpar@877
   842
      ///
deba@459
   843
      ///\code
deba@459
   844
      ///double s=0;
deba@459
   845
      ///for(LpBase::DualExpr::ConstCoeffIt i(e);i!=INVALID;++i)
deba@459
   846
      ///  s+= *i * dual(i);
deba@459
   847
      ///\endcode
deba@459
   848
      class ConstCoeffIt {
deba@459
   849
      private:
deba@459
   850
deba@459
   851
        std::map<int, Value>::const_iterator _it, _end;
deba@459
   852
deba@459
   853
      public:
deba@459
   854
deba@459
   855
        /// Sets the iterator to the first term
alpar@877
   856
deba@459
   857
        /// Sets the iterator to the first term of the expression.
deba@459
   858
        ///
deba@459
   859
        ConstCoeffIt(const DualExpr& e)
deba@459
   860
          : _it(e.comps.begin()), _end(e.comps.end()){}
deba@459
   861
deba@459
   862
        /// Convert the iterator to the row of the term
deba@459
   863
        operator Row() const {
deba@459
   864
          return rowFromId(_it->first);
deba@459
   865
        }
deba@459
   866
deba@459
   867
        /// Returns the coefficient of the term
deba@459
   868
        const Value& operator*() const { return _it->second; }
deba@459
   869
deba@459
   870
        /// Next term
alpar@877
   871
deba@459
   872
        /// Assign the iterator to the next term.
deba@459
   873
        ///
deba@459
   874
        ConstCoeffIt& operator++() { ++_it; return *this; }
deba@459
   875
deba@459
   876
        /// Equality operator
deba@459
   877
        bool operator==(Invalid) const { return _it == _end; }
deba@459
   878
        /// Inequality operator
deba@459
   879
        bool operator!=(Invalid) const { return _it != _end; }
deba@459
   880
      };
deba@458
   881
    };
deba@458
   882
deba@458
   883
deba@459
   884
  protected:
deba@458
   885
deba@459
   886
    class InsertIterator {
deba@459
   887
    private:
deba@459
   888
deba@459
   889
      std::map<int, Value>& _host;
deba@459
   890
      const _solver_bits::VarIndex& _index;
deba@459
   891
deba@458
   892
    public:
deba@458
   893
deba@458
   894
      typedef std::output_iterator_tag iterator_category;
deba@458
   895
      typedef void difference_type;
deba@458
   896
      typedef void value_type;
deba@458
   897
      typedef void reference;
deba@458
   898
      typedef void pointer;
deba@458
   899
deba@459
   900
      InsertIterator(std::map<int, Value>& host,
deba@459
   901
                   const _solver_bits::VarIndex& index)
deba@459
   902
        : _host(host), _index(index) {}
deba@458
   903
deba@459
   904
      InsertIterator& operator=(const std::pair<int, Value>& value) {
deba@459
   905
        typedef std::map<int, Value>::value_type pair_type;
deba@459
   906
        _host.insert(pair_type(_index[value.first], value.second));
deba@458
   907
        return *this;
deba@458
   908
      }
deba@458
   909
deba@459
   910
      InsertIterator& operator*() { return *this; }
deba@459
   911
      InsertIterator& operator++() { return *this; }
deba@459
   912
      InsertIterator operator++(int) { return *this; }
deba@458
   913
deba@458
   914
    };
deba@458
   915
deba@459
   916
    class ExprIterator {
deba@459
   917
    private:
deba@459
   918
      std::map<int, Value>::const_iterator _host_it;
deba@459
   919
      const _solver_bits::VarIndex& _index;
deba@458
   920
    public:
deba@458
   921
deba@459
   922
      typedef std::bidirectional_iterator_tag iterator_category;
deba@459
   923
      typedef std::ptrdiff_t difference_type;
deba@458
   924
      typedef const std::pair<int, Value> value_type;
deba@458
   925
      typedef value_type reference;
deba@459
   926
deba@458
   927
      class pointer {
deba@458
   928
      public:
deba@458
   929
        pointer(value_type& _value) : value(_value) {}
deba@458
   930
        value_type* operator->() { return &value; }
deba@458
   931
      private:
deba@458
   932
        value_type value;
deba@458
   933
      };
deba@458
   934
deba@459
   935
      ExprIterator(const std::map<int, Value>::const_iterator& host_it,
deba@459
   936
                   const _solver_bits::VarIndex& index)
deba@459
   937
        : _host_it(host_it), _index(index) {}
deba@458
   938
deba@458
   939
      reference operator*() {
deba@459
   940
        return std::make_pair(_index(_host_it->first), _host_it->second);
deba@458
   941
      }
deba@458
   942
deba@458
   943
      pointer operator->() {
deba@458
   944
        return pointer(operator*());
deba@458
   945
      }
deba@458
   946
deba@459
   947
      ExprIterator& operator++() { ++_host_it; return *this; }
deba@459
   948
      ExprIterator operator++(int) {
deba@459
   949
        ExprIterator tmp(*this); ++_host_it; return tmp;
deba@458
   950
      }
deba@458
   951
deba@459
   952
      ExprIterator& operator--() { --_host_it; return *this; }
deba@459
   953
      ExprIterator operator--(int) {
deba@459
   954
        ExprIterator tmp(*this); --_host_it; return tmp;
deba@458
   955
      }
deba@458
   956
deba@459
   957
      bool operator==(const ExprIterator& it) const {
deba@459
   958
        return _host_it == it._host_it;
deba@458
   959
      }
deba@458
   960
deba@459
   961
      bool operator!=(const ExprIterator& it) const {
deba@459
   962
        return _host_it != it._host_it;
deba@458
   963
      }
deba@458
   964
deba@458
   965
    };
deba@458
   966
deba@458
   967
  protected:
deba@458
   968
deba@459
   969
    //Abstract virtual functions
deba@458
   970
ggab90@1130
   971
    virtual int _addColId(int col) { return _cols.addIndex(col); }
ggab90@1130
   972
    virtual int _addRowId(int row) { return _rows.addIndex(row); }
deba@458
   973
ggab90@1130
   974
    virtual void _eraseColId(int col) { _cols.eraseIndex(col); }
ggab90@1130
   975
    virtual void _eraseRowId(int row) { _rows.eraseIndex(row); }
deba@458
   976
deba@458
   977
    virtual int _addCol() = 0;
deba@458
   978
    virtual int _addRow() = 0;
deba@458
   979
deba@746
   980
    virtual int _addRow(Value l, ExprIterator b, ExprIterator e, Value u) {
deba@746
   981
      int row = _addRow();
deba@746
   982
      _setRowCoeffs(row, b, e);
deba@746
   983
      _setRowLowerBound(row, l);
deba@746
   984
      _setRowUpperBound(row, u);
deba@746
   985
      return row;
deba@746
   986
    }
deba@746
   987
deba@458
   988
    virtual void _eraseCol(int col) = 0;
deba@458
   989
    virtual void _eraseRow(int row) = 0;
deba@458
   990
deba@459
   991
    virtual void _getColName(int col, std::string& name) const = 0;
deba@459
   992
    virtual void _setColName(int col, const std::string& name) = 0;
deba@458
   993
    virtual int _colByName(const std::string& name) const = 0;
deba@458
   994
deba@459
   995
    virtual void _getRowName(int row, std::string& name) const = 0;
deba@459
   996
    virtual void _setRowName(int row, const std::string& name) = 0;
deba@459
   997
    virtual int _rowByName(const std::string& name) const = 0;
deba@459
   998
deba@459
   999
    virtual void _setRowCoeffs(int i, ExprIterator b, ExprIterator e) = 0;
deba@459
  1000
    virtual void _getRowCoeffs(int i, InsertIterator b) const = 0;
deba@459
  1001
deba@459
  1002
    virtual void _setColCoeffs(int i, ExprIterator b, ExprIterator e) = 0;
deba@459
  1003
    virtual void _getColCoeffs(int i, InsertIterator b) const = 0;
deba@459
  1004
deba@458
  1005
    virtual void _setCoeff(int row, int col, Value value) = 0;
deba@458
  1006
    virtual Value _getCoeff(int row, int col) const = 0;
deba@459
  1007
deba@458
  1008
    virtual void _setColLowerBound(int i, Value value) = 0;
deba@458
  1009
    virtual Value _getColLowerBound(int i) const = 0;
deba@459
  1010
deba@458
  1011
    virtual void _setColUpperBound(int i, Value value) = 0;
deba@458
  1012
    virtual Value _getColUpperBound(int i) const = 0;
deba@459
  1013
deba@459
  1014
    virtual void _setRowLowerBound(int i, Value value) = 0;
deba@459
  1015
    virtual Value _getRowLowerBound(int i) const = 0;
deba@459
  1016
deba@459
  1017
    virtual void _setRowUpperBound(int i, Value value) = 0;
deba@459
  1018
    virtual Value _getRowUpperBound(int i) const = 0;
deba@459
  1019
deba@459
  1020
    virtual void _setObjCoeffs(ExprIterator b, ExprIterator e) = 0;
deba@459
  1021
    virtual void _getObjCoeffs(InsertIterator b) const = 0;
deba@458
  1022
deba@458
  1023
    virtual void _setObjCoeff(int i, Value obj_coef) = 0;
deba@458
  1024
    virtual Value _getObjCoeff(int i) const = 0;
deba@458
  1025
deba@459
  1026
    virtual void _setSense(Sense) = 0;
deba@459
  1027
    virtual Sense _getSense() const = 0;
deba@458
  1028
deba@459
  1029
    virtual void _clear() = 0;
deba@458
  1030
deba@459
  1031
    virtual const char* _solverName() const = 0;
deba@458
  1032
deba@576
  1033
    virtual void _messageLevel(MessageLevel level) = 0;
deba@576
  1034
deba@458
  1035
    //Own protected stuff
deba@458
  1036
deba@458
  1037
    //Constant component of the objective function
deba@458
  1038
    Value obj_const_comp;
deba@458
  1039
ggab90@1130
  1040
    LpBase() : _rows(), _cols(), obj_const_comp(0) {}
deba@459
  1041
deba@458
  1042
  public:
deba@458
  1043
alpar@1076
  1044
    ///Unsupported file format exception
alpar@1063
  1045
    class UnsupportedFormatError : public Exception
alpar@1063
  1046
    {
alpar@1063
  1047
      std::string _format;
alpar@1063
  1048
      mutable std::string _what;
alpar@1063
  1049
    public:
alpar@1063
  1050
      explicit UnsupportedFormatError(std::string format) throw()
alpar@1063
  1051
        : _format(format) { }
alpar@1063
  1052
      virtual ~UnsupportedFormatError() throw() {}
alpar@1063
  1053
      virtual const char* what() const throw() {
alpar@1063
  1054
        try {
alpar@1063
  1055
          _what.clear();
alpar@1063
  1056
          std::ostringstream oss;
alpar@1063
  1057
          oss << "lemon::UnsupportedFormatError: " << _format;
alpar@1063
  1058
          _what = oss.str();
alpar@1063
  1059
        }
alpar@1063
  1060
        catch (...) {}
alpar@1063
  1061
        if (!_what.empty()) return _what.c_str();
alpar@1063
  1062
        else return "lemon::UnsupportedFormatError";
alpar@1063
  1063
      }
alpar@1063
  1064
    };
alpar@1092
  1065
alpar@1063
  1066
  protected:
alpar@1063
  1067
    virtual void _write(std::string, std::string format) const
alpar@1063
  1068
    {
alpar@1063
  1069
      throw UnsupportedFormatError(format);
alpar@1063
  1070
    }
alpar@1092
  1071
alpar@1063
  1072
  public:
alpar@1063
  1073
deba@459
  1074
    /// Virtual destructor
deba@459
  1075
    virtual ~LpBase() {}
deba@458
  1076
deba@459
  1077
    ///Gives back the name of the solver.
deba@459
  1078
    const char* solverName() const {return _solverName();}
deba@458
  1079
kpeter@584
  1080
    ///\name Build Up and Modify the LP
deba@458
  1081
deba@458
  1082
    ///@{
deba@458
  1083
deba@458
  1084
    ///Add a new empty column (i.e a new variable) to the LP
deba@459
  1085
    Col addCol() { Col c; c._id = _addColId(_addCol()); return c;}
deba@458
  1086
deba@459
  1087
    ///\brief Adds several new columns (i.e variables) at once
deba@458
  1088
    ///
deba@459
  1089
    ///This magic function takes a container as its argument and fills
deba@459
  1090
    ///its elements with new columns (i.e. variables)
deba@458
  1091
    ///\param t can be
deba@458
  1092
    ///- a standard STL compatible iterable container with
deba@459
  1093
    ///\ref Col as its \c values_type like
deba@458
  1094
    ///\code
deba@459
  1095
    ///std::vector<LpBase::Col>
deba@459
  1096
    ///std::list<LpBase::Col>
deba@458
  1097
    ///\endcode
deba@458
  1098
    ///- a standard STL compatible iterable container with
deba@459
  1099
    ///\ref Col as its \c mapped_type like
deba@458
  1100
    ///\code
deba@459
  1101
    ///std::map<AnyType,LpBase::Col>
deba@458
  1102
    ///\endcode
deba@458
  1103
    ///- an iterable lemon \ref concepts::WriteMap "write map" like
deba@458
  1104
    ///\code
deba@459
  1105
    ///ListGraph::NodeMap<LpBase::Col>
deba@459
  1106
    ///ListGraph::ArcMap<LpBase::Col>
deba@458
  1107
    ///\endcode
deba@458
  1108
    ///\return The number of the created column.
deba@458
  1109
#ifdef DOXYGEN
deba@458
  1110
    template<class T>
deba@458
  1111
    int addColSet(T &t) { return 0;}
deba@458
  1112
#else
deba@458
  1113
    template<class T>
deba@459
  1114
    typename enable_if<typename T::value_type::LpCol,int>::type
deba@458
  1115
    addColSet(T &t,dummy<0> = 0) {
deba@458
  1116
      int s=0;
deba@458
  1117
      for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addCol();s++;}
deba@458
  1118
      return s;
deba@458
  1119
    }
deba@458
  1120
    template<class T>
deba@459
  1121
    typename enable_if<typename T::value_type::second_type::LpCol,
deba@458
  1122
                       int>::type
deba@458
  1123
    addColSet(T &t,dummy<1> = 1) {
deba@458
  1124
      int s=0;
deba@458
  1125
      for(typename T::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1126
        i->second=addCol();
deba@458
  1127
        s++;
deba@458
  1128
      }
deba@458
  1129
      return s;
deba@458
  1130
    }
deba@458
  1131
    template<class T>
deba@459
  1132
    typename enable_if<typename T::MapIt::Value::LpCol,
deba@458
  1133
                       int>::type
deba@458
  1134
    addColSet(T &t,dummy<2> = 2) {
deba@458
  1135
      int s=0;
deba@458
  1136
      for(typename T::MapIt i(t); i!=INVALID; ++i)
deba@458
  1137
        {
deba@458
  1138
          i.set(addCol());
deba@458
  1139
          s++;
deba@458
  1140
        }
deba@458
  1141
      return s;
deba@458
  1142
    }
deba@458
  1143
#endif
deba@458
  1144
deba@458
  1145
    ///Set a column (i.e a dual constraint) of the LP
deba@458
  1146
deba@458
  1147
    ///\param c is the column to be modified
deba@458
  1148
    ///\param e is a dual linear expression (see \ref DualExpr)
deba@458
  1149
    ///a better one.
deba@459
  1150
    void col(Col c, const DualExpr &e) {
deba@458
  1151
      e.simplify();
ggab90@1130
  1152
      _setColCoeffs(_cols(id(c)), ExprIterator(e.comps.begin(), _rows),
ggab90@1130
  1153
                    ExprIterator(e.comps.end(), _rows));
deba@458
  1154
    }
deba@458
  1155
deba@458
  1156
    ///Get a column (i.e a dual constraint) of the LP
deba@458
  1157
deba@459
  1158
    ///\param c is the column to get
deba@458
  1159
    ///\return the dual expression associated to the column
deba@458
  1160
    DualExpr col(Col c) const {
deba@458
  1161
      DualExpr e;
ggab90@1130
  1162
      _getColCoeffs(_cols(id(c)), InsertIterator(e.comps, _rows));
deba@458
  1163
      return e;
deba@458
  1164
    }
deba@458
  1165
deba@458
  1166
    ///Add a new column to the LP
deba@458
  1167
deba@458
  1168
    ///\param e is a dual linear expression (see \ref DualExpr)
deba@459
  1169
    ///\param o is the corresponding component of the objective
deba@458
  1170
    ///function. It is 0 by default.
deba@458
  1171
    ///\return The created column.
deba@458
  1172
    Col addCol(const DualExpr &e, Value o = 0) {
deba@458
  1173
      Col c=addCol();
deba@458
  1174
      col(c,e);
deba@458
  1175
      objCoeff(c,o);
deba@458
  1176
      return c;
deba@458
  1177
    }
deba@458
  1178
deba@458
  1179
    ///Add a new empty row (i.e a new constraint) to the LP
deba@458
  1180
deba@458
  1181
    ///This function adds a new empty row (i.e a new constraint) to the LP.
deba@458
  1182
    ///\return The created row
deba@459
  1183
    Row addRow() { Row r; r._id = _addRowId(_addRow()); return r;}
deba@458
  1184
deba@459
  1185
    ///\brief Add several new rows (i.e constraints) at once
deba@458
  1186
    ///
deba@459
  1187
    ///This magic function takes a container as its argument and fills
deba@459
  1188
    ///its elements with new row (i.e. variables)
deba@458
  1189
    ///\param t can be
deba@458
  1190
    ///- a standard STL compatible iterable container with
alpar@1144
  1191
    ///\ref LpBase::Row "Row" as its \c values_type like
deba@458
  1192
    ///\code
deba@459
  1193
    ///std::vector<LpBase::Row>
deba@459
  1194
    ///std::list<LpBase::Row>
deba@458
  1195
    ///\endcode
deba@458
  1196
    ///- a standard STL compatible iterable container with
alpar@1144
  1197
    ///\ref LpBase::Row "Row" as its \c mapped_type like
deba@458
  1198
    ///\code
deba@459
  1199
    ///std::map<AnyType,LpBase::Row>
deba@458
  1200
    ///\endcode
deba@458
  1201
    ///- an iterable lemon \ref concepts::WriteMap "write map" like
deba@458
  1202
    ///\code
deba@459
  1203
    ///ListGraph::NodeMap<LpBase::Row>
deba@459
  1204
    ///ListGraph::ArcMap<LpBase::Row>
deba@458
  1205
    ///\endcode
deba@458
  1206
    ///\return The number of rows created.
deba@458
  1207
#ifdef DOXYGEN
deba@458
  1208
    template<class T>
deba@458
  1209
    int addRowSet(T &t) { return 0;}
deba@458
  1210
#else
deba@458
  1211
    template<class T>
deba@459
  1212
    typename enable_if<typename T::value_type::LpRow,int>::type
deba@459
  1213
    addRowSet(T &t, dummy<0> = 0) {
deba@458
  1214
      int s=0;
deba@458
  1215
      for(typename T::iterator i=t.begin();i!=t.end();++i) {*i=addRow();s++;}
deba@458
  1216
      return s;
deba@458
  1217
    }
deba@458
  1218
    template<class T>
deba@459
  1219
    typename enable_if<typename T::value_type::second_type::LpRow, int>::type
deba@459
  1220
    addRowSet(T &t, dummy<1> = 1) {
deba@458
  1221
      int s=0;
deba@458
  1222
      for(typename T::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1223
        i->second=addRow();
deba@458
  1224
        s++;
deba@458
  1225
      }
deba@458
  1226
      return s;
deba@458
  1227
    }
deba@458
  1228
    template<class T>
deba@459
  1229
    typename enable_if<typename T::MapIt::Value::LpRow, int>::type
deba@459
  1230
    addRowSet(T &t, dummy<2> = 2) {
deba@458
  1231
      int s=0;
deba@458
  1232
      for(typename T::MapIt i(t); i!=INVALID; ++i)
deba@458
  1233
        {
deba@458
  1234
          i.set(addRow());
deba@458
  1235
          s++;
deba@458
  1236
        }
deba@458
  1237
      return s;
deba@458
  1238
    }
deba@458
  1239
#endif
deba@458
  1240
deba@458
  1241
    ///Set a row (i.e a constraint) of the LP
deba@458
  1242
deba@458
  1243
    ///\param r is the row to be modified
deba@458
  1244
    ///\param l is lower bound (-\ref INF means no bound)
deba@458
  1245
    ///\param e is a linear expression (see \ref Expr)
deba@458
  1246
    ///\param u is the upper bound (\ref INF means no bound)
deba@458
  1247
    void row(Row r, Value l, const Expr &e, Value u) {
deba@458
  1248
      e.simplify();
ggab90@1130
  1249
      _setRowCoeffs(_rows(id(r)), ExprIterator(e.comps.begin(), _cols),
ggab90@1130
  1250
                    ExprIterator(e.comps.end(), _cols));
ggab90@1130
  1251
      _setRowLowerBound(_rows(id(r)),l - *e);
ggab90@1130
  1252
      _setRowUpperBound(_rows(id(r)),u - *e);
deba@458
  1253
    }
deba@458
  1254
deba@458
  1255
    ///Set a row (i.e a constraint) of the LP
deba@458
  1256
deba@458
  1257
    ///\param r is the row to be modified
deba@458
  1258
    ///\param c is a linear expression (see \ref Constr)
deba@458
  1259
    void row(Row r, const Constr &c) {
deba@458
  1260
      row(r, c.lowerBounded()?c.lowerBound():-INF,
deba@458
  1261
          c.expr(), c.upperBounded()?c.upperBound():INF);
deba@458
  1262
    }
deba@458
  1263
deba@458
  1264
deba@458
  1265
    ///Get a row (i.e a constraint) of the LP
deba@458
  1266
deba@458
  1267
    ///\param r is the row to get
deba@458
  1268
    ///\return the expression associated to the row
deba@458
  1269
    Expr row(Row r) const {
deba@458
  1270
      Expr e;
ggab90@1130
  1271
      _getRowCoeffs(_rows(id(r)), InsertIterator(e.comps, _cols));
deba@458
  1272
      return e;
deba@458
  1273
    }
deba@458
  1274
deba@458
  1275
    ///Add a new row (i.e a new constraint) to the LP
deba@458
  1276
deba@458
  1277
    ///\param l is the lower bound (-\ref INF means no bound)
deba@458
  1278
    ///\param e is a linear expression (see \ref Expr)
deba@458
  1279
    ///\param u is the upper bound (\ref INF means no bound)
deba@458
  1280
    ///\return The created row.
deba@458
  1281
    Row addRow(Value l,const Expr &e, Value u) {
deba@746
  1282
      Row r;
deba@746
  1283
      e.simplify();
ggab90@1130
  1284
      r._id = _addRowId(_addRow(l - *e, ExprIterator(e.comps.begin(), _cols),
ggab90@1130
  1285
                                ExprIterator(e.comps.end(), _cols), u - *e));
deba@458
  1286
      return r;
deba@458
  1287
    }
deba@458
  1288
deba@458
  1289
    ///Add a new row (i.e a new constraint) to the LP
deba@458
  1290
deba@458
  1291
    ///\param c is a linear expression (see \ref Constr)
deba@458
  1292
    ///\return The created row.
deba@458
  1293
    Row addRow(const Constr &c) {
deba@746
  1294
      Row r;
deba@746
  1295
      c.expr().simplify();
alpar@877
  1296
      r._id = _addRowId(_addRow(c.lowerBounded()?c.lowerBound()-*c.expr():-INF,
ggab90@1130
  1297
                                ExprIterator(c.expr().comps.begin(), _cols),
ggab90@1130
  1298
                                ExprIterator(c.expr().comps.end(), _cols),
deba@834
  1299
                                c.upperBounded()?c.upperBound()-*c.expr():INF));
deba@458
  1300
      return r;
deba@458
  1301
    }
deba@459
  1302
    ///Erase a column (i.e a variable) from the LP
deba@458
  1303
deba@459
  1304
    ///\param c is the column to be deleted
deba@459
  1305
    void erase(Col c) {
ggab90@1130
  1306
      _eraseCol(_cols(id(c)));
ggab90@1130
  1307
      _eraseColId(_cols(id(c)));
deba@458
  1308
    }
deba@459
  1309
    ///Erase a row (i.e a constraint) from the LP
deba@458
  1310
deba@458
  1311
    ///\param r is the row to be deleted
deba@459
  1312
    void erase(Row r) {
ggab90@1130
  1313
      _eraseRow(_rows(id(r)));
ggab90@1130
  1314
      _eraseRowId(_rows(id(r)));
deba@458
  1315
    }
deba@458
  1316
deba@458
  1317
    /// Get the name of a column
deba@458
  1318
deba@459
  1319
    ///\param c is the coresponding column
deba@458
  1320
    ///\return The name of the colunm
deba@458
  1321
    std::string colName(Col c) const {
deba@458
  1322
      std::string name;
ggab90@1130
  1323
      _getColName(_cols(id(c)), name);
deba@458
  1324
      return name;
deba@458
  1325
    }
deba@458
  1326
deba@458
  1327
    /// Set the name of a column
deba@458
  1328
deba@459
  1329
    ///\param c is the coresponding column
deba@458
  1330
    ///\param name The name to be given
deba@458
  1331
    void colName(Col c, const std::string& name) {
ggab90@1130
  1332
      _setColName(_cols(id(c)), name);
deba@458
  1333
    }
deba@458
  1334
deba@458
  1335
    /// Get the column by its name
deba@458
  1336
deba@458
  1337
    ///\param name The name of the column
deba@458
  1338
    ///\return the proper column or \c INVALID
deba@458
  1339
    Col colByName(const std::string& name) const {
deba@458
  1340
      int k = _colByName(name);
ggab90@1130
  1341
      return k != -1 ? Col(_cols[k]) : Col(INVALID);
deba@459
  1342
    }
deba@459
  1343
deba@459
  1344
    /// Get the name of a row
deba@459
  1345
deba@459
  1346
    ///\param r is the coresponding row
deba@459
  1347
    ///\return The name of the row
deba@459
  1348
    std::string rowName(Row r) const {
deba@459
  1349
      std::string name;
ggab90@1130
  1350
      _getRowName(_rows(id(r)), name);
deba@459
  1351
      return name;
deba@459
  1352
    }
deba@459
  1353
deba@459
  1354
    /// Set the name of a row
deba@459
  1355
deba@459
  1356
    ///\param r is the coresponding row
deba@459
  1357
    ///\param name The name to be given
deba@459
  1358
    void rowName(Row r, const std::string& name) {
ggab90@1130
  1359
      _setRowName(_rows(id(r)), name);
deba@459
  1360
    }
deba@459
  1361
deba@459
  1362
    /// Get the row by its name
deba@459
  1363
deba@459
  1364
    ///\param name The name of the row
deba@459
  1365
    ///\return the proper row or \c INVALID
deba@459
  1366
    Row rowByName(const std::string& name) const {
deba@459
  1367
      int k = _rowByName(name);
ggab90@1130
  1368
      return k != -1 ? Row(_rows[k]) : Row(INVALID);
deba@458
  1369
    }
deba@458
  1370
deba@458
  1371
    /// Set an element of the coefficient matrix of the LP
deba@458
  1372
deba@458
  1373
    ///\param r is the row of the element to be modified
deba@459
  1374
    ///\param c is the column of the element to be modified
deba@458
  1375
    ///\param val is the new value of the coefficient
deba@458
  1376
    void coeff(Row r, Col c, Value val) {
ggab90@1130
  1377
      _setCoeff(_rows(id(r)),_cols(id(c)), val);
deba@458
  1378
    }
deba@458
  1379
deba@458
  1380
    /// Get an element of the coefficient matrix of the LP
deba@458
  1381
deba@459
  1382
    ///\param r is the row of the element
deba@459
  1383
    ///\param c is the column of the element
deba@458
  1384
    ///\return the corresponding coefficient
deba@458
  1385
    Value coeff(Row r, Col c) const {
ggab90@1130
  1386
      return _getCoeff(_rows(id(r)),_cols(id(c)));
deba@458
  1387
    }
deba@458
  1388
deba@458
  1389
    /// Set the lower bound of a column (i.e a variable)
deba@458
  1390
deba@458
  1391
    /// The lower bound of a variable (column) has to be given by an
deba@458
  1392
    /// extended number of type Value, i.e. a finite number of type
deba@458
  1393
    /// Value or -\ref INF.
deba@458
  1394
    void colLowerBound(Col c, Value value) {
ggab90@1130
  1395
      _setColLowerBound(_cols(id(c)),value);
deba@458
  1396
    }
deba@458
  1397
deba@458
  1398
    /// Get the lower bound of a column (i.e a variable)
deba@458
  1399
deba@459
  1400
    /// This function returns the lower bound for column (variable) \c c
deba@458
  1401
    /// (this might be -\ref INF as well).
deba@459
  1402
    ///\return The lower bound for column \c c
deba@458
  1403
    Value colLowerBound(Col c) const {
ggab90@1130
  1404
      return _getColLowerBound(_cols(id(c)));
deba@458
  1405
    }
deba@458
  1406
deba@458
  1407
    ///\brief Set the lower bound of  several columns
deba@459
  1408
    ///(i.e variables) at once
deba@458
  1409
    ///
deba@458
  1410
    ///This magic function takes a container as its argument
deba@458
  1411
    ///and applies the function on all of its elements.
deba@459
  1412
    ///The lower bound of a variable (column) has to be given by an
deba@459
  1413
    ///extended number of type Value, i.e. a finite number of type
deba@459
  1414
    ///Value or -\ref INF.
deba@458
  1415
#ifdef DOXYGEN
deba@458
  1416
    template<class T>
deba@458
  1417
    void colLowerBound(T &t, Value value) { return 0;}
deba@458
  1418
#else
deba@458
  1419
    template<class T>
deba@459
  1420
    typename enable_if<typename T::value_type::LpCol,void>::type
deba@458
  1421
    colLowerBound(T &t, Value value,dummy<0> = 0) {
deba@458
  1422
      for(typename T::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1423
        colLowerBound(*i, value);
deba@458
  1424
      }
deba@458
  1425
    }
deba@458
  1426
    template<class T>
deba@459
  1427
    typename enable_if<typename T::value_type::second_type::LpCol,
deba@458
  1428
                       void>::type
deba@458
  1429
    colLowerBound(T &t, Value value,dummy<1> = 1) {
deba@458
  1430
      for(typename T::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1431
        colLowerBound(i->second, value);
deba@458
  1432
      }
deba@458
  1433
    }
deba@458
  1434
    template<class T>
deba@459
  1435
    typename enable_if<typename T::MapIt::Value::LpCol,
deba@458
  1436
                       void>::type
deba@458
  1437
    colLowerBound(T &t, Value value,dummy<2> = 2) {
deba@458
  1438
      for(typename T::MapIt i(t); i!=INVALID; ++i){
deba@458
  1439
        colLowerBound(*i, value);
deba@458
  1440
      }
deba@458
  1441
    }
deba@458
  1442
#endif
deba@458
  1443
deba@458
  1444
    /// Set the upper bound of a column (i.e a variable)
deba@458
  1445
deba@458
  1446
    /// The upper bound of a variable (column) has to be given by an
deba@458
  1447
    /// extended number of type Value, i.e. a finite number of type
deba@458
  1448
    /// Value or \ref INF.
deba@458
  1449
    void colUpperBound(Col c, Value value) {
ggab90@1130
  1450
      _setColUpperBound(_cols(id(c)),value);
deba@458
  1451
    };
deba@458
  1452
deba@458
  1453
    /// Get the upper bound of a column (i.e a variable)
deba@458
  1454
deba@459
  1455
    /// This function returns the upper bound for column (variable) \c c
deba@458
  1456
    /// (this might be \ref INF as well).
deba@459
  1457
    /// \return The upper bound for column \c c
deba@458
  1458
    Value colUpperBound(Col c) const {
ggab90@1130
  1459
      return _getColUpperBound(_cols(id(c)));
deba@458
  1460
    }
deba@458
  1461
deba@458
  1462
    ///\brief Set the upper bound of  several columns
deba@459
  1463
    ///(i.e variables) at once
deba@458
  1464
    ///
deba@458
  1465
    ///This magic function takes a container as its argument
deba@458
  1466
    ///and applies the function on all of its elements.
deba@459
  1467
    ///The upper bound of a variable (column) has to be given by an
deba@459
  1468
    ///extended number of type Value, i.e. a finite number of type
deba@459
  1469
    ///Value or \ref INF.
deba@458
  1470
#ifdef DOXYGEN
deba@458
  1471
    template<class T>
deba@458
  1472
    void colUpperBound(T &t, Value value) { return 0;}
deba@458
  1473
#else
tapolcai@490
  1474
    template<class T1>
tapolcai@490
  1475
    typename enable_if<typename T1::value_type::LpCol,void>::type
tapolcai@490
  1476
    colUpperBound(T1 &t, Value value,dummy<0> = 0) {
tapolcai@490
  1477
      for(typename T1::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1478
        colUpperBound(*i, value);
deba@458
  1479
      }
deba@458
  1480
    }
tapolcai@490
  1481
    template<class T1>
tapolcai@490
  1482
    typename enable_if<typename T1::value_type::second_type::LpCol,
deba@458
  1483
                       void>::type
tapolcai@490
  1484
    colUpperBound(T1 &t, Value value,dummy<1> = 1) {
tapolcai@490
  1485
      for(typename T1::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1486
        colUpperBound(i->second, value);
deba@458
  1487
      }
deba@458
  1488
    }
tapolcai@490
  1489
    template<class T1>
tapolcai@490
  1490
    typename enable_if<typename T1::MapIt::Value::LpCol,
deba@458
  1491
                       void>::type
tapolcai@490
  1492
    colUpperBound(T1 &t, Value value,dummy<2> = 2) {
tapolcai@490
  1493
      for(typename T1::MapIt i(t); i!=INVALID; ++i){
deba@458
  1494
        colUpperBound(*i, value);
deba@458
  1495
      }
deba@458
  1496
    }
deba@458
  1497
#endif
deba@458
  1498
deba@458
  1499
    /// Set the lower and the upper bounds of a column (i.e a variable)
deba@458
  1500
deba@458
  1501
    /// The lower and the upper bounds of
deba@458
  1502
    /// a variable (column) have to be given by an
deba@458
  1503
    /// extended number of type Value, i.e. a finite number of type
deba@458
  1504
    /// Value, -\ref INF or \ref INF.
deba@458
  1505
    void colBounds(Col c, Value lower, Value upper) {
ggab90@1130
  1506
      _setColLowerBound(_cols(id(c)),lower);
ggab90@1130
  1507
      _setColUpperBound(_cols(id(c)),upper);
deba@458
  1508
    }
deba@458
  1509
deba@458
  1510
    ///\brief Set the lower and the upper bound of several columns
deba@459
  1511
    ///(i.e variables) at once
deba@458
  1512
    ///
deba@458
  1513
    ///This magic function takes a container as its argument
deba@458
  1514
    ///and applies the function on all of its elements.
deba@458
  1515
    /// The lower and the upper bounds of
deba@458
  1516
    /// a variable (column) have to be given by an
deba@458
  1517
    /// extended number of type Value, i.e. a finite number of type
deba@458
  1518
    /// Value, -\ref INF or \ref INF.
deba@458
  1519
#ifdef DOXYGEN
deba@458
  1520
    template<class T>
deba@458
  1521
    void colBounds(T &t, Value lower, Value upper) { return 0;}
deba@458
  1522
#else
tapolcai@490
  1523
    template<class T2>
tapolcai@490
  1524
    typename enable_if<typename T2::value_type::LpCol,void>::type
tapolcai@490
  1525
    colBounds(T2 &t, Value lower, Value upper,dummy<0> = 0) {
tapolcai@490
  1526
      for(typename T2::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1527
        colBounds(*i, lower, upper);
deba@458
  1528
      }
deba@458
  1529
    }
tapolcai@490
  1530
    template<class T2>
tapolcai@490
  1531
    typename enable_if<typename T2::value_type::second_type::LpCol, void>::type
tapolcai@490
  1532
    colBounds(T2 &t, Value lower, Value upper,dummy<1> = 1) {
tapolcai@490
  1533
      for(typename T2::iterator i=t.begin();i!=t.end();++i) {
deba@458
  1534
        colBounds(i->second, lower, upper);
deba@458
  1535
      }
deba@458
  1536
    }
tapolcai@490
  1537
    template<class T2>
tapolcai@490
  1538
    typename enable_if<typename T2::MapIt::Value::LpCol, void>::type
tapolcai@490
  1539
    colBounds(T2 &t, Value lower, Value upper,dummy<2> = 2) {
tapolcai@490
  1540
      for(typename T2::MapIt i(t); i!=INVALID; ++i){
deba@458
  1541
        colBounds(*i, lower, upper);
deba@458
  1542
      }
deba@458
  1543
    }
deba@458
  1544
#endif
deba@458
  1545
deba@459
  1546
    /// Set the lower bound of a row (i.e a constraint)
deba@458
  1547
deba@459
  1548
    /// The lower bound of a constraint (row) has to be given by an
deba@459
  1549
    /// extended number of type Value, i.e. a finite number of type
deba@459
  1550
    /// Value or -\ref INF.
deba@459
  1551
    void rowLowerBound(Row r, Value value) {
ggab90@1130
  1552
      _setRowLowerBound(_rows(id(r)),value);
deba@458
  1553
    }
deba@458
  1554
deba@459
  1555
    /// Get the lower bound of a row (i.e a constraint)
deba@458
  1556
deba@459
  1557
    /// This function returns the lower bound for row (constraint) \c c
deba@459
  1558
    /// (this might be -\ref INF as well).
deba@459
  1559
    ///\return The lower bound for row \c r
deba@459
  1560
    Value rowLowerBound(Row r) const {
ggab90@1130
  1561
      return _getRowLowerBound(_rows(id(r)));
deba@459
  1562
    }
deba@459
  1563
deba@459
  1564
    /// Set the upper bound of a row (i.e a constraint)
deba@459
  1565
deba@459
  1566
    /// The upper bound of a constraint (row) has to be given by an
deba@459
  1567
    /// extended number of type Value, i.e. a finite number of type
deba@459
  1568
    /// Value or -\ref INF.
deba@459
  1569
    void rowUpperBound(Row r, Value value) {
ggab90@1130
  1570
      _setRowUpperBound(_rows(id(r)),value);
deba@459
  1571
    }
deba@459
  1572
deba@459
  1573
    /// Get the upper bound of a row (i.e a constraint)
deba@459
  1574
deba@459
  1575
    /// This function returns the upper bound for row (constraint) \c c
deba@459
  1576
    /// (this might be -\ref INF as well).
deba@459
  1577
    ///\return The upper bound for row \c r
deba@459
  1578
    Value rowUpperBound(Row r) const {
ggab90@1130
  1579
      return _getRowUpperBound(_rows(id(r)));
deba@458
  1580
    }
deba@458
  1581
deba@458
  1582
    ///Set an element of the objective function
ggab90@1130
  1583
    void objCoeff(Col c, Value v) {_setObjCoeff(_cols(id(c)),v); };
deba@458
  1584
deba@458
  1585
    ///Get an element of the objective function
ggab90@1130
  1586
    Value objCoeff(Col c) const { return _getObjCoeff(_cols(id(c))); };
deba@458
  1587
deba@458
  1588
    ///Set the objective function
deba@458
  1589
deba@458
  1590
    ///\param e is a linear expression of type \ref Expr.
deba@459
  1591
    ///
deba@459
  1592
    void obj(const Expr& e) {
ggab90@1130
  1593
      _setObjCoeffs(ExprIterator(e.comps.begin(), _cols),
ggab90@1130
  1594
                    ExprIterator(e.comps.end(), _cols));
deba@459
  1595
      obj_const_comp = *e;
deba@458
  1596
    }
deba@458
  1597
deba@458
  1598
    ///Get the objective function
deba@458
  1599
deba@459
  1600
    ///\return the objective function as a linear expression of type
deba@459
  1601
    ///Expr.
deba@458
  1602
    Expr obj() const {
deba@458
  1603
      Expr e;
ggab90@1130
  1604
      _getObjCoeffs(InsertIterator(e.comps, _cols));
deba@459
  1605
      *e = obj_const_comp;
deba@458
  1606
      return e;
deba@458
  1607
    }
deba@458
  1608
deba@458
  1609
deba@459
  1610
    ///Set the direction of optimization
deba@459
  1611
    void sense(Sense sense) { _setSense(sense); }
deba@458
  1612
deba@459
  1613
    ///Query the direction of the optimization
deba@459
  1614
    Sense sense() const {return _getSense(); }
deba@458
  1615
deba@459
  1616
    ///Set the sense to maximization
deba@459
  1617
    void max() { _setSense(MAX); }
deba@459
  1618
deba@459
  1619
    ///Set the sense to maximization
deba@459
  1620
    void min() { _setSense(MIN); }
deba@459
  1621
alpar@1063
  1622
    ///Clear the problem
ggab90@1130
  1623
    void clear() { _clear(); _rows.clear(); _cols.clear(); }
deba@458
  1624
alpar@1063
  1625
    /// Set the message level of the solver
deba@576
  1626
    void messageLevel(MessageLevel level) { _messageLevel(level); }
deba@576
  1627
alpar@1063
  1628
    /// Write the problem to a file in the given format
alpar@1063
  1629
alpar@1063
  1630
    /// This function writes the problem to a file in the given format.
alpar@1063
  1631
    /// Different solver backends may support different formats.
alpar@1063
  1632
    /// Trying to write in an unsupported format will trigger
alpar@1063
  1633
    /// \ref UnsupportedFormatError. For the supported formats,
alpar@1063
  1634
    /// visit the documentation of the base class of the related backends
alpar@1063
  1635
    /// (\ref CplexBase, \ref GlpkBase etc.)
alpar@1063
  1636
    /// \param file The file path
alpar@1063
  1637
    /// \param format The output file format.
alpar@1063
  1638
    void write(std::string file, std::string format = "MPS") const
alpar@1063
  1639
    {
alpar@1063
  1640
      _write(file.c_str(),format.c_str());
alpar@1063
  1641
    }
alpar@1063
  1642
deba@458
  1643
    ///@}
deba@458
  1644
deba@459
  1645
  };
deba@459
  1646
deba@459
  1647
  /// Addition
deba@459
  1648
deba@459
  1649
  ///\relates LpBase::Expr
deba@459
  1650
  ///
deba@459
  1651
  inline LpBase::Expr operator+(const LpBase::Expr &a, const LpBase::Expr &b) {
deba@459
  1652
    LpBase::Expr tmp(a);
deba@459
  1653
    tmp+=b;
deba@459
  1654
    return tmp;
deba@459
  1655
  }
deba@459
  1656
  ///Substraction
deba@459
  1657
deba@459
  1658
  ///\relates LpBase::Expr
deba@459
  1659
  ///
deba@459
  1660
  inline LpBase::Expr operator-(const LpBase::Expr &a, const LpBase::Expr &b) {
deba@459
  1661
    LpBase::Expr tmp(a);
deba@459
  1662
    tmp-=b;
deba@459
  1663
    return tmp;
deba@459
  1664
  }
deba@459
  1665
  ///Multiply with constant
deba@459
  1666
deba@459
  1667
  ///\relates LpBase::Expr
deba@459
  1668
  ///
deba@459
  1669
  inline LpBase::Expr operator*(const LpBase::Expr &a, const LpBase::Value &b) {
deba@459
  1670
    LpBase::Expr tmp(a);
deba@459
  1671
    tmp*=b;
deba@459
  1672
    return tmp;
deba@459
  1673
  }
deba@459
  1674
deba@459
  1675
  ///Multiply with constant
deba@459
  1676
deba@459
  1677
  ///\relates LpBase::Expr
deba@459
  1678
  ///
deba@459
  1679
  inline LpBase::Expr operator*(const LpBase::Value &a, const LpBase::Expr &b) {
deba@459
  1680
    LpBase::Expr tmp(b);
deba@459
  1681
    tmp*=a;
deba@459
  1682
    return tmp;
deba@459
  1683
  }
deba@459
  1684
  ///Divide with constant
deba@459
  1685
deba@459
  1686
  ///\relates LpBase::Expr
deba@459
  1687
  ///
deba@459
  1688
  inline LpBase::Expr operator/(const LpBase::Expr &a, const LpBase::Value &b) {
deba@459
  1689
    LpBase::Expr tmp(a);
deba@459
  1690
    tmp/=b;
deba@459
  1691
    return tmp;
deba@459
  1692
  }
deba@459
  1693
deba@459
  1694
  ///Create constraint
deba@459
  1695
deba@459
  1696
  ///\relates LpBase::Constr
deba@459
  1697
  ///
deba@459
  1698
  inline LpBase::Constr operator<=(const LpBase::Expr &e,
deba@459
  1699
                                   const LpBase::Expr &f) {
retvari@957
  1700
    return LpBase::Constr(0, f - e, LpBase::NaN);
deba@459
  1701
  }
deba@459
  1702
deba@459
  1703
  ///Create constraint
deba@459
  1704
deba@459
  1705
  ///\relates LpBase::Constr
deba@459
  1706
  ///
deba@459
  1707
  inline LpBase::Constr operator<=(const LpBase::Value &e,
deba@459
  1708
                                   const LpBase::Expr &f) {
deba@459
  1709
    return LpBase::Constr(e, f, LpBase::NaN);
deba@459
  1710
  }
deba@459
  1711
deba@459
  1712
  ///Create constraint
deba@459
  1713
deba@459
  1714
  ///\relates LpBase::Constr
deba@459
  1715
  ///
deba@459
  1716
  inline LpBase::Constr operator<=(const LpBase::Expr &e,
deba@459
  1717
                                   const LpBase::Value &f) {
retvari@957
  1718
    return LpBase::Constr(LpBase::NaN, e, f);
deba@459
  1719
  }
deba@459
  1720
deba@459
  1721
  ///Create constraint
deba@459
  1722
deba@459
  1723
  ///\relates LpBase::Constr
deba@459
  1724
  ///
deba@459
  1725
  inline LpBase::Constr operator>=(const LpBase::Expr &e,
deba@459
  1726
                                   const LpBase::Expr &f) {
retvari@957
  1727
    return LpBase::Constr(0, e - f, LpBase::NaN);
deba@459
  1728
  }
deba@459
  1729
deba@459
  1730
deba@459
  1731
  ///Create constraint
deba@459
  1732
deba@459
  1733
  ///\relates LpBase::Constr
deba@459
  1734
  ///
deba@459
  1735
  inline LpBase::Constr operator>=(const LpBase::Value &e,
deba@459
  1736
                                   const LpBase::Expr &f) {
deba@459
  1737
    return LpBase::Constr(LpBase::NaN, f, e);
deba@459
  1738
  }
deba@459
  1739
deba@459
  1740
deba@459
  1741
  ///Create constraint
deba@459
  1742
deba@459
  1743
  ///\relates LpBase::Constr
deba@459
  1744
  ///
deba@459
  1745
  inline LpBase::Constr operator>=(const LpBase::Expr &e,
deba@459
  1746
                                   const LpBase::Value &f) {
retvari@957
  1747
    return LpBase::Constr(f, e, LpBase::NaN);
deba@459
  1748
  }
deba@459
  1749
deba@459
  1750
  ///Create constraint
deba@459
  1751
deba@459
  1752
  ///\relates LpBase::Constr
deba@459
  1753
  ///
deba@459
  1754
  inline LpBase::Constr operator==(const LpBase::Expr &e,
deba@459
  1755
                                   const LpBase::Value &f) {
deba@459
  1756
    return LpBase::Constr(f, e, f);
deba@459
  1757
  }
deba@459
  1758
deba@459
  1759
  ///Create constraint
deba@459
  1760
deba@459
  1761
  ///\relates LpBase::Constr
deba@459
  1762
  ///
deba@459
  1763
  inline LpBase::Constr operator==(const LpBase::Expr &e,
deba@459
  1764
                                   const LpBase::Expr &f) {
deba@459
  1765
    return LpBase::Constr(0, f - e, 0);
deba@459
  1766
  }
deba@459
  1767
deba@459
  1768
  ///Create constraint
deba@459
  1769
deba@459
  1770
  ///\relates LpBase::Constr
deba@459
  1771
  ///
deba@459
  1772
  inline LpBase::Constr operator<=(const LpBase::Value &n,
deba@459
  1773
                                   const LpBase::Constr &c) {
deba@459
  1774
    LpBase::Constr tmp(c);
alpar@487
  1775
    LEMON_ASSERT(isNaN(tmp.lowerBound()), "Wrong LP constraint");
deba@459
  1776
    tmp.lowerBound()=n;
deba@459
  1777
    return tmp;
deba@459
  1778
  }
deba@459
  1779
  ///Create constraint
deba@459
  1780
deba@459
  1781
  ///\relates LpBase::Constr
deba@459
  1782
  ///
deba@459
  1783
  inline LpBase::Constr operator<=(const LpBase::Constr &c,
deba@459
  1784
                                   const LpBase::Value &n)
deba@459
  1785
  {
deba@459
  1786
    LpBase::Constr tmp(c);
alpar@487
  1787
    LEMON_ASSERT(isNaN(tmp.upperBound()), "Wrong LP constraint");
deba@459
  1788
    tmp.upperBound()=n;
deba@459
  1789
    return tmp;
deba@459
  1790
  }
deba@459
  1791
deba@459
  1792
  ///Create constraint
deba@459
  1793
deba@459
  1794
  ///\relates LpBase::Constr
deba@459
  1795
  ///
deba@459
  1796
  inline LpBase::Constr operator>=(const LpBase::Value &n,
deba@459
  1797
                                   const LpBase::Constr &c) {
deba@459
  1798
    LpBase::Constr tmp(c);
alpar@487
  1799
    LEMON_ASSERT(isNaN(tmp.upperBound()), "Wrong LP constraint");
deba@459
  1800
    tmp.upperBound()=n;
deba@459
  1801
    return tmp;
deba@459
  1802
  }
deba@459
  1803
  ///Create constraint
deba@459
  1804
deba@459
  1805
  ///\relates LpBase::Constr
deba@459
  1806
  ///
deba@459
  1807
  inline LpBase::Constr operator>=(const LpBase::Constr &c,
deba@459
  1808
                                   const LpBase::Value &n)
deba@459
  1809
  {
deba@459
  1810
    LpBase::Constr tmp(c);
alpar@487
  1811
    LEMON_ASSERT(isNaN(tmp.lowerBound()), "Wrong LP constraint");
deba@459
  1812
    tmp.lowerBound()=n;
deba@459
  1813
    return tmp;
deba@459
  1814
  }
deba@459
  1815
deba@459
  1816
  ///Addition
deba@459
  1817
deba@459
  1818
  ///\relates LpBase::DualExpr
deba@459
  1819
  ///
deba@459
  1820
  inline LpBase::DualExpr operator+(const LpBase::DualExpr &a,
deba@459
  1821
                                    const LpBase::DualExpr &b) {
deba@459
  1822
    LpBase::DualExpr tmp(a);
deba@459
  1823
    tmp+=b;
deba@459
  1824
    return tmp;
deba@459
  1825
  }
deba@459
  1826
  ///Substraction
deba@459
  1827
deba@459
  1828
  ///\relates LpBase::DualExpr
deba@459
  1829
  ///
deba@459
  1830
  inline LpBase::DualExpr operator-(const LpBase::DualExpr &a,
deba@459
  1831
                                    const LpBase::DualExpr &b) {
deba@459
  1832
    LpBase::DualExpr tmp(a);
deba@459
  1833
    tmp-=b;
deba@459
  1834
    return tmp;
deba@459
  1835
  }
deba@459
  1836
  ///Multiply with constant
deba@459
  1837
deba@459
  1838
  ///\relates LpBase::DualExpr
deba@459
  1839
  ///
deba@459
  1840
  inline LpBase::DualExpr operator*(const LpBase::DualExpr &a,
deba@459
  1841
                                    const LpBase::Value &b) {
deba@459
  1842
    LpBase::DualExpr tmp(a);
deba@459
  1843
    tmp*=b;
deba@459
  1844
    return tmp;
deba@459
  1845
  }
deba@459
  1846
deba@459
  1847
  ///Multiply with constant
deba@459
  1848
deba@459
  1849
  ///\relates LpBase::DualExpr
deba@459
  1850
  ///
deba@459
  1851
  inline LpBase::DualExpr operator*(const LpBase::Value &a,
deba@459
  1852
                                    const LpBase::DualExpr &b) {
deba@459
  1853
    LpBase::DualExpr tmp(b);
deba@459
  1854
    tmp*=a;
deba@459
  1855
    return tmp;
deba@459
  1856
  }
deba@459
  1857
  ///Divide with constant
deba@459
  1858
deba@459
  1859
  ///\relates LpBase::DualExpr
deba@459
  1860
  ///
deba@459
  1861
  inline LpBase::DualExpr operator/(const LpBase::DualExpr &a,
deba@459
  1862
                                    const LpBase::Value &b) {
deba@459
  1863
    LpBase::DualExpr tmp(a);
deba@459
  1864
    tmp/=b;
deba@459
  1865
    return tmp;
deba@459
  1866
  }
deba@459
  1867
deba@459
  1868
  /// \ingroup lp_group
deba@459
  1869
  ///
deba@459
  1870
  /// \brief Common base class for LP solvers
deba@459
  1871
  ///
deba@459
  1872
  /// This class is an abstract base class for LP solvers. This class
deba@459
  1873
  /// provides a full interface for set and modify an LP problem,
deba@459
  1874
  /// solve it and retrieve the solution. You can use one of the
deba@459
  1875
  /// descendants as a concrete implementation, or the \c Lp
deba@459
  1876
  /// default LP solver. However, if you would like to handle LP
deba@459
  1877
  /// solvers as reference or pointer in a generic way, you can use
deba@459
  1878
  /// this class directly.
deba@459
  1879
  class LpSolver : virtual public LpBase {
deba@459
  1880
  public:
deba@459
  1881
deba@459
  1882
    /// The problem types for primal and dual problems
deba@459
  1883
    enum ProblemType {
kpeter@584
  1884
      /// = 0. Feasible solution hasn't been found (but may exist).
deba@459
  1885
      UNDEFINED = 0,
kpeter@584
  1886
      /// = 1. The problem has no feasible solution.
deba@459
  1887
      INFEASIBLE = 1,
kpeter@584
  1888
      /// = 2. Feasible solution found.
deba@459
  1889
      FEASIBLE = 2,
kpeter@584
  1890
      /// = 3. Optimal solution exists and found.
deba@459
  1891
      OPTIMAL = 3,
kpeter@584
  1892
      /// = 4. The cost function is unbounded.
deba@459
  1893
      UNBOUNDED = 4
deba@459
  1894
    };
deba@459
  1895
deba@459
  1896
    ///The basis status of variables
deba@459
  1897
    enum VarStatus {
deba@459
  1898
      /// The variable is in the basis
alpar@877
  1899
      BASIC,
deba@459
  1900
      /// The variable is free, but not basic
deba@459
  1901
      FREE,
alpar@877
  1902
      /// The variable has active lower bound
deba@459
  1903
      LOWER,
deba@459
  1904
      /// The variable has active upper bound
deba@459
  1905
      UPPER,
deba@459
  1906
      /// The variable is non-basic and fixed
deba@459
  1907
      FIXED
deba@459
  1908
    };
deba@459
  1909
deba@459
  1910
  protected:
deba@459
  1911
deba@459
  1912
    virtual SolveExitStatus _solve() = 0;
deba@459
  1913
deba@459
  1914
    virtual Value _getPrimal(int i) const = 0;
deba@459
  1915
    virtual Value _getDual(int i) const = 0;
deba@459
  1916
deba@459
  1917
    virtual Value _getPrimalRay(int i) const = 0;
deba@459
  1918
    virtual Value _getDualRay(int i) const = 0;
deba@459
  1919
deba@459
  1920
    virtual Value _getPrimalValue() const = 0;
deba@459
  1921
deba@459
  1922
    virtual VarStatus _getColStatus(int i) const = 0;
deba@459
  1923
    virtual VarStatus _getRowStatus(int i) const = 0;
deba@459
  1924
deba@459
  1925
    virtual ProblemType _getPrimalType() const = 0;
deba@459
  1926
    virtual ProblemType _getDualType() const = 0;
deba@459
  1927
deba@459
  1928
  public:
deba@458
  1929
alpar@540
  1930
    ///Allocate a new LP problem instance
alpar@540
  1931
    virtual LpSolver* newSolver() const = 0;
alpar@540
  1932
    ///Make a copy of the LP problem
alpar@540
  1933
    virtual LpSolver* cloneSolver() const = 0;
alpar@540
  1934
deba@458
  1935
    ///\name Solve the LP
deba@458
  1936
deba@458
  1937
    ///@{
deba@458
  1938
deba@458
  1939
    ///\e Solve the LP problem at hand
deba@458
  1940
    ///
deba@458
  1941
    ///\return The result of the optimization procedure. Possible
deba@458
  1942
    ///values and their meanings can be found in the documentation of
deba@458
  1943
    ///\ref SolveExitStatus.
deba@458
  1944
    SolveExitStatus solve() { return _solve(); }
deba@458
  1945
deba@458
  1946
    ///@}
deba@458
  1947
kpeter@584
  1948
    ///\name Obtain the Solution
deba@458
  1949
deba@458
  1950
    ///@{
deba@458
  1951
deba@459
  1952
    /// The type of the primal problem
deba@459
  1953
    ProblemType primalType() const {
deba@459
  1954
      return _getPrimalType();
deba@458
  1955
    }
deba@458
  1956
deba@459
  1957
    /// The type of the dual problem
deba@459
  1958
    ProblemType dualType() const {
deba@459
  1959
      return _getDualType();
deba@458
  1960
    }
deba@458
  1961
deba@459
  1962
    /// Return the primal value of the column
deba@459
  1963
deba@459
  1964
    /// Return the primal value of the column.
deba@459
  1965
    /// \pre The problem is solved.
ggab90@1130
  1966
    Value primal(Col c) const { return _getPrimal(_cols(id(c))); }
deba@459
  1967
deba@459
  1968
    /// Return the primal value of the expression
deba@459
  1969
deba@459
  1970
    /// Return the primal value of the expression, i.e. the dot
deba@459
  1971
    /// product of the primal solution and the expression.
deba@459
  1972
    /// \pre The problem is solved.
deba@459
  1973
    Value primal(const Expr& e) const {
deba@459
  1974
      double res = *e;
deba@459
  1975
      for (Expr::ConstCoeffIt c(e); c != INVALID; ++c) {
deba@459
  1976
        res += *c * primal(c);
deba@459
  1977
      }
deba@459
  1978
      return res;
deba@458
  1979
    }
deba@459
  1980
    /// Returns a component of the primal ray
alpar@877
  1981
deba@459
  1982
    /// The primal ray is solution of the modified primal problem,
deba@459
  1983
    /// where we change each finite bound to 0, and we looking for a
deba@459
  1984
    /// negative objective value in case of minimization, and positive
deba@459
  1985
    /// objective value for maximization. If there is such solution,
deba@459
  1986
    /// that proofs the unsolvability of the dual problem, and if a
deba@459
  1987
    /// feasible primal solution exists, then the unboundness of
deba@459
  1988
    /// primal problem.
deba@459
  1989
    ///
deba@459
  1990
    /// \pre The problem is solved and the dual problem is infeasible.
deba@459
  1991
    /// \note Some solvers does not provide primal ray calculation
deba@459
  1992
    /// functions.
ggab90@1130
  1993
    Value primalRay(Col c) const { return _getPrimalRay(_cols(id(c))); }
deba@458
  1994
deba@459
  1995
    /// Return the dual value of the row
deba@459
  1996
deba@459
  1997
    /// Return the dual value of the row.
deba@459
  1998
    /// \pre The problem is solved.
ggab90@1130
  1999
    Value dual(Row r) const { return _getDual(_rows(id(r))); }
deba@459
  2000
deba@459
  2001
    /// Return the dual value of the dual expression
deba@459
  2002
deba@459
  2003
    /// Return the dual value of the dual expression, i.e. the dot
deba@459
  2004
    /// product of the dual solution and the dual expression.
deba@459
  2005
    /// \pre The problem is solved.
deba@459
  2006
    Value dual(const DualExpr& e) const {
deba@459
  2007
      double res = 0.0;
deba@459
  2008
      for (DualExpr::ConstCoeffIt r(e); r != INVALID; ++r) {
deba@459
  2009
        res += *r * dual(r);
deba@458
  2010
      }
deba@458
  2011
      return res;
deba@458
  2012
    }
deba@458
  2013
deba@459
  2014
    /// Returns a component of the dual ray
alpar@877
  2015
deba@459
  2016
    /// The dual ray is solution of the modified primal problem, where
deba@459
  2017
    /// we change each finite bound to 0 (i.e. the objective function
deba@459
  2018
    /// coefficients in the primal problem), and we looking for a
deba@459
  2019
    /// ositive objective value. If there is such solution, that
deba@459
  2020
    /// proofs the unsolvability of the primal problem, and if a
deba@459
  2021
    /// feasible dual solution exists, then the unboundness of
deba@459
  2022
    /// dual problem.
deba@459
  2023
    ///
deba@459
  2024
    /// \pre The problem is solved and the primal problem is infeasible.
deba@459
  2025
    /// \note Some solvers does not provide dual ray calculation
deba@459
  2026
    /// functions.
ggab90@1130
  2027
    Value dualRay(Row r) const { return _getDualRay(_rows(id(r))); }
deba@458
  2028
deba@459
  2029
    /// Return the basis status of the column
deba@458
  2030
deba@459
  2031
    /// \see VarStatus
ggab90@1130
  2032
    VarStatus colStatus(Col c) const { return _getColStatus(_cols(id(c))); }
deba@459
  2033
deba@459
  2034
    /// Return the basis status of the row
deba@459
  2035
deba@459
  2036
    /// \see VarStatus
ggab90@1130
  2037
    VarStatus rowStatus(Row r) const { return _getRowStatus(_rows(id(r))); }
deba@459
  2038
deba@459
  2039
    ///The value of the objective function
deba@458
  2040
deba@458
  2041
    ///\return
deba@458
  2042
    ///- \ref INF or -\ref INF means either infeasibility or unboundedness
deba@458
  2043
    /// of the primal problem, depending on whether we minimize or maximize.
deba@458
  2044
    ///- \ref NaN if no primal solution is found.
deba@458
  2045
    ///- The (finite) objective value if an optimal solution is found.
deba@459
  2046
    Value primal() const { return _getPrimalValue()+obj_const_comp;}
deba@458
  2047
    ///@}
deba@458
  2048
deba@459
  2049
  protected:
deba@459
  2050
deba@458
  2051
  };
deba@458
  2052
deba@458
  2053
deba@458
  2054
  /// \ingroup lp_group
deba@458
  2055
  ///
deba@458
  2056
  /// \brief Common base class for MIP solvers
deba@459
  2057
  ///
deba@459
  2058
  /// This class is an abstract base class for MIP solvers. This class
deba@459
  2059
  /// provides a full interface for set and modify an MIP problem,
deba@459
  2060
  /// solve it and retrieve the solution. You can use one of the
deba@459
  2061
  /// descendants as a concrete implementation, or the \c Lp
deba@459
  2062
  /// default MIP solver. However, if you would like to handle MIP
deba@459
  2063
  /// solvers as reference or pointer in a generic way, you can use
deba@459
  2064
  /// this class directly.
deba@459
  2065
  class MipSolver : virtual public LpBase {
deba@458
  2066
  public:
deba@458
  2067
deba@459
  2068
    /// The problem types for MIP problems
deba@459
  2069
    enum ProblemType {
kpeter@584
  2070
      /// = 0. Feasible solution hasn't been found (but may exist).
deba@459
  2071
      UNDEFINED = 0,
kpeter@584
  2072
      /// = 1. The problem has no feasible solution.
deba@459
  2073
      INFEASIBLE = 1,
kpeter@584
  2074
      /// = 2. Feasible solution found.
deba@459
  2075
      FEASIBLE = 2,
kpeter@584
  2076
      /// = 3. Optimal solution exists and found.
deba@459
  2077
      OPTIMAL = 3,
kpeter@584
  2078
      /// = 4. The cost function is unbounded.
kpeter@584
  2079
      ///The Mip or at least the relaxed problem is unbounded.
deba@459
  2080
      UNBOUNDED = 4
deba@459
  2081
    };
deba@459
  2082
alpar@540
  2083
    ///Allocate a new MIP problem instance
alpar@540
  2084
    virtual MipSolver* newSolver() const = 0;
alpar@540
  2085
    ///Make a copy of the MIP problem
alpar@540
  2086
    virtual MipSolver* cloneSolver() const = 0;
alpar@540
  2087
deba@459
  2088
    ///\name Solve the MIP
deba@459
  2089
deba@459
  2090
    ///@{
deba@459
  2091
deba@459
  2092
    /// Solve the MIP problem at hand
deba@459
  2093
    ///
deba@459
  2094
    ///\return The result of the optimization procedure. Possible
deba@459
  2095
    ///values and their meanings can be found in the documentation of
deba@459
  2096
    ///\ref SolveExitStatus.
deba@459
  2097
    SolveExitStatus solve() { return _solve(); }
deba@459
  2098
deba@459
  2099
    ///@}
deba@459
  2100
kpeter@584
  2101
    ///\name Set Column Type
deba@459
  2102
    ///@{
deba@459
  2103
deba@459
  2104
    ///Possible variable (column) types (e.g. real, integer, binary etc.)
deba@458
  2105
    enum ColTypes {
kpeter@584
  2106
      /// = 0. Continuous variable (default).
deba@458
  2107
      REAL = 0,
kpeter@584
  2108
      /// = 1. Integer variable.
deba@459
  2109
      INTEGER = 1
deba@458
  2110
    };
deba@458
  2111
deba@459
  2112
    ///Sets the type of the given column to the given type
deba@459
  2113
deba@459
  2114
    ///Sets the type of the given column to the given type.
deba@458
  2115
    ///
deba@458
  2116
    void colType(Col c, ColTypes col_type) {
ggab90@1130
  2117
      _setColType(_cols(id(c)),col_type);
deba@458
  2118
    }
deba@458
  2119
deba@458
  2120
    ///Gives back the type of the column.
deba@459
  2121
deba@459
  2122
    ///Gives back the type of the column.
deba@458
  2123
    ///
deba@458
  2124
    ColTypes colType(Col c) const {
ggab90@1130
  2125
      return _getColType(_cols(id(c)));
deba@459
  2126
    }
deba@459
  2127
    ///@}
deba@459
  2128
kpeter@584
  2129
    ///\name Obtain the Solution
deba@459
  2130
deba@459
  2131
    ///@{
deba@459
  2132
deba@459
  2133
    /// The type of the MIP problem
deba@459
  2134
    ProblemType type() const {
deba@459
  2135
      return _getType();
deba@458
  2136
    }
deba@458
  2137
deba@459
  2138
    /// Return the value of the row in the solution
deba@459
  2139
deba@459
  2140
    ///  Return the value of the row in the solution.
deba@459
  2141
    /// \pre The problem is solved.
ggab90@1130
  2142
    Value sol(Col c) const { return _getSol(_cols(id(c))); }
deba@459
  2143
deba@459
  2144
    /// Return the value of the expression in the solution
deba@459
  2145
deba@459
  2146
    /// Return the value of the expression in the solution, i.e. the
deba@459
  2147
    /// dot product of the solution and the expression.
deba@459
  2148
    /// \pre The problem is solved.
deba@459
  2149
    Value sol(const Expr& e) const {
deba@459
  2150
      double res = *e;
deba@459
  2151
      for (Expr::ConstCoeffIt c(e); c != INVALID; ++c) {
deba@459
  2152
        res += *c * sol(c);
deba@459
  2153
      }
deba@459
  2154
      return res;
deba@458
  2155
    }
deba@459
  2156
    ///The value of the objective function
alpar@877
  2157
deba@459
  2158
    ///\return
deba@459
  2159
    ///- \ref INF or -\ref INF means either infeasibility or unboundedness
deba@459
  2160
    /// of the problem, depending on whether we minimize or maximize.
deba@459
  2161
    ///- \ref NaN if no primal solution is found.
deba@459
  2162
    ///- The (finite) objective value if an optimal solution is found.
deba@459
  2163
    Value solValue() const { return _getSolValue()+obj_const_comp;}
deba@459
  2164
    ///@}
deba@458
  2165
deba@458
  2166
  protected:
deba@458
  2167
deba@459
  2168
    virtual SolveExitStatus _solve() = 0;
deba@459
  2169
    virtual ColTypes _getColType(int col) const = 0;
deba@459
  2170
    virtual void _setColType(int col, ColTypes col_type) = 0;
deba@459
  2171
    virtual ProblemType _getType() const = 0;
deba@459
  2172
    virtual Value _getSol(int i) const = 0;
deba@459
  2173
    virtual Value _getSolValue() const = 0;
deba@458
  2174
deba@458
  2175
  };
deba@458
  2176
deba@458
  2177
deba@458
  2178
deba@458
  2179
} //namespace lemon
deba@458
  2180
deba@458
  2181
#endif //LEMON_LP_BASE_H