source: stable/0.1/Couenne/src/problem/CouenneSolverInterface.cpp @ 113

Last change on this file since 113 was 113, checked in by pbelotti, 13 years ago

skip restore unused if there are none (see bug pointed out by C. D'Ambrosio). Restored precedence in if+assignment+OR in impliedBounds-exprSum while (hopefully) avoiding warnings Stefan mentioned.

File size: 10.1 KB
Line 
1/*
2 * Name:    CouenneSolverInterface.cpp
3 * Authors: Pietro Belotti, Carnegie Mellon University
4 *          Andreas Waechter, IBM Corp.
5 * Purpose: Implementation of the OsiSolverInterface::resolve () method
6 *
7 * (C) Carnegie-Mellon University, 2006-09.
8 * This file is licensed under the Common Public License (CPL)
9 */
10
11#include "OsiClpSolverInterface.hpp"
12#include "CouenneProblem.hpp"
13#include "CouenneSolverInterface.hpp"
14#include "CglTreeInfo.hpp"
15
16/// constructor
17CouenneSolverInterface::CouenneSolverInterface (CouenneCutGenerator *cg /*= NULL*/):
18
19  OsiClpSolverInterface(),
20  cutgen_ (cg),
21  knowInfeasible_(false),
22  knowOptimal_(false) {
23
24  // prevents from running OsiClpSolverInterface::tightenBounds()
25  specialOptions_ = specialOptions_ | 262144; 
26}
27
28/// copy constructor
29CouenneSolverInterface::CouenneSolverInterface (const CouenneSolverInterface &src):
30
31  OsiSolverInterface (src),
32  OsiClpSolverInterface (src),
33  cutgen_ (src.cutgen_),
34  knowInfeasible_ (src.knowInfeasible_),
35  knowOptimal_ (src.knowOptimal_) {}
36
37/// Destructor
38CouenneSolverInterface::~CouenneSolverInterface () {
39  //  if (cutgen_)
40  //    delete cutgen_;
41}
42
43
44/// Solve initial LP relaxation
45void CouenneSolverInterface::initialSolve () {
46  /*printf ("------------------------------------- INITIAL SOLVE\n");
47  for (int i=0; i<getNumCols(); i++)
48    printf ("%4d. %20.5g [%20.5g %20.5g]\n",
49            i, getColSolution () [i],
50            getColLower () [i],
51            getColUpper () [i]);
52
53            cutgen_ -> Problem () -> print ();*/
54
55  knowInfeasible_ = 
56  knowOptimal_    = false;
57
58  OsiClpSolverInterface::initialSolve ();
59  //writeLp ("initialLP");
60
61  // some originals may be unused due to their zero multiplicity (that
62  // happens when they are duplicates), restore their value
63  if (cutgen_ -> Problem () -> nUnusedOriginals () > 0) {
64    CouNumber *x = new CouNumber [getNumCols ()];
65    CoinCopyN (getColSolution (), getNumCols (), x);
66    cutgen_ -> Problem () -> restoreUnusedOriginals (x);
67    setColSolution (x);
68    delete [] x;
69  }
70
71  /*
72    printf ("------------------------------------- INITSOLV\n");
73  for (int i=0; i<cutgen_ -> Problem () -> nOrigVars (); i++)
74    //if (cutgen_ -> Problem () -> Var (i) -> Multiplicity () <= 0)
75      {
76      printf ("%4d. %20.5g [%20.5g %20.5g]   ",
77              i, getColSolution () [i],
78              getColLower () [i],
79              getColUpper () [i]);
80      cutgen_ -> Problem () -> Var (i) -> print ();
81      printf ("\n");
82    }
83  */
84}
85
86bool CouenneSolverInterface::isProvenPrimalInfeasible() const {
87  return knowInfeasible_ || OsiClpSolverInterface::isProvenPrimalInfeasible();
88}
89
90bool CouenneSolverInterface::isProvenOptimal() const {
91  return knowOptimal_ || OsiClpSolverInterface::isProvenOptimal();
92}
93
94
95/// Defined in Couenne/src/convex/generateCuts.cpp
96void sparse2dense (int, t_chg_bounds *, int *&, int &);
97
98
99/// Resolve an LP relaxation after problem modification
100void CouenneSolverInterface::resolve () {
101  /*printf ("------------------------------------- RESOLVE\n");
102  for (int i=0; i<getNumCols(); i++)
103    printf ("%4d. %20.5g [%20.5g %20.5g]\n",
104            i, getColSolution () [i],
105            getColLower () [i],
106            getColUpper () [i]);*/
107
108  // CUT THIS! Some about-to-be-resolved problems have variables with
109  // lower +inf. That is, between CbcModel::initialSolve() and
110  // CbcModel::resolve(). I couldn't spot where in Couenne this
111  // happens.
112
113  ////////////////////////////////////// Cut {
114  /*const double
115    *lb = getColLower (),
116    *ub = getColUpper ();
117
118  for (int i=getNumCols(); i--;) {
119    if (lb [i] >  COUENNE_INFINITY)
120      setColLower (i, cutgen_ -> Problem () -> Lb (i));
121    //setColLower (i, -COIN_DBL_MAX);//cutgen_ -> Problem () -> Lb (i));
122    if (ub [i] < -COUENNE_INFINITY)
123      setColUpper (i, cutgen_ -> Problem () -> Ub (i));
124    //setColUpper (i,  COIN_DBL_MAX);//cutgen_ -> Problem () -> Ub (i));
125    }*/
126  ////////////////////////////////////// Cut }
127
128  static int count = -1;
129  char filename [30];
130
131  // save problem to be loaded later
132  if (cutgen_ && (cutgen_ -> check_lp ())) {
133    count++;
134    sprintf (filename, "resolve_%d", count);
135    writeMps (filename);
136  }
137
138  knowInfeasible_ = false;
139  knowOptimal_    = false;
140
141  const CoinWarmStart *ws = NULL;
142
143  if (cutgen_ && (cutgen_ -> check_lp ()))
144    ws = getWarmStart ();
145
146  //deleteScaleFactors ();
147
148  // re-solve problem
149  OsiClpSolverInterface::resolve ();
150
151  // some originals may be unused due to their zero multiplicity (that
152  // happens when they are duplicates), restore their value
153  if (cutgen_ -> Problem () -> nUnusedOriginals () > 0) {
154    CouNumber *x = new CouNumber [getNumCols ()];
155    CoinCopyN (getColSolution (), getNumCols (), x);
156    cutgen_ -> Problem () -> restoreUnusedOriginals (x);
157    setColSolution (x);
158    delete [] x;
159  }
160
161  //cutgen_ -> Problem () -> restoreUnusedOriginals (this);
162
163  // check LP independently
164  if (cutgen_ && (cutgen_ -> check_lp ())) {
165
166    OsiSolverInterface
167      *nsi = new OsiClpSolverInterface,
168      *csi = clone ();
169
170    sprintf (filename, "resolve_%d.mps", count);
171    nsi -> readMps (filename);
172
173    nsi -> messageHandler () -> setLogLevel (0);
174    nsi -> setWarmStart (ws);
175
176    nsi -> initialSolve ();
177
178    if ((nsi -> isProvenOptimal () && isProvenOptimal ()) ||
179        (!(nsi -> isProvenOptimal ()) && !isProvenOptimal ())) {
180
181      if (nsi -> isProvenOptimal () &&
182          (fabs (nsi -> getObjValue () - getObjValue ()) / 
183           (1. + fabs (nsi -> getObjValue ()) + fabs (getObjValue ())) > 1e-2))
184
185        printf ("Warning: discrepancy between saved %g and current %g [%g], file %s\n", 
186                nsi -> getObjValue (),  getObjValue (),
187                nsi -> getObjValue () - getObjValue (),
188                filename);
189    }
190
191    csi -> messageHandler () -> setLogLevel (0);
192    csi -> setWarmStart (ws);
193
194    csi -> initialSolve ();
195
196    if (((csi -> isProvenOptimal () && isProvenOptimal ())) ||
197        (!(csi -> isProvenOptimal ()) && !isProvenOptimal ())) {
198
199      if (csi -> isProvenOptimal () &&
200          (fabs (csi -> getObjValue () - getObjValue ()) / 
201           (1. + fabs (csi -> getObjValue ()) + fabs (getObjValue ())) > 1e-2))
202
203        printf ("Warning: discrepancy between cloned %g and current %g [%g]\n", 
204                csi -> getObjValue (),  getObjValue (),
205                csi -> getObjValue () - getObjValue ());
206    }
207
208    delete nsi;
209    delete csi;
210   
211    delete ws;
212
213    //else printf ("Warning: discrepancy between statuses %s -- %s feasible\n",
214    //filename, isProvenOptimal () ? "current" : "saved");
215  }
216}
217
218
219/// Create a hot start snapshot of the optimization process.
220void CouenneSolverInterface::markHotStart () {
221  //printf(">>>> markHotStart\n");
222  // Using OsiClpSolverInterface doesn't work yet...
223  OsiSolverInterface::markHotStart ();
224}
225
226
227/// Delete the hot start snapshot.
228void CouenneSolverInterface::unmarkHotStart() {
229  //printf("<<<< unmarkHotStart\n");
230  OsiSolverInterface::unmarkHotStart();
231}
232
233
234
235/// Optimize starting from the hot start snapshot.
236void CouenneSolverInterface::solveFromHotStart() {
237
238  //OsiClpSolverInterface::solveFromHotStart ();
239
240  //#if 0
241  knowInfeasible_ = false;
242  knowOptimal_ = false;
243
244  /*
245  const int ncols = cutgen_ -> Problem () -> nVars ();
246
247  cutgen_ -> Problem () -> domain () -> push
248    (cutgen_ -> Problem () -> nVars (),
249     getColSolution (),
250     getColLower    (),
251     getColUpper    ());
252
253  // This vector contains variables whose bounds have changed due to
254  // branching, reduced cost fixing, or bound tightening below. To be
255  // used with malloc/realloc/free
256
257  t_chg_bounds *chg_bds = new t_chg_bounds [ncols];
258
259  OsiCuts cs;
260
261  Bonmin::BabInfo *babInfo = dynamic_cast <Bonmin::BabInfo *> (getAuxiliaryInfo ());
262
263  if (cutgen_ -> Problem () -> doFBBT () &&
264      (!(cutgen_ -> Problem () -> boundTightening (chg_bds, babInfo)))) {
265
266#ifdef DEBUG
267    printf ("#### BT says infeasible before re-solve\n");
268#endif
269    knowInfeasible_ = true;
270    cutgen_ -> Problem () -> domain () -> pop ();
271    return;
272  }
273
274  int *changed = NULL, nchanged;
275  sparse2dense (ncols, chg_bds, changed, nchanged);
276
277  // change tightened bounds through OsiCuts
278  if (nchanged)
279    cutgen_ -> genColCuts (*this, cs, nchanged, changed);
280
281  const int nRowsBeforeRowCuts = getNumRows();
282  //printf("NumRows before getRowCuts = %d\n", getNumRows());
283  cutgen_ -> genRowCuts (*this, cs, nchanged, changed, //CglTreeInfo(),
284                         chg_bds, false);
285
286  // Now go through the list of cuts and apply the column cuts
287  // directly as changes on bounds
288  while(cs.sizeColCuts()) {
289    const OsiColCut& ccut = cs.colCut(0);
290    const CoinPackedVector& lbs = ccut.lbs();
291    int nele = lbs.getNumElements();
292    const int* idxs = lbs.getIndices();
293    const double* eles = lbs.getElements();
294    const double* bnds = getColLower();
295    for (int i=0; i<nele; i++) {
296      if (bnds[*idxs] < *eles) {
297        //printf ("setcolLower %d %g\n", *idxs, *eles);
298        setColLower(*idxs,*eles);
299      }
300      idxs++;
301      eles++;
302    }
303    const CoinPackedVector& ubs = ccut.ubs();
304    nele = ubs.getNumElements();
305    idxs = ubs.getIndices();
306    eles = ubs.getElements();
307    bnds = getColUpper();
308    for (int i=0; i<nele; i++) {
309      if (bnds[*idxs] > *eles) {
310        //printf ("setcolUpper %d %g\n", *idxs, *eles);
311        setColUpper(*idxs,*eles);
312      }
313      idxs++;
314      eles++;
315    }
316    cs.eraseColCut(0);
317  }
318
319  applyCuts (cs);
320
321  const int nRowsAfterRowCuts = getNumRows();
322  //printf("NumRows after applyCuts = %d\n", getNumRows());
323  */
324
325  resolve();
326
327  // some originals may be unused due to their zero multiplicity (that
328  // happens when they are duplicates), restore their value
329  if (cutgen_ -> Problem () -> nUnusedOriginals () > 0) {
330    CouNumber *x = new CouNumber [getNumCols ()];
331    CoinCopyN (getColSolution (), getNumCols (), x);
332    cutgen_ -> Problem () -> restoreUnusedOriginals (x);
333    setColSolution (x);
334    delete [] x;
335  }
336
337  if (isProvenPrimalInfeasible ()) knowInfeasible_ = true;
338  if (isProvenOptimal ())          knowOptimal_    = true;
339
340  //printf("obj value = %e\n",getObjValue());
341
342  // now undo the row cuts
343  /*
344  int nrowsdel = nRowsAfterRowCuts-nRowsBeforeRowCuts;
345  int* rowsdel = new int[nrowsdel];
346  for(int i=0; i<nrowsdel; i++) {
347    rowsdel[i] = nRowsBeforeRowCuts+i;
348  }
349  deleteRows(nrowsdel, rowsdel);
350  delete [] rowsdel;
351  //printf("NumRows after deleting = %d\n", getNumRows());
352
353  cutgen_ -> Problem () -> domain () -> pop ();
354  */
355  //#endif
356}
Note: See TracBrowser for help on using the repository browser.