source: trunk/Bonmin/src/Algorithms/OaGenerators/BonOACutGenerator2.cpp @ 1422

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

Remove debug output

File size: 10.2 KB
Line 
1// (C) Copyright Carnegie Mellon University 2005, 2006
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 "BonOACutGenerator2.hpp"
11#include "BonminConfig.h"
12
13#include "OsiClpSolverInterface.hpp"
14
15#include "CbcModel.hpp"
16#include "BonCbcLpStrategy.hpp"
17#ifdef COIN_HAS_CPX
18#include "OsiCpxSolverInterface.hpp"
19#endif
20#include "OsiAuxInfo.hpp"
21#include "BonSolverHelp.hpp"
22
23
24namespace Bonmin
25{
26
27/// Constructor with basic setup
28  OACutGenerator2::OACutGenerator2(BabSetupBase & b):
29      OaDecompositionBase(b, true, false)
30  {
31    int ivalue;
32    std::string bonmin="bonmin.";
33    std::string prefix = (b.prefix() == bonmin) ? "" : b.prefix();
34    prefix += "oa_decomposition.";
35    b.options()->GetEnumValue("milp_solver",ivalue,prefix);
36    if (ivalue <= 0) {//uses cbc
37      //nothing to do?
38    }
39    else if (ivalue == 1) {
40      int nodeS, nStrong, nTrust, mig, mir, probe, cover;
41      b.options()->GetEnumValue("node_comparison",nodeS,prefix);
42      b.options()->GetIntegerValue("number_strong_branch",nStrong,prefix);
43      b.options()->GetIntegerValue("number_before_trust", nTrust,prefix);
44      b.options()->GetIntegerValue("Gomory_cuts", mig,prefix);
45      //b.options()->GetIntegerValue("probing_cuts",probe,prefix);
46      b.options()->GetIntegerValue("mir_cuts",mir,prefix);
47      b.options()->GetIntegerValue("cover_cuts",cover,prefix);
48     
49      CbcStrategy * strategy =
50        new CbcOaStrategy(mig, probe, mir, cover, nTrust,
51            nStrong, nodeS, parameters_.cbcIntegerTolerance_, parameters_.subMilpLogLevel_);
52      setStrategy(*strategy);
53      delete strategy;
54
55    }
56    else if (ivalue == 2) {
57#ifdef COIN_HAS_CPX
58      OsiCpxSolverInterface * cpxSolver = new OsiCpxSolverInterface;
59      b.nonlinearSolver()->extractLinearRelaxation(*cpxSolver);
60      assignLpInterface(cpxSolver);
61#else
62      std::cerr << "You have set an option to use CPLEX as the milp\n"
63      << "subsolver in oa decomposition. However, apparently\n"
64      << "CPLEX is not configured to be used in bonmin.\n"
65      << "See the manual for configuring CPLEX\n";
66      throw -1;
67#endif
68    }
69
70    double oaTime;
71    b.options()->GetNumericValue("time_limit",oaTime,prefix);
72    parameter().maxLocalSearchTime_ =
73    std::min(b.getDoubleParameter(BabSetupBase::MaxTime), oaTime);
74    parameter().localSearchNodeLimit_ = 1000000;
75    parameter().maxLocalSearch_ = 100000;
76    parameter().maxLocalSearchPerNode_ = 10000;
77  }
78  OACutGenerator2::~OACutGenerator2()
79  {}
80
81  /// virutal method to decide if local search is performed
82  bool
83  OACutGenerator2::doLocalSearch(BabInfo * babInfo) const
84  {
85    return (nLocalSearch_<parameters_.maxLocalSearch_ &&
86        parameters_.localSearchNodeLimit_ > 0 &&
87        CoinCpuTime() - timeBegin_ < parameters_.maxLocalSearchTime_ &&
88        numSols_ < parameters_.maxSols_);
89  }
90  /// virtual method which performs the OA algorithm by modifying lp and nlp.
91  double
92  OACutGenerator2::performOa(OsiCuts &cs,
93      solverManip &lpManip,
94      SubMipSolver * &subMip,
95      BabInfo * babInfo,
96      double & cutoff) const
97  {
98
99    bool interuptOnLimit = false;
100    double lastPeriodicLog = CoinCpuTime();
101
102    //const int numcols = nlp_->getNumCols();
103
104
105    bool isInteger = false;
106
107    OsiSolverInterface * lp = lpManip.si();
108    OsiBranchingInformation info(lp, false);
109    bool milpOptimal = 1;
110
111
112    double milpBound = -COIN_DBL_MAX;
113    bool milpFeasible = 1;
114    bool feasible = 1;
115
116    if (subMip)//Perform a local search
117    {
118      subMip->find_good_sol(cutoff, parameters_.subMilpLogLevel_,
119          (parameters_.maxLocalSearchTime_ + timeBegin_ - CoinCpuTime()) /* time limit */,
120          parameters_.localSearchNodeLimit_);
121      milpBound = subMip->lowBound();
122      milpOptimal = subMip->optimal();
123
124      feasible = milpBound < cutoff;
125      milpFeasible = feasible;
126      isInteger = (subMip->getLastSolution() != NULL);
127      nLocalSearch_++;
128
129      if (milpOptimal)
130        handler_->message(SOLVED_LOCAL_SEARCH, messages_)<<subMip->nodeCount()<<subMip->iterationCount()<<CoinMessageEol;
131      else
132      {
133        handler_->message(LOCAL_SEARCH_ABORT, messages_)<<subMip->nodeCount()<<subMip->iterationCount()<<CoinMessageEol;
134      }
135    }
136    int numberPasses = 0;
137
138#ifdef OA_DEBUG
139    bool foundSolution = 0;
140#endif
141    double * nlpSol = NULL;
142
143    while (isInteger && feasible ) {
144      numberPasses++;
145      //after a prescribed elapsed time give some information to user
146      double time = CoinCpuTime();
147      if (time - lastPeriodicLog > parameters_.logFrequency_) {
148        double lb = (subMip == NULL) ?lp->getObjValue() : subMip->getLowerBound();
149        handler_->message(PERIODIC_MSG,messages_)
150        <<time - timeBegin_<<cutoff
151        <<lb
152        <<CoinMessageEol;
153        lastPeriodicLog = CoinCpuTime();
154      }
155
156
157      //setup the nlp
158      int numberCutsBefore = cs.sizeRowCuts();
159
160      //Fix the variable which have to be fixed, after having saved the bounds
161      const double * colsol = subMip == NULL ? lp->getColSolution():
162          subMip->getLastSolution();
163      info.solution_ = colsol;
164
165      fixIntegers(*nlp_,info, parameters_.cbcIntegerTolerance_,objects_, nObjects_);
166
167      nlp_->resolve();
168      if (post_nlp_solve(babInfo, cutoff)) {
169        //nlp solved and feasible
170        // Update the cutoff
171        cutoff = nlp_->getObjValue() *(1 - parameters_.cbcCutoffIncrement_);
172        // Update the lp solver cutoff
173        lp->setDblParam(OsiDualObjectiveLimit, cutoff);
174        numSols_++;
175      }
176
177      nlpSol = const_cast<double *>(nlp_->getColSolution());
178
179      // Get the cuts outer approximation at the current point
180      const double * toCut = (parameter().addOnlyViolated_)?
181          colsol:NULL;
182      nlp_->getOuterApproximation(cs, nlpSol, 1, toCut,
183                                  parameter().global_);
184
185      int numberCuts = cs.sizeRowCuts() - numberCutsBefore;
186      assert(numberCuts);
187      installCuts(*lp, cs, numberCuts);
188
189      lp->resolve();
190
191      double objvalue = lp->getObjValue();
192      //milpBound = max(milpBound, lp->getObjValue());
193      feasible = (lp->isProvenOptimal() &&
194          !lp->isDualObjectiveLimitReached() && (objvalue<cutoff)) ;
195      //if value of integers are unchanged then we have to get out
196      bool changed = !feasible;//if lp is infeasible we don't have to check anything
197      info.solution_ = lp->getColSolution();
198      if (!changed)
199        changed = isDifferentOnIntegers(*nlp_, objects_, nObjects_,
200                                        0.1,
201                                        nlp_->getColSolution(), lp->getColSolution());
202      if (changed) {
203
204        isInteger = integerFeasible(*lp, info, parameters_.cbcIntegerTolerance_,
205                                     objects_, nObjects_);
206      }
207      else {
208        isInteger = 0;
209        //        if(!fixed)//fathom on bounds
210        milpBound = 1e200;
211      }
212#ifdef OA_DEBUG
213      printf("Obj value after cuts %g %d rows\n",lp->getObjValue(),
214          numberCuts) ;
215#endif
216      //check time
217      if (CoinCpuTime() - timeBegin_ > parameters_.maxLocalSearchTime_)
218        break;
219      //do we perform a new local search ?
220      if (feasible && !isInteger &&
221          nLocalSearch_ < parameters_.maxLocalSearch_ &&
222          numberPasses < parameters_.maxLocalSearchPerNode_ &&
223          parameters_.localSearchNodeLimit_ > 0 && numSols_ < parameters_.maxSols_) {
224
225        /** do we have a subMip? if not create a new one. */
226        if (subMip == NULL) subMip = new SubMipSolver(lp, parameters_.strategy());
227
228        nLocalSearch_++;
229
230        subMip->find_good_sol(cutoff, parameters_.subMilpLogLevel_,
231            parameters_.maxLocalSearchTime_ + timeBegin_ - CoinCpuTime(),
232            parameters_.localSearchNodeLimit_);
233
234        milpBound = subMip->lowBound();
235
236        if (subMip->optimal())
237          handler_->message(SOLVED_LOCAL_SEARCH, messages_)<<subMip->nodeCount()<<subMip->iterationCount()<<CoinMessageEol;
238        else
239          handler_->message(LOCAL_SEARCH_ABORT, messages_)<<subMip->nodeCount()<<subMip->iterationCount()<<CoinMessageEol;
240
241
242        colsol = const_cast<double *> (subMip->getLastSolution());
243        isInteger = (colsol != 0);
244
245        feasible =  (milpBound < cutoff);
246
247        if (feasible && isInteger) {
248          bool changed = false;
249          changed = isDifferentOnIntegers(*nlp_, objects_, nObjects_,
250                                          0.1,
251                                          nlp_->getColSolution(), colsol);
252          //solution problem is solved
253          if (!changed) {
254            feasible = 0;
255            milpBound = 1e50;
256          }
257          milpFeasible = feasible;
258        }
259        if (subMip->optimal())
260          milpOptimal = 1;
261        else {
262          milpOptimal = 0;
263        }
264
265        if (milpBound < cutoff)
266          handler_->message(UPDATE_LB, messages_)<<milpBound<<CoinCpuTime() - timeBegin_<<CoinMessageEol;
267        else {
268          milpBound = 1e50;
269          feasible = 0;
270          handler_->message(OASUCCESS, messages_)<<CoinCpuTime() - timeBegin_ <<CoinMessageEol;
271        }
272      }/** endif localSearch*/
273      else if (subMip!=NULL) {
274        delete subMip;
275        subMip = NULL;
276      }
277    }
278
279#ifdef OA_DEBUG
280    debug_.printEndOfProcedureDebugMessage(cs, foundSolution, cutoff, milpBound, isInteger, feasible, std::cout);
281#endif
282    return milpBound;
283  }
284
285  /** Register OA options.*/
286  void
287  OACutGenerator2::registerOptions(Ipopt::SmartPtr<Bonmin::RegisteredOptions> roptions)
288  {
289    roptions->SetRegisteringCategory("Options for OA decomposition", RegisteredOptions::BonminCategory);
290    roptions->AddStringOption2("oa_decomposition", "If yes do initial OA decomposition",
291                               "no",
292                               "no","",
293                               "yes","",
294                               "");
295    roptions->AddBoundedIntegerOption("oa_log_level",
296        "specify OA iterations log level.",
297        0,2,1,
298        "Set the level of output of OA decomposition solver : "
299        "0 - none, 1 - normal, 2 - verbose"
300                                     );
301
302    roptions->AddLowerBoundedNumberOption("oa_log_frequency",
303        "display an update on lower and upper bounds in OA every n seconds",
304        0.,1.,100.,
305        "");
306  }
307}/* End namespace Bonmin. */
308
Note: See TracBrowser for help on using the repository browser.