source: stable/0.99/Bonmin/src/Interfaces/BonOsiTMINLPInterface.hpp @ 1255

Last change on this file since 1255 was 1255, checked in by pbonami, 11 years ago

merge with trunk (revision r1171)

File size: 41.2 KB
Line 
1// (C) Copyright International Business Machines Corporation, Carnegie Mellon University 2004, 2007
2// All Rights Reserved.
3// This code is published under the Common Public License.
4//
5// Authors :
6// Pierre Bonami, Carnegie Mellon University,
7// Carl D. Laird, Carnegie Mellon University,
8// Andreas Waechter, International Business Machines Corporation
9//
10// Date : 12/01/2004
11
12
13#ifndef OsiTMINLPInterface_H
14#define OsiTMINLPInterface_H
15
16#define INT_BIAS 0e-8
17
18#include <string>
19#include <iostream>
20
21#include "OsiSolverInterface.hpp"
22#include "CoinWarmStartBasis.hpp"
23
24#include "BonTMINLP.hpp"
25#include "BonTMINLP2TNLP.hpp"
26#include "BonTNLP2FPNLP.hpp"
27#include "BonTNLPSolver.hpp"
28#include "BonCutStrengthener.hpp"
29//#include "BonRegisteredOptions.hpp"
30
31namespace Bonmin {
32  class RegisteredOptions;
33  class StrongBranchingSolver;
34
35  /** Solvers for solving nonlinear programs.*/
36  enum Solver{
37    EIpopt=0 /** <a href="http://projects.coin-or.org/Ipopt"> Ipopt </a> interior point algorithm.*/,
38    EFilterSQP /** <a href="http://www-unix.mcs.anl.gov/~leyffer/solvers.html"> filterSQP </a> Sequential Quadratic Programming algorithm.*/,
39    EAll/** Use all solvers.*/
40  };
41/**
42   This is class provides an Osi interface for a Mixed Integer Linear Program
43   expressed as a TMINLP
44   (so that we can use it for example as the continuous solver in Cbc).
45*/
46
47class OsiTMINLPInterface : public OsiSolverInterface
48{
49  friend class BonminParam;
50
51public:
52
53  //#############################################################################
54
55  /**Error class to throw exceptions from OsiTMINLPInterface.
56   * Inherited from CoinError, we just want to have a different class to be able to catch
57   * errors thrown by OsiTMINLPInterface.
58  */
59class SimpleError : public CoinError
60  {
61  private:
62    SimpleError();
63
64  public:
65    ///Alternate constructor using strings
66    SimpleError(std::string message,
67        std::string methodName,
68        std::string f = std::string(),
69        int l = -1)
70        :
71        CoinError(message,methodName,std::string("OsiTMINLPInterface"), f, l)
72    {}
73  }
74  ;
75
76#ifdef __LINE__
77#define SimpleError(x, y) SimpleError((x), (y), __FILE__, __LINE__)
78#endif
79
80  // Error when problem is not solved
81  TNLPSolver::UnsolvedError * newUnsolvedError(int num, Ipopt::SmartPtr<TMINLP2TNLP> problem, std::string name){
82    return app_->newUnsolvedError(num, problem, name);
83  }
84  //#############################################################################
85
86
87  /** Type of the messages specifically outputed by OsiTMINLPInterface.*/
88  enum MessagesTypes{
89    SOLUTION_FOUND/**found a feasible solution*/,
90    INFEASIBLE_SOLUTION_FOUND/**found an infeasible problem*/,
91    UNSOLVED_PROBLEM_FOUND/**found an unsolved problem*/,
92    WARNING_RESOLVING /** Warn that a problem is resolved*/,
93    WARN_SUCCESS_WS/** Problem not solved with warm start but solved without*/,
94    WARN_SUCCESS_RANDOM/** Subproblem not solve with warm start but solved with random point*/,
95    WARN_CONTINUING_ON_FAILURE/** a failure occured but is continuing*/,
96    SUSPECT_PROBLEM/** Output the number of the problem.*/,
97    SUSPECT_PROBLEM2/** Output the number of the problem.*/,
98    IPOPT_SUMMARY /** Output summary statistics on Ipopt solution.*/,
99    BETTER_SOL /** Found a better solution with random values.*/,
100    LOG_HEAD/** Head of "civilized" log.*/,
101    LOG_FIRST_LINE/** First line (first solve) of log.*/,
102    LOG_LINE/**standard line (retry solving) of log.*/,
103    ALTERNATE_OBJECTIVE/** Recomputed integer feasible with alternate objective function*/,
104    WARN_RESOLVE_BEFORE_INITIAL_SOLVE /** resolve() has been called but there
105                                              was no previous call to initialSolve().
106                                         */,
107    ERROR_NO_TNLPSOLVER /** Trying to access non-existent TNLPSolver*/,
108    WARNING_NON_CONVEX_OA /** Warn that there are equality or ranged constraints and OA may works bad.*/,
109    SOLVER_DISAGREE_STATUS /** Different solver gives different status for problem.*/,
110    SOLVER_DISAGREE_VALUE /** Different solver gives different optimal value for problem.*/,
111    OSITMINLPINTERFACE_DUMMY_END
112  };
113
114  //#############################################################################
115
116
117  /** Messages outputed by an OsiTMINLPInterface. */
118class Messages : public CoinMessages
119  {
120  public:
121    /// Constructor
122    Messages();
123  };
124
125
126  //#############################################################################
127
128
129  /**@name Constructors and destructors */
130  //@{
131  /// Default Constructor
132  OsiTMINLPInterface();
133
134  /** Facilitator to initialize interface. */
135  void initialize(Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions,
136                  Ipopt::SmartPtr<Ipopt::OptionsList> options,
137                  Ipopt::SmartPtr<Ipopt::Journalist> journalist_,
138                  Ipopt::SmartPtr<TMINLP> tminlp);
139
140  /** Set the model to be solved by interface.*/
141  void setModel(Ipopt::SmartPtr<TMINLP> tminlp);
142  /** Set the solver to be used by interface.*/
143  void setSolver(Ipopt::SmartPtr<TNLPSolver> app);
144  /** Sets the TMINLP2TNLP to be used by the interface.*/
145  void use(Ipopt::SmartPtr<TMINLP2TNLP> tminlp2tnlp){
146     problem_ = tminlp2tnlp;}
147  /** Copy constructor
148  */
149  OsiTMINLPInterface (const OsiTMINLPInterface &);
150
151  /** Virtual copy constructor */
152  OsiSolverInterface * clone(bool copyData = true) const;
153
154  /// Assignment operator
155  OsiTMINLPInterface & operator=(const OsiTMINLPInterface& rhs);
156
157  /// Destructor
158  virtual ~OsiTMINLPInterface ();
159
160
161  /// Read parameter file
162  void readOptionFile(const std::string & fileName);
163
164  /// Retrieve OsiTMINLPApplication option list
165  const Ipopt::SmartPtr<Ipopt::OptionsList> options() const;
166  /// Retrieve OsiTMINLPApplication option list
167  Ipopt::SmartPtr<Ipopt::OptionsList> options();
168
169  //---------------------------------------------------------------------------
170  /**@name Solve methods */
171  //@{
172  /// Solve initial continuous relaxation
173  virtual void initialSolve();
174
175  /** Resolve the continuous relaxation after problem modification.
176      initialSolve may or may not have been called before this is called. In
177      any case, this must solve the problem, and speed the process up if it
178      can reuse any remnants of data that might exist from a previous solve.
179   */
180  virtual void resolve();
181
182  /** Resolve the problem with different random starting points to try to find
183      a better solution (only makes sense for a non-convex problem.*/
184  virtual void resolveForCost(int numretry, bool keepWs);
185
186  /** Method to be called when a problem has failed to be solved. Will try
187      to resolve it with different settings.
188  */
189  virtual void resolveForRobustness(int numretry);
190
191  /// Nescessary for compatibility with OsiSolverInterface but does nothing.
192  virtual void branchAndBound()
193  {
194    throw SimpleError("Function not implemented for OsiTMINLPInterface","branchAndBound()");
195  }
196  //@}
197
198
199
200  //---------------------------------------------------------------------------
201  ///@name Methods returning info on how the solution process terminated
202  //@{
203  /// Are there a numerical difficulties?
204  virtual bool isAbandoned() const;
205  /// Is optimality proven?
206  virtual bool isProvenOptimal() const;
207  /// Is primal infeasiblity proven?
208  virtual bool isProvenPrimalInfeasible() const;
209  /// Is dual infeasiblity proven?
210  virtual bool isProvenDualInfeasible() const;
211  /// Is the given primal objective limit reached?
212  virtual bool isPrimalObjectiveLimitReached() const;
213  /// Is the given dual objective limit reached?
214  virtual bool isDualObjectiveLimitReached() const;
215  /// Iteration limit reached?
216  virtual bool isIterationLimitReached() const;
217
218  ///Warn solver that branch-and-bound is continuing after a failure
219  void continuingOnAFailure()
220  {
221    hasContinuedAfterNlpFailure_ = true;
222  }
223  /// Did we continue on a failure
224  bool hasContinuedOnAFailure()
225  {
226    return hasContinuedAfterNlpFailure_;
227  }
228  /// tell to ignore the failures (don't throw, don't fathom, don't report)
229  void ignoreFailures()
230  {
231    pretendFailIsInfeasible_ = 2;
232  }
233  /// Force current solution to be infeasible
234  void forceInfeasible()
235  {
236    problem_->set_obj_value(1e200);
237  }
238  /// Force current solution to be branched on (make it fractionnal with small objective)
239  void forceBranchable()
240  {
241    problem_->set_obj_value(-1e200);
242    problem_->force_fractionnal_sol();
243  }
244  //@}
245
246
247  //---------------------------------------------------------------------------
248  /**@name Parameter set/get methods
249
250     The set methods return true if the parameter was set to the given value,
251     false otherwise. There can be various reasons for failure: the given
252     parameter is not applicable for the solver (e.g., refactorization
253     frequency for the clp algorithm), the parameter is not yet implemented
254     for the solver or simply the value of the parameter is out of the range
255     the solver accepts. If a parameter setting call returns false check the
256     details of your solver.
257
258     The get methods return true if the given parameter is applicable for the
259     solver and is implemented. In this case the value of the parameter is
260     returned in the second argument. Otherwise they return false.
261  */
262  //@{
263  // Set an integer parameter
264  bool setIntParam(OsiIntParam key, int value);
265  // Set an double parameter
266  bool setDblParam(OsiDblParam key, double value);
267  // Set a string parameter
268  bool setStrParam(OsiStrParam key, const std::string & value);
269  // Get an integer parameter
270  bool getIntParam(OsiIntParam key, int& value) const;
271  // Get an double parameter
272  bool getDblParam(OsiDblParam key, double& value) const;
273  // Get a string parameter
274  bool getStrParam(OsiStrParam key, std::string& value) const;
275
276  // Get the push values for starting point
277  inline double getPushFact() const
278  {
279    return pushValue_;
280  }
281
282  //@}
283
284
285  //---------------------------------------------------------------------------
286  /**@name Problem information methods
287
288     These methods call the solver's query routines to return
289     information about the problem referred to by the current object.
290     Querying a problem that has no data associated with it result in
291     zeros for the number of rows and columns, and NULL pointers from
292     the methods that return vectors.
293
294     Const pointers returned from any data-query method are valid as
295     long as the data is unchanged and the solver is not called.
296  */
297  //@{
298  /**@name Methods related to querying the input data */
299  //@{
300  /// Get number of columns
301  virtual int getNumCols() const;
302
303  /// Get number of rows
304  virtual int getNumRows() const;
305
306  ///get name of variables
307  const OsiSolverInterface::OsiNameVec& getVarNames() ;
308  /// Get pointer to array[getNumCols()] of column lower bounds
309  virtual const double * getColLower() const;
310
311  /// Get pointer to array[getNumCols()] of column upper bounds
312  virtual const double * getColUpper() const;
313
314  /** Get pointer to array[getNumRows()] of row constraint senses.
315      <ul>
316      <li>'L': <= constraint
317      <li>'E': =  constraint
318      <li>'G': >= constraint
319      <li>'R': ranged constraint
320      <li>'N': free constraint
321      </ul>
322  */
323  virtual const char * getRowSense() const;
324
325  /** Get pointer to array[getNumRows()] of rows right-hand sides
326      <ul>
327      <li> if rowsense()[i] == 'L' then rhs()[i] == rowupper()[i]
328      <li> if rowsense()[i] == 'G' then rhs()[i] == rowlower()[i]
329      <li> if rowsense()[i] == 'R' then rhs()[i] == rowupper()[i]
330      <li> if rowsense()[i] == 'N' then rhs()[i] == 0.0
331      </ul>
332  */
333  virtual const double * getRightHandSide() const;
334
335  /** Get pointer to array[getNumRows()] of row ranges.
336      <ul>
337      <li> if rowsense()[i] == 'R' then
338      rowrange()[i] == rowupper()[i] - rowlower()[i]
339      <li> if rowsense()[i] != 'R' then
340      rowrange()[i] is 0.0
341      </ul>
342  */
343  virtual const double * getRowRange() const;
344
345  /// Get pointer to array[getNumRows()] of row lower bounds
346  virtual const double * getRowLower() const;
347
348  /// Get pointer to array[getNumRows()] of row upper bounds
349  virtual const double * getRowUpper() const;
350
351  /** Get objective function sense (1 for min (default), -1 for max)
352   * Always minimizes */
353  virtual double getObjSense() const
354  {
355    return 1;
356  }
357
358  /// Return true if column is continuous
359  virtual bool isContinuous(int colNumber) const;
360
361  /// Return true if column is binary
362  virtual bool isBinary(int columnNumber) const;
363
364  /** Return true if column is integer.
365      Note: This function returns true if the the column
366      is binary or a general integer.
367  */
368  virtual bool isInteger(int columnNumber) const;
369
370  /// Return true if column is general integer
371  virtual bool isIntegerNonBinary(int columnNumber) const;
372
373  /// Return true if column is binary and not fixed at either bound
374  virtual bool isFreeBinary(int columnNumber) const;
375
376  /// Get solver's value for infinity
377  virtual double getInfinity() const;
378
379  ///Get priorities on integer variables.
380  const int * getPriorities() const
381  {
382    const TMINLP::BranchingInfo * branch = tminlp_->branchingInfo();
383    if(branch)
384      return branch->priorities;
385    else return NULL;
386  }
387  ///get prefered branching directions
388  const int * getBranchingDirections() const
389  {
390    const TMINLP::BranchingInfo * branch = tminlp_->branchingInfo();
391    if(branch)
392      return branch->branchingDirections;
393    else return NULL;
394  }
395  const double * getUpPsCosts() const
396  {
397    const TMINLP::BranchingInfo * branch = tminlp_->branchingInfo();
398    if(branch)
399    return branch->upPsCosts;
400    else return NULL;
401  }
402  const double * getDownPsCosts() const
403  {
404    const TMINLP::BranchingInfo * branch = tminlp_->branchingInfo();
405    if(branch)
406    return branch->downPsCosts;
407    else return NULL;
408  }
409
410
411  //@}
412
413  /**@name Methods related to querying the solution */
414  //@{
415  /// Get pointer to array[getNumCols()] of primal solution vector
416  virtual const double * getColSolution() const;
417
418  /// Get pointer to array[getNumRows()] of dual prices
419  virtual const double * getRowPrice() const;
420
421  /// Get a pointer to array[getNumCols()] of reduced costs
422  virtual const double * getReducedCost() const;
423
424  /** Get pointer to array[getNumRows()] of row activity levels (constraint
425      matrix times the solution vector */
426  virtual const double * getRowActivity() const;
427
428
429  /** Get how many iterations it took to solve the problem (whatever
430      "iteration" mean to the solver.
431      * \todo Figure out what it could mean for Ipopt.
432      */
433  virtual int getIterationCount() const;
434
435  /** get total number of calls to solve.*/
436  int nCallOptimizeTNLP()
437  {
438    return nCallOptimizeTNLP_;
439  }
440  /** get total time taken to solve NLP's. */
441  double totalNlpSolveTime()
442  {
443    return totalNlpSolveTime_;
444  }
445  /** get total number of iterations */
446  int totalIterations()
447  {
448    return totalIterations_;
449  }
450
451
452  //@}
453  //-------------------------------------------------------------------------
454  /**@name Methods to modify the objective, bounds, and solution
455  */
456  //@{
457
458  /** Set a single column lower bound.
459      Use -getInfinity() for -infinity. */
460  virtual void setColLower( int elementIndex, double elementValue );
461
462  /** Set a single column upper bound.
463      Use getInfinity() for infinity. */
464  virtual void setColUpper( int elementIndex, double elementValue );
465
466  /** Set the lower bounds for all columns
467      array [getNumCols()] is an array of values for the objective.
468  */
469  virtual void setColLower(const double * array);
470
471  /** Set the upper bounds for all columns
472      array [getNumCols()] is an array of values for the objective.
473  */
474  virtual void setColUpper(const double * array);
475
476
477  /** Set a single row lower bound.
478      Use -getInfinity() for -infinity. */
479  virtual void setRowLower( int elementIndex, double elementValue );
480
481  /** Set a single row upper bound.
482      Use getInfinity() for infinity. */
483  virtual void setRowUpper( int elementIndex, double elementValue );
484
485  /** Set the type of a single row */
486  virtual void setRowType(int index, char sense, double rightHandSide,
487      double range);
488
489
490  /** \brief Set the objective function sense (disabled).
491   * (1 for min (default), -1 for max)
492   \todo Make it work.
493   \bug Can not treat maximisation problems. */
494  virtual void setObjSense(double s);
495
496  /** Set the primal solution variable values
497      Set the values for the starting point.
498      \warning getColSolution will never return this vector (unless it is optimal).
499  */
500  virtual void setColSolution(const double *colsol);
501
502  /** Set dual solution variable values.
503      set the values for the starting point.
504      \warning getRowPrice will never return this vector (unless it is optimal).
505  */
506  virtual void setRowPrice(const double * rowprice);
507
508  //@}
509
510
511  //---------------------------------------------------------------------------
512  /**@name WarmStart related methods (those should really do nothing for the moment)*/
513  //@{
514
515  /*! \brief Get an empty warm start object
516
517  This routine returns an empty CoinWarmStartBasis object. Its purpose is
518  to provide a way to give a client a warm start basis object of the
519  appropriate type, which can resized and modified as desired.
520  */
521  virtual CoinWarmStart *getEmptyWarmStart () const;
522
523  /** Get warmstarting information */
524  virtual CoinWarmStart* getWarmStart() const;
525
526  /** Set warmstarting information. Return true/false depending on whether
527      the warmstart information was accepted or not. */
528  virtual bool setWarmStart(const CoinWarmStart* warmstart);
529
530  void setExposeWarmStart(bool value) {
531    exposeWarmStart_ = value;
532  }
533
534  bool getExposeWarmStart() {
535    return exposeWarmStart_;
536  }
537
538  void randomStartingPoint();
539
540  //Returns true if a basis is available
541  virtual bool basisIsAvailable() const
542  {
543    // Throw an exception
544    throw SimpleError("Needs coding for this interface", "basisIsAvailable");
545  }
546
547
548  //@}
549
550  //-------------------------------------------------------------------------
551  /**@name Methods to set variable type */
552  //@{
553  /** Set the index-th variable to be a continuous variable */
554  virtual void setContinuous(int index);
555  /** Set the index-th variable to be an integer variable */
556  virtual void setInteger(int index);
557  //@}
558
559  //Set numIterationSuspect_
560  void setNumIterationSuspect(int value)
561  {
562    numIterationSuspect_ = value;
563  }
564
565  /**@name Dummy functions
566   * Functions which have to be implemented in an OsiSolverInterface,
567   * but which do not do anything (but throwing exceptions) here in the case of a
568   * minlp solved using an nlp solver for continuous relaxations */
569  //@{
570
571  /** Cbc will understand that no matrix exsits if return -1
572  */
573  virtual int getNumElements() const
574  {
575    return -1;
576  }
577
578
579  /** This returns the objective function gradient at the current
580   *  point.  It seems to be required for Cbc's pseudo cost
581   *  initialization
582  */
583  virtual const double * getObjCoefficients() const;
584
585  /** We have to keep this but it will return NULL.
586   */
587  virtual const CoinPackedMatrix * getMatrixByRow() const
588  {
589      return NULL;
590  }
591
592
593  /** We have to keep this but it will return NULL.
594   */
595  virtual const CoinPackedMatrix * getMatrixByCol() const
596  {
597      return NULL;
598  }
599
600  /** We have to keep this but it will throw an error.
601  */
602  virtual void setObjCoeff( int elementIndex, double elementValue )
603  {
604    throw SimpleError("OsiTMINLPInterface does not implement this function.",
605        "setObjCoeff");
606  }
607
608  /** We have to keep this but it will throw an error.
609  */
610  virtual void addCol(const CoinPackedVectorBase& vec,
611      const double collb, const double colub,
612      const double obj)
613  {
614    throw SimpleError("OsiTMINLPInterface does not implement this function.",
615        "addCol");
616  }
617  /** We have to keep this but it will throw an error.
618  */
619  virtual void deleteCols(const int num, const int * colIndices)
620  {
621    throw SimpleError("OsiTMINLPInterface does not implement this function.",
622        "deleteCols");
623  }
624
625  /** We have to keep this but it will throw an error.
626  */
627  virtual void addRow(const CoinPackedVectorBase& vec,
628      const double rowlb, const double rowub)
629  {
630    throw SimpleError("OsiTMINLPInterface does not implement this function.",
631        "addRow");
632  }
633  /** We have to keep this but it will throw an error.
634  */
635  virtual void addRow(const CoinPackedVectorBase& vec,
636      const char rowsen, const double rowrhs,
637      const double rowrng)
638  {
639    throw SimpleError("OsiTMINLPInterface model does not implement this function.",
640        "addRow");
641  }
642  /** We have to keep this but it will throw an error.
643  */
644  virtual void deleteRows(const int num, const int * rowIndices)
645  {
646    if(num)
647      freeCachedRowRim();
648     problem_->removeCuts(num, rowIndices);
649  }
650
651
652  /** We have to keep this but it will throw an error
653  */
654  virtual void loadProblem(const CoinPackedMatrix& matrix,
655      const double* collb, const double* colub,
656      const double* obj,
657      const double* rowlb, const double* rowub)
658  {
659    throw SimpleError("OsiTMINLPInterface does not implement this function.",
660        "loadProblem");
661  }
662
663
664  /** We have to keep this but it will throw an error.
665  */
666  virtual void assignProblem(CoinPackedMatrix*& matrix,
667      double*& collb, double*& colub, double*& obj,
668      double*& rowlb, double*& rowub)
669  {
670    throw SimpleError("OsiTMINLPInterface does not implement this function.",
671        "assignProblem");
672  }
673
674  /** We have to keep this but it will throw an error.
675  */
676  virtual void loadProblem(const CoinPackedMatrix& matrix,
677      const double* collb, const double* colub,
678      const double* obj,
679      const char* rowsen, const double* rowrhs,
680      const double* rowrng)
681  {
682    throw SimpleError("OsiTMINLPInterface does not implement this function.",
683        "loadProblem");
684  }
685
686  /** We have to keep this but it will throw an error.
687  */
688  virtual void assignProblem(CoinPackedMatrix*& matrix,
689      double*& collb, double*& colub, double*& obj,
690      char*& rowsen, double*& rowrhs,
691      double*& rowrng)
692  {
693    throw SimpleError("OsiTMINLPInterface does not implement this function.",
694        "assignProblem");
695  }
696
697
698  /** We have to keep this but it will throw an error.
699  */
700  virtual void loadProblem(const int numcols, const int numrows,
701      const int* start, const int* index,
702      const double* value,
703      const double* collb, const double* colub,
704      const double* obj,
705      const double* rowlb, const double* rowub)
706  {
707    throw SimpleError("OsiTMINLPInterface does not implement this function.",
708        "loadProblem");
709  }
710
711  /** We have to keep this but it will throw an error.
712  */
713  virtual void loadProblem(const int numcols, const int numrows,
714      const int* start, const int* index,
715      const double* value,
716      const double* collb, const double* colub,
717      const double* obj,
718      const char* rowsen, const double* rowrhs,
719      const double* rowrng)
720  {
721    throw SimpleError("OsiTMINLPInterface model does not implement this function.",
722        "loadProblem");
723  }
724
725  /** We have to keep this but it will throw an error.
726  */
727  virtual int readMps(const char *filename,
728      const char *extension = "mps")
729  {
730    throw SimpleError("OsiTMINLPInterface does not implement this function.",
731        "readMps");
732  }
733
734
735  /** We have to keep this but it will throw an error.
736  */
737  virtual void writeMps(const char *filename,
738      const char *extension = "mps",
739      double objSense=0.0) const
740  {
741    throw SimpleError("OsiTMINLPInterface does not implement this function.",
742        "writeMps");
743  }
744
745  /** Throws an error */
746  virtual std::vector<double*> getDualRays(int maxNumRays) const
747  {
748    throw SimpleError("OsiTMINLPInterface does not implement this function.",
749        "getDualRays");
750  }
751
752  /** Throws an error */
753  virtual std::vector<double*> getPrimalRays(int maxNumRays) const
754  {
755    throw SimpleError("OsiTMINLPInterface does not implement this function.",
756        "getPrimalRays");
757  }
758
759  //@}
760
761
762 
763  //---------------------------------------------------------------------------
764
765
766
767  /**@name Control of Ipopt output
768   */
769  //@{
770  void turnOffSolverOutput(){
771  app_->turnOffOutput();}
772  void turnOnSolverOutput(){
773  app_->turnOnOutput();}
774  //@}
775
776  /**@name Sets and Getss */
777  //@{
778  /// Get objective function value (can't use default)
779  virtual double getObjValue() const;
780
781  //@}
782
783  /** get pointer to the TMINLP2TNLP adapter */
784  const TMINLP2TNLP * problem() const
785  {
786    return GetRawPtr(problem_);
787  }
788
789  TMINLP2TNLP * problem()
790  {
791    return GetRawPtr(problem_);
792  }
793
794  const TMINLP * model() const
795  {
796    return GetRawPtr(tminlp_);
797  }
798 
799  Bonmin::TMINLP * model()
800  {
801    return GetRawPtr(tminlp_);
802  }
803 
804  const Bonmin::TNLPSolver * solver() const
805  {
806    return GetRawPtr(app_);
807  } 
808 
809  TNLPSolver * solver()
810  {
811    return GetRawPtr(app_);
812  } 
813  /** \name Methods to build outer approximations */
814  //@{
815  /** \name Methods to build outer approximations */
816  //@{
817  /** \brief Extract a linear relaxation of the MINLP.
818   * Use user-provided point to build first-order outer-approximation constraints at the optimum.
819   * And put it in an OsiSolverInterface.
820   */
821  virtual void extractLinearRelaxation(OsiSolverInterface &si, const double *x, 
822                                       bool getObj = 1);
823
824  /** \brief Extract a linear relaxation of the MINLP.
825   * Solve the continuous relaxation and takes first-order outer-approximation constraints at the optimum.
826   * The put everything in an OsiSolverInterface.
827   */
828  virtual void extractLinearRelaxation(OsiSolverInterface &si, bool getObj = 1,
829                                       bool solveNlp = 1){
830     if(solveNlp)
831       initialSolve();
832     extractLinearRelaxation(si, getColSolution(), getObj); 
833     if(solveNlp){
834        app_->enableWarmStart();
835        setColSolution(problem()->x_sol());
836        setRowPrice(problem()->duals_sol());
837     }
838   }
839
840  /** Get the outer approximation constraints at the current optimal point.
841      If x2 is different from NULL only add cuts violated by x2.
842   (Only get outer-approximations of nonlinear constraints of the problem.)*/
843  void getOuterApproximation(OsiCuts &cs, bool getObj, const double * x2, bool global)
844{
845  getOuterApproximation(cs, getColSolution(), getObj, x2, global);
846}
847
848  /** Get the outer approximation constraints at provided point.
849      If x2 is different from NULL only add cuts violated by x2.
850   (Only get outer-approximations of nonlinear constraints of the problem.)*/
851  void getOuterApproximation(OsiCuts &cs, const double * x, bool getObj, const double * x2, bool global){
852    getOuterApproximation(cs, x, getObj, x2, 0., global);}
853
854  /** Get the outer approximation constraints at provided point.
855      If x2 is different from NULL only add cuts violated by x2 by more than delta.
856   (Only get outer-approximations of nonlinear constraints of the problem.)*/
857  virtual void getOuterApproximation(OsiCuts &cs, const double * x, bool getObj, const double * x2,
858                                     double theta, bool global);
859
860 /** Get the outer approximation at provided point for given constraint. */
861  virtual void getConstraintOuterApproximation(OsiCuts & cs, int constraintNumber,
862                                               const double * x, 
863                                               const double * x2, bool global);
864
865 /** Get the outer approximation at current optimal point for given constraint. */
866  void getConstraintOuterApproximation(OsiCuts & cs, int constraintNumber,
867                                       const double * x2, bool global){
868     getConstraintOuterApproximation(cs, constraintNumber, getColSolution(),x2,global);
869  }
870
871  /** Get the Benders cut at provided point with provided multipliers.*/
872  void getBendersCut(OsiCuts &cs, const double * x, const double *lambda, bool getObj = 1);
873
874  /** solve the problem of finding the closest point to x_bar in the subspace of coordinates given by ind
875   * (i.e., \f$ min \sum\limits_{i=1}^n (x_{ind[i]} -\overline{x}_i)^2 \f$ ,
876   * and get the corresponding outer-approximation constraints.
877      (Only get outer-approximations of nonlinear constraints of the problem.)
878   * \return Distance between feasibility set and x
879   * \param n number of element in arrays x and ind
880   * \param ind indices of the coordinate*/
881  double getFeasibilityOuterApproximation(int n, const double * x_bar,const int *ind, OsiCuts &cs, bool addOnlyViolated, bool global);
882
883  /** Given a point x_bar this solves the problem of finding the point which minimize a convex
884    *combination between the distance to  x_bar and the original objective function f(x):
885   * \f$ min a * (\sum\limits_{i=1}^n  ||x_{ind[i]} -\overline{x}_i)||_L) + (1 - a)* s *f(x) \f$
886   * \return Distance between feasibility set a x_bar on components in ind
887   * \param n number of elements in array x_bar and ind
888   * \param s scaling of the original objective.
889   * \param a Combination to take between feasibility and original objective (must be between 0 and 1).
890   * \param L L-norm to use (can be either 1 or 2).
891   */
892  double solveFeasibilityProblem(int n, const double * x_bar, const int* ind, double a, double s, int L);
893
894  /** Given a point x_bar this solves the problem of finding the point which minimize
895    * the distance to x_bar while satisfying the additional cutoff constraint:
896   * \f$ min \sum\limits_{i=1}^n  ||x_{ind[i]} -\overline{x}_i)||_L$
897   * \return Distance between feasibility set a x_bar on components in ind
898   * \param n number of elements in array x_bar and ind
899   * \param L L-norm to use (can be either 1 or 2).
900   * \param cutoff objective function value of a known integer feasible solution
901   */
902  double solveFeasibilityProblem(int n, const double * x_bar, const int* ind, int L, double cutoff);
903
904  /** Given a point x_bar setup feasibility problem and switch so that every call to initialSolve or resolve will
905      solve it.*/
906  void switchToFeasibilityProblem(int n, const double * x_bar, const int* ind, double a, double s, int L);
907
908  /** Given a point x_bar setup feasibility problem and switch so that every call to initialSolve or resolve will
909      solve it. This is to be used in the local branching heuristic */
910  void switchToFeasibilityProblem(int n, const double * x_bar, const int* ind,
911                                  double rhs_local_branching_constraint);
912
913  /** switch back to solving original problem.*/
914  void switchToOriginalProblem();
915  //@}
916
917  /** \name output for OA cut generation
918       \todo All OA code here should be moved to a separate class sometime.*/
919  //@{
920  /** OA Messages types.*/
921  enum OaMessagesTypes {
922    CUT_NOT_VIOLATED_ENOUGH = 0/** Says that one cut has been generarted, where from, which is the violation.*/,
923    VIOLATED_OA_CUT_GENERATED/** Cut is not violated enough, give violation.*/,
924    OA_CUT_GENERATED/** Print the cut which has been generated.*/,
925    OA_MESSAGES_DUMMY_END/** Dummy end.*/};
926  /** Class to store OA Messages.*/
927  class OaMessages :public CoinMessages{
928     public:
929     /** Default constructor.*/
930     OaMessages();
931  };
932  /** Like a CoinMessageHandler but can print a cut also.*/
933  class OaMessageHandler : public CoinMessageHandler{
934    public:
935    /** Default constructor.*/
936    OaMessageHandler():CoinMessageHandler(){
937    }
938    /** Constructor to put to file pointer (fp won't be closed).*/
939    OaMessageHandler(FILE * fp):CoinMessageHandler(fp){
940    }
941    /** Destructor.*/
942    virtual ~OaMessageHandler(){
943    }
944    /** Copy constructor.*/
945    OaMessageHandler(const OaMessageHandler &other):
946    CoinMessageHandler(other){}
947    /** Constructor from a regular CoinMessageHandler.*/
948    OaMessageHandler(const CoinMessageHandler &other):
949    CoinMessageHandler(other){}
950    /** Assignment operator.*/
951    OaMessageHandler & operator=(const OaMessageHandler &rhs){
952       CoinMessageHandler::operator=(rhs);
953       return *this;}
954    /** Virtual copy */
955    virtual CoinMessageHandler* clone() const{
956      return new OaMessageHandler(*this);}
957    /** print an OsiRowCut.*/
958    void print(OsiRowCut &row);
959  };
960  void setOaMessageHandler(const CoinMessageHandler &handler){
961    delete oaHandler_;
962    oaHandler_ = new OaMessageHandler(handler);
963  }
964  //@}
965
966    //-----------------------------------------------------------------------
967    /** Apply a collection of cuts.
968    */
969    virtual ApplyCutsReturnCode applyCuts(const OsiCuts & cs,
970                                          double effectivenessLb = 0.0){
971       freeCachedRowRim();
972      problem_->addCuts(cs);
973      ApplyCutsReturnCode rc;
974      return rc;}
975
976   /** Add a collection of linear cuts to problem formulation.*/
977  virtual void applyRowCuts(int numberCuts, const OsiRowCut * cuts);
978
979
980  /** Add a collection of linear cuts to the problem formulation */
981  virtual void applyRowCuts(int numberCuts, const OsiRowCut ** cuts)
982  {
983    if(numberCuts)
984      freeCachedRowRim();
985    problem_->addCuts(numberCuts, cuts);
986  }
987
988 /** Get infinity norm of constraint violation for x. Put into
989     obj the objective value of x.*/
990 double getConstraintsViolation(const double * x, double & obj);
991
992  /** Get infinity norm of constraint violation for x and error in objective
993      value where obj is the estimated objective value of x.*/
994  double getNonLinearitiesViolation(const double *x, const double obj);
995
996//---------------------------------------------------------------------------
997
998  void extractInterfaceParams();
999
1000
1001  /** To set some application specific defaults. */
1002  virtual void setAppDefaultOptions(Ipopt::SmartPtr<Ipopt::OptionsList> Options);
1003
1004  /** Register all possible options to Bonmin */
1005  static void registerOptions (Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions);
1006 
1007  Ipopt::SmartPtr<Bonmin::RegisteredOptions> regOptions(){
1008    if(IsValid(app_))
1009      return app_->roptions();
1010    else
1011      return NULL;
1012  }
1013
1014  /** @name Methods related to strong branching */
1015  //@{
1016  /// Set the strong branching solver
1017  void SetStrongBrachingSolver(Ipopt::SmartPtr<StrongBranchingSolver> strong_branching_solver);
1018  /// Create a hot start snapshot of the optimization process.  In our
1019  /// case, we initialize the StrongBrachingSolver.
1020  virtual void markHotStart();
1021  /// Optimize starting from the hot start snapshot. In our case, we
1022  /// call the StrongBranchingSolver to give us an approximate
1023  /// solution for the current state of the bounds
1024  virtual void solveFromHotStart();
1025  /// Delete the hot start snapshot. In our case we deactivate the
1026  /// StrongBrachingSolver.
1027  virtual void unmarkHotStart();
1028  //@}
1029
1030protected:
1031 
1032  //@}
1033
1034  enum RandomGenerationType{
1035    uniform =0, perturb=1, perturb_suffix=2};
1036  /// Initialize data structures for storing the jacobian
1037  int initializeJacobianArrays();
1038
1039  ///@name Virtual callbacks for application specific stuff
1040  //@{
1041  virtual std::string  appName()
1042  {
1043    return "bonmin";
1044  }
1045  //@}
1046  ///@name Protected methods
1047  //@{
1048
1049  /** Call Ipopt to solve or resolve the problem and check for errors.*/
1050  void solveAndCheckErrors(bool doResolve, bool throwOnFailure,
1051      const char * whereFrom);
1052
1053
1054  /** Add a linear cut to the problem formulation.
1055  */
1056  virtual void applyRowCut( const OsiRowCut & rc )
1057  {
1058    const OsiRowCut * cut = &rc;
1059    problem_->addCuts(1, &cut);
1060  }
1061  /** We have to keep this but it will throw an error.
1062  */
1063  virtual void applyColCut( const OsiColCut & cc )
1064  {
1065    throw SimpleError("Ipopt model does not implement this function.",
1066        "applyColCut");
1067  }
1068
1069//  /** Read the name of the variables in an ampl .col file. */
1070//  void readVarNames() const;
1071
1072  //@}
1073
1074  /**@name Model and solver */
1075  //@{
1076  /** TMINLP model.*/
1077  Ipopt::SmartPtr<TMINLP> tminlp_;
1078  /** Adapter for a MINLP to a NLP */
1079  Ipopt::SmartPtr<TMINLP2TNLP> problem_;
1080  /** Problem currently optimized (may be problem_ or feasibilityProblem_)*/
1081  Ipopt::SmartPtr<Ipopt::TNLP> problem_to_optimize_;
1082  /** Is true if and only if in feasibility mode.*/
1083  bool feasibility_mode_;
1084  /** Solver for a TMINLP. */
1085  Ipopt::SmartPtr<TNLPSolver> app_;
1086
1087  /** Alternate solvers for TMINLP.*/
1088  std::list<Ipopt::SmartPtr<TNLPSolver> > debug_apps_;
1089  /** Do we use the other solvers?*/
1090  bool testOthers_;
1091  //@}
1092
1093  /** Warmstart information for reoptimization */
1094  CoinWarmStart* warmstart_;
1095
1096  /**@name Cached information on the problem */
1097  //@{
1098  /** Free cached data relative to variables */
1099  void freeCachedColRim();
1100  /** Free cached data relative to constraints */
1101  void freeCachedRowRim();
1102  /** Free all cached data*/
1103  void freeCachedData();
1104  /** Extract rowsense_ vector rhs_ vector and rowrange_ vector from the lower and upper bounds
1105   *  on the constraints */
1106  void extractSenseRhsAndRange() const;
1107  /// Pointer to dense vector of row sense indicators
1108  mutable char    *rowsense_;
1109
1110  /// Pointer to dense vector of row right-hand side values
1111  mutable double  *rhs_;
1112
1113  /// Pointer to dense vector of slack upper bounds for range constraints (undefined for non-range rows)
1114  mutable double  *rowrange_;
1115  /** Pointer to dense vector of reduced costs
1116      \warning Always 0. with Ipopt*/
1117  mutable double  *reducedCosts_;
1118  /** DualObjectiveLimit is used to store the cutoff in Cbc*/
1119  double OsiDualObjectiveLimit_;
1120  /** does the file variable names exists (will check automatically).*/
1121  mutable bool hasVarNamesFile_;
1122  //@}
1123  /// number of time NLP has been solved
1124  int nCallOptimizeTNLP_;
1125  /// Total solution time of NLP
1126  double totalNlpSolveTime_;
1127  /// toatal number of iterations
1128  int totalIterations_;
1129  /// max radius for random point
1130  double maxRandomRadius_;
1131  /// Method to pick a random starting point.
1132  int randomGenerationType_;
1133  /// Maximum perturbation value
1134  double max_perturbation_;
1135  /// Ipopt value for pushing initial point inside the bounds
1136  double pushValue_;
1137  /// Number of times problem will be resolved in initialSolve (root node)
1138  int numRetryInitial_;
1139  /// Number of times problem will be resolved in resolve
1140  int numRetryResolve_;
1141  /// Number of times infeasible problem will be resolved.
1142  int numRetryInfeasibles_;
1143  /// Number of times problem will be resolved in case of a failure
1144  int numRetryUnsolved_;
1145  /** Messages specific to an OsiTMINLPInterface. */
1146  Messages messages_;
1147  /** If not 0 when a problem is not solved (failed to be solved)
1148      will pretend that it is infeasible. If == 1 will care
1149      (i.e. record the fact issue messages to user), if ==2 don't care (somebody else will) */
1150  int pretendFailIsInfeasible_;
1151  /** did we ever continue optimization ignoring a failure. */
1152  bool hasContinuedAfterNlpFailure_;
1153  /** number iterations above which a problem is considered suspect (-1 is considered \f$+ \infty \f$).
1154        If in a call to solve a problem takes more than that number of iterations it will be outputed to files.*/
1155  int numIterationSuspect_ ;
1156  /** Has problem been optimized since last change (include setColSolution).
1157     If yes getColSolution will return Ipopt point, otherwise will return
1158     initial point.*/
1159  bool hasBeenOptimized_;
1160  /** A fake objective function (all variables to 1) to please Cbc
1161      pseudo costs initialization.  AW: I changed this, it will now be
1162      the objective gradient at current point. */
1163  mutable double * obj_;
1164  /** flag to say wether options have been printed or not.*/
1165  static bool hasPrintedOptions;
1166
1167  /** Adapter for TNLP to a feasibility problem */
1168  Ipopt::SmartPtr<TNLP2FPNLP> feasibilityProblem_;
1169
1170
1171  /** \name Arrays to store Jacobian matrix */
1172  //@{
1173  /** Row indices.*/
1174  int * jRow_;
1175  /** Column indices.*/
1176  int * jCol_;
1177  /** Values */
1178  double * jValues_;
1179  /** Number of elements.*/
1180  int nnz_jac;
1181  //@}
1182
1183  ///Store the types of the constraints (linear and nonlinear).
1184  Ipopt::TNLP::LinearityType * constTypes_;
1185  /** Number of nonlinear constraint
1186   */
1187  int nNonLinear_;
1188  /** Value for small non-zero element which we will try to remove cleanly in OA cuts.*/
1189  double tiny_;
1190  /** Value for small non-zero element which we will take the risk to ignore in OA cuts.*/
1191  double veryTiny_;
1192  /** Value for infinity. */
1193  double infty_;
1194  /** status of last optimization. */
1195  TNLPSolver::ReturnStatus optimizationStatus_;
1196  /** Flag indicating if the warm start methods actually do something.*/
1197  bool exposeWarmStart_;
1198  /** Is it the first solve (for random starting point at root options).*/
1199  bool firstSolve_;
1200  /** Object for strengthening cuts */
1201  SmartPtr<CutStrengthener> cutStrengthener_;
1202
1203  /** \name output for OA cut generation
1204       \todo All OA code here should be moved to a separate class sometime.*/
1205  //@{
1206  /** OA Messages.*/
1207  OaMessages oaMessages_;
1208  /** OA Message handler. */
1209  OaMessageHandler * oaHandler_;
1210  //@}
1211protected:
1212  /** Facilitator to create an application. */
1213  void createApplication(Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions,
1214                         Ipopt::SmartPtr<Ipopt::OptionsList> options,
1215                         Ipopt::SmartPtr<Ipopt::Journalist> journalist);
1216  ///Constructor without model only for derived classes
1217  OsiTMINLPInterface(Ipopt::SmartPtr<TNLPSolver> app);
1218
1219  /** Internal set warm start.*/
1220  bool internal_setWarmStart(const CoinWarmStart* ws);
1221
1222  /** internal get warm start.*/
1223  CoinWarmStart* internal_getWarmStart() const; 
1224private:
1225  /** solver to be used for all strong branching solves */
1226  SmartPtr<StrongBranchingSolver> strong_branching_solver_;
1227  /** status of last optimization before hot start was marked. */
1228  TNLPSolver::ReturnStatus optimizationStatusBeforeHotStart_;
1229static const char * OPT_SYMB;
1230static const char * FAILED_SYMB;
1231static const char * INFEAS_SYMB;
1232static const char * UNBOUND_SYMB;
1233  /** Get status as a char * for log.*/
1234  const char * statusAsString(TNLPSolver::ReturnStatus r){
1235    if(r == TNLPSolver::solvedOptimal || r == TNLPSolver::solvedOptimalTol){
1236      return OPT_SYMB;} 
1237    else if(r == TNLPSolver::provenInfeasible){
1238      return INFEAS_SYMB;}
1239    else if(r == TNLPSolver::unbounded){
1240      return UNBOUND_SYMB;}
1241    else return FAILED_SYMB;
1242  }
1243  const char * statusAsString(){
1244    return statusAsString(optimizationStatus_);}
1245};
1246}
1247#endif
Note: See TracBrowser for help on using the repository browser.