source: stable/0.2/Couenne/src/problem/CouenneProblem.hpp @ 159

Last change on this file since 159 was 159, checked in by pbelotti, 11 years ago

created new stable branch 0.2 from trunk (rev. 157)

File size: 18.9 KB
Line 
1/* $Id: CouenneProblem.hpp 155 2009-06-16 20:19:39Z pbelotti $
2 *
3 * Name:    CouenneProblem.hpp
4 * Author:  Pietro Belotti
5 * Purpose: define the class CouenneProblem
6 *
7 * (C) Carnegie-Mellon University, 2006-09.
8 * This file is licensed under the Common Public License (CPL)
9 */
10
11#ifndef COUENNE_PROBLEM_HPP
12#define COUENNE_PROBLEM_HPP
13
14#include <vector>
15
16#include "OsiRowCut.hpp"
17#include "BonBabInfos.hpp"
18
19#include "CouenneTypes.hpp"
20#include "expression.hpp"
21#include "exprAux.hpp"
22#include "CouenneProblemElem.hpp"
23#include "CouenneObject.hpp"
24#include "CouenneJournalist.hpp"
25#include "domain.hpp"
26
27struct ASL;
28struct expr;
29
30class DepGraph;
31class CouenneCutGenerator;
32class CouenneSolverInterface;
33class quadElem;
34class LinMap;
35class QuadMap;
36
37// default tolerance for checking feasibility (and integrality) of NLP solutions
38const CouNumber feas_tolerance_default = 1e-5;
39
40
41/** Class for MINLP problems with symbolic information
42 *
43 *  It is read from an AMPL .nl file and contains variables, AMPL's
44 *  "defined variables" (aka common expressions), objective(s), and
45 *  constraints in the form of expression's. Changes throughout the
46 *  program occur in standardization.
47 */
48
49class CouenneProblem {
50
51  /// Class for storing a global cutoff for a CouenneProblem and all
52  /// its clones
53  class GlobalCutOff {
54  private:
55    GlobalCutOff(const GlobalCutOff&);
56    double cutoff_;
57  public:
58    GlobalCutOff ()         : cutoff_ (COIN_DBL_MAX) {}
59    GlobalCutOff (double c) : cutoff_ (c) {}
60    ~GlobalCutOff() {}
61    inline void setCutOff (double cutoff) {cutoff_ = cutoff;}
62    inline double getCutOff() const {return cutoff_;}
63  };
64
65  /// structure to record fixed, non-fixed, and continuous variables
66  enum fixType {UNFIXED, FIXED, CONTINUOUS};
67
68 protected:
69
70  /// problem name
71  std::string problemName_;
72
73  std::vector <exprVar           *> variables_;   ///< Variables (original, auxiliary, and defined)
74  std::vector <CouenneObjective  *> objectives_;  ///< Objectives
75  std::vector <CouenneConstraint *> constraints_; ///< Constraints
76
77  /// AMPL's common expressions (read from AMPL through structures cexps and cexps1)
78  std::vector <expression *> commonexprs_; 
79
80  mutable Domain domain_; ///< current point and bounds;
81
82  /// Expression map for comparison in standardization and to count
83  /// occurrences of an auxiliary
84  std::set <exprAux *, compExpr> *auxSet_;
85
86  /// Number of elements in the x_, lb_, ub_ arrays
87  mutable int curnvars_;
88
89  /// Number of discrete variables
90  int nIntVars_;
91
92  /// Best solution known to be loaded from file -- for testing purposes
93  mutable CouNumber *optimum_;
94
95  /// Best known objective function
96  CouNumber bestObj_;
97
98  /// Indices of variables appearing in products (used for SDP cuts)
99  int *quadIndex_;
100
101  /// Variables that have commuted to auxiliary
102  bool *commuted_;
103
104  /// numbering of variables. No variable xi with associated pi(i)
105  /// greater than pi(j) should be evaluated before variable xj
106  int *numbering_;
107
108  /// Number of "defined variables" (aka "common expressions")
109  int ndefined_;
110
111  /// Dependence (acyclic) graph: shows dependence of all auxiliary
112  /// variables on one another and on original variables. Used to
113  /// create a numbering of all variables for evaluation and bound
114  /// tightening (reverse order for implied bounds)
115  DepGraph *graph_;
116
117  /// Number of original variables
118  int nOrigVars_;
119
120  /// Number of original constraints (disregarding those that turned
121  /// into auxiliary variable definition)
122  int nOrigCons_;
123
124  /// Number of original integer variables
125  int nOrigIntVars_;
126
127  /// Pointer to a global cutoff object
128  mutable GlobalCutOff* pcutoff_;
129
130  /// flag indicating if this class is creator of global cutoff object
131  mutable bool created_pcutoff_;
132
133  bool doFBBT_;  ///< do Feasibility-based bound tightening
134  bool doRCBT_;  ///< do reduced cost      bound tightening
135  bool doOBBT_;  ///< do Optimality-based  bound tightening
136  bool doABT_;   ///< do Aggressive        bound tightening
137
138  int logObbtLev_;   ///< frequency of Optimality-based bound tightening
139  int logAbtLev_;    ///< frequency of Aggressive       bound tightening
140
141  /// SmartPointer to the Journalist
142  JnlstPtr jnlst_;
143
144  /// window around known optimum (for testing purposes)
145  CouNumber opt_window_;
146
147  /// Use quadratic expressions?
148  bool useQuadratic_;
149
150  /// feasibility tolerance (to be used in checkNLP)
151  CouNumber feas_tolerance_;
152
153  /// inverse dependence structure: for each variable x give set of
154  /// auxiliary variables (or better, their indices) whose expression
155  /// depends on x
156  std::vector <std::set <int> > dependence_;
157
158  /// vector of pointer to CouenneObjects. Used by CouenneVarObjects
159  /// when finding all objects related to (having as argument) a
160  /// single variable
161  std::vector <CouenneObject> objects_;
162
163  /// each element is true if variable is integer and, if auxiliary,
164  /// depends on no integer
165  mutable int *integerRank_;
166
167  /// numberInRank_ [i] is the number of integer variables in rank i
168  mutable std::vector <int> numberInRank_;
169
170  /// maximum cpu time
171  double maxCpuTime_;
172
173  /// options
174  Bonmin::BabSetupBase *bonBase_;
175
176#ifdef COIN_HAS_ASL
177  /// AMPL structure pointer (temporary --- looking forward to embedding into OS...)
178  ASL *asl_;
179#endif
180
181  /// some originals may be unused due to their zero multiplicity
182  /// (that happens when they are duplicates). This array keeps track
183  /// of their indices and is sorted by evaluation order
184  int *unusedOriginalsIndices_;
185
186  /// number of unused originals
187  int nUnusedOriginals_;
188
189 public:
190
191  CouenneProblem  (ASL * = NULL,
192                   Bonmin::BabSetupBase *base = NULL,
193                   JnlstPtr jnlst = NULL);  ///< Constructor
194  CouenneProblem  (const CouenneProblem &); ///< Copy constructor
195  ~CouenneProblem ();                       ///< Destructor
196
197  /// Clone method (for use within CouenneCutGenerator::clone)
198  CouenneProblem *clone () const
199  {return new CouenneProblem (*this);}
200
201  int nObjs     () const {return (int) objectives_.   size ();} ///< Get number of objectives
202  int nCons     () const {return (int) constraints_.  size ();} ///< Get number of constraints
203  int nOrigCons () const {return nOrigCons_;}                   ///< Get number of original constraints
204
205  inline int nOrigVars    () const {return nOrigVars_;}                ///< Number of orig. variables
206  inline int nDefVars     () const {return ndefined_;}                 ///< Number of def'd variables
207  inline int nOrigIntVars () const {return nOrigIntVars_;}             ///< Number of original integers
208  inline int nIntVars     () const {return nIntVars_;}                 ///< Number of integer variables
209  inline int nVars        () const {return (int) variables_. size ();} ///< Total number of variables
210
211  /// get evaluation order index
212  inline int evalOrder (int i) const
213  {return numbering_ [i];}
214
215  /// get evaluation order vector (numbering_)
216  inline int *evalVector ()
217  {return numbering_;}
218
219  // get elements from vectors
220  inline CouenneConstraint *Con (int i) const {return constraints_ [i];} ///< i-th constraint
221  inline CouenneObjective  *Obj (int i) const {return objectives_  [i];} ///< i-th objective
222
223  /// Return pointer to i-th variable
224  inline exprVar *Var   (int i) const 
225  {return variables_ [i];}
226
227  /// Return vector of variables (symbolic representation)
228  inline std::vector <exprVar *> &Variables () 
229  {return variables_;}
230
231  /// Return pointer to set for comparisons
232  inline std::set <exprAux *, compExpr> *& AuxSet () 
233  {return auxSet_;}
234
235  /// Return pointer to dependence graph
236  inline DepGraph *getDepGraph () 
237  {return graph_;}
238
239  /// return current point & bounds
240  inline Domain *domain () const
241  {return &domain_;}
242
243  // Get and set current variable and bounds
244  inline CouNumber   &X     (int i) const {return domain_.x   (i);} ///< \f$x_i\f$
245  inline CouNumber   &Lb    (int i) const {return domain_.lb  (i);} ///< lower bound on \f$x_i\f$
246  inline CouNumber   &Ub    (int i) const {return domain_.ub  (i);} ///< upper bound on\f$x_i\f$
247
248  // get and set current variable and bounds
249  inline CouNumber  *X     () const {return domain_.();} ///< Return vector of variables
250  inline CouNumber  *Lb    () const {return domain_.lb ();} ///< Return vector of lower bounds
251  inline CouNumber  *Ub    () const {return domain_.ub ();} ///< Return vector of upper bounds
252
253  // get optimal solution and objective value
254  CouNumber  *&bestSol () const {return optimum_;} ///< Best known solution (read from file)
255  CouNumber    bestObj () const {return bestObj_;} ///< Objective of best known solution
256
257  /// Get vector of commuted variables
258  bool *&Commuted () 
259  {return commuted_;}
260
261  /// Add (non linear) objective function
262  void addObjective     (expression *, const std::string &);
263
264  // Add (non linear) "=", ">=", "<=", and range constraints
265  void addEQConstraint  (expression *, expression *); ///< Add equality constraint \f$ h(x) = b\f$
266  void addGEConstraint  (expression *, expression *); ///< Add \f$\ge\f$ constraint, \f$h(x)\ge b\f$
267  void addLEConstraint  (expression *, expression *); ///< Add \f$\le\f$ constraint, \f$h(x)\le b\f$
268  void addRNGConstraint (expression *, expression *, 
269                         expression *);               ///< Add range constraint, \f$a\le h(x)\le b\f$
270
271  /// Add original variable.
272  ///
273  /// @param isint if true, this variable is integer, otherwise it is
274  /// continuous
275  expression *addVariable (bool isint = false, Domain *d = NULL);
276
277  /// Add auxiliary variable and associate it with expression given as
278  /// argument (used in standardization)
279  exprAux *addAuxiliary (expression *);
280
281  /// preprocess problem in order to extract linear relaxations etc.
282  void reformulate ();
283
284  /// Break problem's nonlinear constraints in simple expressions to
285  /// be convexified later. Return true if problem looks feasible,
286  /// false if proven infeasible.
287  bool standardize ();
288
289  /// Display current representation of problem: objective, linear and
290  /// nonlinear constraints, and auxiliary variables.
291  void print (std::ostream & = std::cout);
292
293#ifdef COIN_HAS_ASL
294  /// Read problem from .nl file using the Ampl Solver Library (ASL)
295  int readnl (const struct ASL *);
296
297  /// Generate a Couenne expression from an ASL expression
298  expression *nl2e (struct expr *, const ASL *asl);
299#endif
300
301  // bound tightening parameters
302  bool doFBBT () const {return doFBBT_;} ///< shall we do Feasibility Based Bound Tightening?
303  bool doRCBT () const {return doRCBT_;} ///< shall we do reduced cost      Bound Tightening?
304  bool doOBBT () const {return doOBBT_;} ///< shall we do Optimality  Based Bound Tightening?
305  bool doABT  () const {return doABT_;}  ///< shall we do Aggressive        Bound Tightening?
306
307  int  logObbtLev () const {return logObbtLev_;} ///< How often shall we do OBBT?
308  int  logAbtLev  () const {return logAbtLev_;}  ///< How often shall we do ABT?
309
310  /// Write nonlinear problem to a .mod file (with lots of defined
311  /// variables)
312  ///
313  /// @param fname Name of the .mod file to be written
314  ///
315  /// @param aux controls the use of auxiliaries. If true, a problem
316  /// is written with auxiliary variables written with their
317  /// associated expression, i.e. \f$w_i = h_i(x,y,w)\f$ and bounds
318  /// \f$l_i \le w_i \le u_i\f$, while if false these constraints are
319  /// written in the form \f$l_i \le h_i (x,y) \le u_i\f$.
320  ///
321  /// Note: if used before standardization, writes original AMPL formulation
322  void writeAMPL (const std::string &fname, bool aux);
323
324  /// Write nonlinear problem to a .gms file
325  ///
326  /// @param fname Name of the .gams file to be written.
327  void writeGAMS (const std::string &fname);
328
329  /// Initialize auxiliary variables and their bounds from original
330  /// variables
331  //void initAuxs (const CouNumber *, const CouNumber *, const CouNumber *);
332  void initAuxs () const;
333
334  /// Get auxiliary variables from original variables
335  void getAuxs (CouNumber *) const;
336
337  /// tighten bounds using propagation, implied bounds and reduced costs
338  bool boundTightening (t_chg_bounds *, 
339                        Bonmin::BabInfo * = NULL) const;
340
341  /// core of the bound tightening procedure
342  bool btCore (t_chg_bounds *chg_bds) const;
343
344  /// Optimality Based Bound Tightening
345  int obbt (const CouenneCutGenerator *cg,
346            const OsiSolverInterface &csi,
347            OsiCuts &cs,
348            const CglTreeInfo &info,
349            Bonmin::BabInfo * babInfo,
350            t_chg_bounds *chg_bds);
351
352  /// aggressive bound tightening. Fake bounds in order to cut
353  /// portions of the solution space by fathoming on
354  /// bounds/infeasibility
355  bool aggressiveBT (Bonmin::OsiTMINLPInterface *nlp,
356                     t_chg_bounds *, 
357                     Bonmin::BabInfo * = NULL) const;
358
359  /// procedure to strengthen variable bounds. Return false if problem
360  /// turns out to be infeasible with given bounds, true otherwise.
361  int redCostBT (const OsiSolverInterface *psi,
362                 t_chg_bounds *chg_bds) const;
363
364  /// "Forward" bound tightening, that is, propagate bound of variable
365  /// \f$x\f$ in an expression \f$w = f(x)\f$ to the bounds of \f$w\f$.
366  int tightenBounds (t_chg_bounds *) const;
367
368  /// "Backward" bound tightening, aka implied bounds.
369  int impliedBounds (t_chg_bounds *) const;
370
371  /// Look for quadratic terms to be used with SDP cuts
372  void fillQuadIndices ();
373
374  /// Fill vector with coefficients of objective function
375  void fillObjCoeff (double *&);
376
377  /// Replace all occurrences of original variable with new aux given
378  /// as argument
379  void auxiliarize (exprVar *, exprVar * = NULL);
380
381  /// Set cutoff
382  void setCutOff (CouNumber cutoff) const;
383
384  /// Set cutoff
385  CouNumber getCutOff () const
386  {return pcutoff_ -> getCutOff ();}
387
388  /// Make cutoff known to the problem
389  void installCutOff () const;
390
391  /// Provide Journalist
392  ConstJnlstPtr Jnlst() const 
393  {return ConstPtr(jnlst_);}
394
395  /// Check if solution is MINLP feasible
396  bool checkNLP (const double *solution, double &obj, bool recompute = false) const;
397
398  /// generate integer NLP point Y starting from fractional solution
399  /// using bound tightening
400  int getIntegerCandidate (const double *xFrac, double *xInt, double *lb, double *ub) const;
401
402  /// Read best known solution from file given in argument
403  bool readOptimum (std::string *fname = NULL);
404
405  /// Read cutoff value (for testing purposes)
406  void readCutoff (const std::string &fname);
407
408  /// Add list of options to be read from file
409  static void registerOptions (Ipopt::SmartPtr <Bonmin::RegisteredOptions> roptions);
410
411  /// standardization of linear exprOp's
412  exprAux *linStandardize (bool addAux, 
413                           CouNumber c0, 
414                           LinMap  &lmap,
415                           QuadMap &qmap);
416
417  /// split a constraint w - f(x) = c into w's index (it is returned)
418  /// and rest = f(x) + c
419  int splitAux (CouNumber, expression *, expression *&, bool *);
420
421  /// translates pair (indices, coefficients) into vector with pointers to variables
422  void indcoe2vector (int *indexL,
423                      CouNumber *coeff,
424                      std::vector <std::pair <exprVar *, CouNumber> > &lcoeff);
425
426  /// translates triplet (indicesI, indicesJ, coefficients) into vector with pointers to variables
427  void indcoe2vector (int *indexI,
428                      int *indexJ,
429                      CouNumber *coeff,
430                      std::vector <quadElem> &qcoeff);
431
432  /// given (expression *) element of sum, returns (coe,ind0,ind1)
433  /// depending on element:
434  ///
435  /// 1) a * x_i ^ 2   ---> (a,i,?)   return COU_EXPRPOW
436  /// 2) a * x_i       ---> (a,i,?)   return COU_EXPRVAR
437  /// 3) a * x_i * x_j ---> (a,i,j)   return COU_EXPRMUL
438  /// 4) a             ---> (a,?,?)   return COU_EXPRCONST
439  ///
440  /// x_i and/or x_j may come from standardizing other (linear or
441  /// quadratic operator) sub-expressions
442  void decomposeTerm (expression *term,
443                      CouNumber initCoe,
444                      CouNumber &c0,
445                      LinMap  &lmap,
446                      QuadMap &qmap);
447
448  /// return problem name
449  const std::string &problemName () const
450  {return problemName_;}
451
452  /// return inverse dependence structure
453  const std::vector <std::set <int> > &Dependence () const
454  {return dependence_;}
455
456  /// return object vector
457  const std::vector <CouenneObject> &Objects () const
458  {return objects_;}
459
460  /// find SOS constraints in problem
461  int findSOS (OsiSolverInterface *solver, OsiObject ** objects);
462
463  /// set maximum CPU time
464  inline void setMaxCpuTime (double time)
465  {maxCpuTime_ = time;}
466
467  /// return maximum CPU time
468  inline double getMaxCpuTime () const
469  {return maxCpuTime_;}
470
471  /// save CouenneBase
472  inline void setBase (Bonmin::BabSetupBase *base) {
473    bonBase_ = base;
474    jnlst_   = base -> journalist ();
475  }
476
477  /// Some originals may be unused due to their zero multiplicity
478  /// (that happens when they are duplicates). This procedure creates
479  /// a structure for quickly checking and restoring their value after
480  /// solving.
481  void createUnusedOriginals ();
482
483  /// Some originals may be unused due to their zero multiplicity (that
484  /// happens when they are duplicates). This procedure restores their
485  /// value after solving
486  void restoreUnusedOriginals (CouNumber * = NULL) const;
487
488  /// return indices of neglected redundant variables
489  int *unusedOriginalsIndices () {return unusedOriginalsIndices_;}
490
491  /// number of unused originals
492  int nUnusedOriginals ()        {return nUnusedOriginals_;}
493
494protected:
495
496  /// single fake tightening. Return
497  ///
498  /// -1   if infeasible
499  ///  0   if no improvement
500  /// +1   if improved
501  int fake_tighten (char direction,  ///< 0: left, 1: right
502                    int index,       ///< index of the variable tested
503                    const double *X, ///< point round which tightening is done
504                    CouNumber *olb,  ///< cur. lower bound
505                    CouNumber *oub,  ///< cur. upper bound
506                    t_chg_bounds *chg_bds,
507                    t_chg_bounds *f_chg) const;
508
509  /// Optimality Based Bound Tightening -- inner loop
510  int obbtInner (CouenneSolverInterface *, 
511                 OsiCuts &,
512                 t_chg_bounds *, 
513                 Bonmin::BabInfo *) const;
514
515  int obbt_iter (CouenneSolverInterface *csi, 
516                 t_chg_bounds *chg_bds, 
517                 const CoinWarmStart *warmstart, 
518                 Bonmin::BabInfo *babInfo,
519                 double *objcoe,
520                 int sense, 
521                 int index) const;
522
523  int call_iter (CouenneSolverInterface *csi, 
524                 t_chg_bounds *chg_bds, 
525                 const CoinWarmStart *warmstart, 
526                 Bonmin::BabInfo *babInfo,
527                 double *objcoe,
528                 enum nodeType type,
529                 int sense) const;
530
531  /// analyze sparsity of potential exprQuad/exprGroup and change
532  /// linear/quadratic maps accordingly, if necessary by adding new
533  /// auxiliary variables and including them in the linear map
534  void analyzeSparsity (CouNumber, 
535                        LinMap &,
536                        QuadMap &);
537
538  /// re-organizes multiplication and stores indices (and exponents) of
539  /// its variables
540  void flattenMul (expression *mul, 
541                   CouNumber &coe, 
542                   std::map <int, CouNumber> &indices);
543
544  /// clear all spurious variables pointers not referring to the variables_ vector
545  void realign ();
546
547  /// fill dependence_ structure
548  void fillDependence (Bonmin::BabSetupBase *base);
549
550  /// fill freeIntegers_ array
551  void fillIntegerRank () const;
552
553  /// Test fixing of an integer variable (used in getIntegerCandidate())
554  int testIntFix (int index, 
555                  CouNumber xFrac, 
556                  enum fixType *fixed,
557                  CouNumber *xInt,
558                  CouNumber *dualL, CouNumber *dualR,
559                  CouNumber *olb,   CouNumber *oub,
560                  bool patient) const;
561};
562
563#endif
Note: See TracBrowser for help on using the repository browser.