source: branches/devel/Bonmin/src/OaInterface/BonOACutGenerator.cpp @ 62

Last change on this file since 62 was 62, checked in by pbonami, 13 years ago

astyled the devel branch

  • Property svn:eol-style set to native
  • Property svn:keywords set to "Author Date Id Revision"
File size: 5.6 KB
Line 
1// (C) Copyright Carnegie Mellon University 2005
2// All Rights Reserved.
3// This code is published under the Common Public License.
4//
5// Authors :
6// P. Bonami, Carnegie Mellon University
7//
8// Date :  05/26/2005
9
10#include "BonOACutGenerator.hpp"
11#include "OsiAuxInfo.hpp"
12
13namespace Bonmin
14{
15/// Default constructor
16  OACutGenerator::OACutGenerator(OsiTMINLPInterface * si,
17      int maxDepth)
18      :
19      CglCutGenerator(),
20      nlp_(si),
21      maxDepth_(maxDepth),
22      nSolve_(0)
23  {
24    handler_ = new CoinMessageHandler();
25    handler_ -> setLogLevel(1);
26    messages_ = OaMessages();
27  }
28/// Assign an OsiTMINLPInterface
29  void
30  OACutGenerator::assignInterface(OsiTMINLPInterface * si)
31
32  {
33    nlp_ = si;
34  }
35  static int nCalls = 0;
36/// cut generation method
37  void
38  OACutGenerator::generateCuts( const OsiSolverInterface & si, OsiCuts & cs,
39      const CglTreeInfo info) const
40  {
41    if (nlp_ == NULL) {
42      std::cerr<<"Error in cut generator for outer approximation no ipopt NLP assigned"<<std::endl;
43      throw -1;
44    }
45
46    int numcols = nlp_->getNumCols();
47
48    //Get the continuous solution
49    //const double *colsol = si.getColSolution();
50    //Check for integer feasibility
51    nCalls++;
52    if (info.level > maxDepth_ || info.pass > 0 || !info.inTree  )
53      return;
54    //Fix the variable which have to be fixed, after having saved the bounds
55    double * saveColLb = new double[numcols];
56    double * saveColUb = new double[numcols];
57    CoinCopyN(nlp_->getColLower(), numcols , saveColLb);
58    CoinCopyN(nlp_->getColUpper(), numcols , saveColUb);
59    for (int i = 0; i < numcols ; i++) {
60      if (nlp_->isInteger(i)) {
61        nlp_->setColBounds(i,si.getColLower()[i],si.getColUpper()[i]);
62      }
63    }
64
65#if 0
66//Add tight cuts at LP optimum
67    int numberCuts = si.getNumRows() - nlp_->getNumRows();
68    const OsiRowCut ** cuts = new const OsiRowCut*[numberCuts];
69    int begin = nlp_->getNumRows();
70    numberCuts = 0;
71    int end = si.getNumRows();
72    const double * rowLower = si.getRowLower();
73    const double * rowUpper = si.getRowUpper();
74    const CoinPackedMatrix * mat = si.getMatrixByRow();
75    const CoinBigIndex * starts = mat->getVectorStarts();
76    const int * lengths = mat->getVectorLengths();
77    const double * elements = mat->getElements();
78    const int * indices = mat->getIndices();
79
80    for (int i = begin ; i < end ; i++, numberCuts++) {
81      bool nnzExists=false;
82      for (int k = starts[i] ; k < starts[i]+lengths[i] ; k++) {
83        if (indices[k] == nlp_->getNumCols()) {
84          nnzExists = true;
85          char sign = (elements[k]>0.)?'+':'-';
86          char type='<';
87          if (rowLower[i]>-1e20) type='>';
88          std::cout<<"Non zero with sign: "<<sign<<", type: "<<type<<std::endl;
89        }
90      }
91      if (nnzExists) {
92        numberCuts--;
93        continue;
94      }
95      else
96        std::cout<<"No nonzero element"<<std::endl;
97      int * indsCopy = CoinCopyOfArray(&indices[starts[i]], lengths[i]);
98      double * elemsCopy = CoinCopyOfArray(&elements[starts[i]], lengths[i]);
99      cuts[numberCuts] = new OsiRowCut(rowLower[i], rowUpper[i], lengths[i], lengths[i],
100          indsCopy, elemsCopy);
101    }
102    nlp_->applyRowCuts(numberCuts,cuts);
103#endif
104    //Now solve the NLP get the cuts, reset bounds and get out
105
106    //  nlp_->turnOnIpoptOutput();
107    nSolve_++;
108    nlp_->resolve();
109    nlp_->getOuterApproximation(cs);
110    if (nlp_->isProvenOptimal()) {
111      handler_->message(LP_ERROR,messages_)
112      <<nlp_->getObjValue()-si.getObjValue()<<CoinMessageEol;
113      bool feasible = 1;
114      const double * colsol2 = nlp_->getColSolution();
115      for (int i = 0 ; i < numcols && feasible; i++) {
116        if (nlp_->isInteger(i)) {
117          if (fabs(colsol2[i] - floor(colsol2[i] + 0.5) ) > 1e-07)
118            feasible = 0;
119        }
120      }
121      if (feasible ) {
122#if 1
123        // Also store into solver
124        OsiAuxInfo * auxInfo = si.getAuxiliaryInfo();
125        OsiBabSolver * auxiliaryInfo = dynamic_cast<OsiBabSolver *> (auxInfo);
126        if (auxiliaryInfo) {
127          double * lpSolution = new double[numcols + 1];
128          CoinCopyN(colsol2, numcols, lpSolution);
129          lpSolution[numcols] = nlp_->getObjValue();
130          auxiliaryInfo->setSolution(lpSolution, numcols + 1, lpSolution[numcols]);
131          delete [] lpSolution;
132        }
133        else
134          printf("No auxiliary info in nlp solve!\n");
135#endif
136
137      }
138    }
139    else if (nlp_->isAbandoned() || nlp_->isIterationLimitReached()) {
140      std::cerr<<"Unsolved NLP ... exit"<<std::endl;
141      //    nlp_->turnOnIpoptOutput();
142      nlp_->resolve();
143      throw -1;
144    }
145    else {
146      //   std::cout<<"Infeasible NLP => Infeasible node"<<std::endl;
147      //       //Add an infeasibility local constraint
148      //       CoinPackedVector v;
149      //       double rhs = 1.;
150      //       for(int i = 0; i < numcols ; i++)
151      //        {
152      //          if(nlp_->isInteger(i) && (si.getColUpper()[i] - si.getColLower()[i] < 0.9))
153      //            {
154      //              double value = floor(colsol[i] + 0.5);
155      //              assert(fabs(colsol[i]-value)<1e-8 && value >= -1e-08 && value <= 1+ 1e08);
156      //              v.insert(i, -(2*value - 1));
157      //              rhs -= value;
158      //            }
159      //        }
160      //       OsiRowCut cut;
161      //       cut.setRow(v);
162      //       cut.setLb(rhs);
163      //       cut.setUb(1e300);
164      //       cut.setGloballyValid();
165      //       cs.insert(cut);
166    }
167    for (int i = 0; i < numcols ; i++) {
168      if (nlp_->isInteger(i)) {
169        nlp_->setColBounds(i,saveColLb[i],saveColUb[i]);
170      }
171    }
172#if 0
173    nlp_->deleteLastRows(numberCuts);
174#endif
175    delete [] saveColLb;
176    delete [] saveColUb;
177  }
178}
Note: See TracBrowser for help on using the repository browser.