source: tags/release_2002-11-05/PresolveSingleton.cpp

Last change on this file was 50, checked in by ladanyi, 17 years ago

devel-1 merged into HEAD

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 8.9 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4#include <stdio.h>
5#include <math.h>
6#include <strings.h>
7
8#include "CoinHelperFunctions.hpp"
9#include "PresolveMatrix.hpp"
10
11#include "PresolveEmpty.hpp"    // for DROP_COL/DROP_ROW
12#include "PresolveFixed.hpp"
13#include "PresolveSingleton.hpp"
14#include "ClpMessage.hpp"
15
16/*
17 * Transfers singleton row bound information to the corresponding column bounds.
18 * What I refer to as a row singleton would be called a doubleton
19 * in the paper, since my terminology doesn't refer to the slacks.
20 * In terms of the paper, we transfer the bounds of the slack onto
21 * the variable (vii) and then "substitute" the slack out of the problem
22 * (which is a noop).
23 */
24const PresolveAction *
25slack_doubleton_action::presolve(PresolveMatrix *prob,
26                                 const PresolveAction *next,
27                                 bool & notFinished)
28{
29  double *colels        = prob->colels_;
30  int *hrow             = prob->hrow_;
31  CoinBigIndex *mcstrt          = prob->mcstrt_;
32  int *hincol           = prob->hincol_;
33  int ncols             = prob->ncols_;
34
35  double *clo           = prob->clo_;
36  double *cup           = prob->cup_;
37
38  double *rowels        = prob->rowels_;
39  const int *hcol       = prob->hcol_;
40  const CoinBigIndex *mrstrt    = prob->mrstrt_;
41  int *hinrow           = prob->hinrow_;
42  //  int nrows         = prob->nrows_;
43
44  double *rlo           = prob->rlo_;
45  double *rup           = prob->rup_;
46
47  // If rowstat exists then all do
48  unsigned char *rowstat        = prob->rowstat_;
49  //  double *acts      = prob->acts_;
50  double * sol = prob->sol_;
51  //  unsigned char * colstat = prob->colstat_;
52
53  //  const char *integerType = prob->integerType_;
54
55  const double ztolzb   = prob->ztolzb_;
56
57  int numberLook = prob->numberRowsToDo_;
58  int iLook;
59  int * look = prob->rowsToDo_;
60
61  action actions[MAX_SLACK_DOUBLETONS];
62  int nactions = 0;
63  notFinished = false;
64
65  int * fixed_cols = new int [ncols];
66  int nfixed_cols       = 0;
67
68  // this loop can apparently create new singleton rows;
69  // I haven't bothered to detect this.
70  // wasfor (int irow=0; irow<nrows; irow++)
71  for (iLook=0;iLook<numberLook;iLook++) {
72    int irow = look[iLook];
73    if (hinrow[irow] == 1) {
74      int jcol = hcol[mrstrt[irow]];
75      double coeff = rowels[mrstrt[irow]];
76      double lo = rlo[irow];
77      double up = rup[irow];
78      double acoeff = fabs(coeff);
79      //      const bool singleton_col = (hincol[jcol] == 1);
80
81      if (acoeff < ZTOLDP)
82        continue;
83
84      // don't bother with fixed cols
85      if (fabs(cup[jcol] - clo[jcol]) < ztolzb)
86        continue;
87
88      {
89        // put column on stack of things to do next time
90        prob->addCol(jcol);
91        action *s = &actions[nactions];
92        nactions++;
93
94        s->col = jcol;
95        s->clo = clo[jcol];
96        s->cup = cup[jcol];
97
98        s->row = irow;
99        s->rlo = rlo[irow];
100        s->rup = rup[irow];
101
102        s->coeff = coeff;
103      }
104
105      if (coeff < 0.0) {
106        swap(lo, up);
107        lo = -lo;
108        up = -up;
109      }
110     
111      if (lo <= -PRESOLVE_INF)
112        lo = -PRESOLVE_INF;
113      else {
114        lo /= acoeff;
115        if (lo <= -PRESOLVE_INF)
116          lo = -PRESOLVE_INF;
117      }
118
119      if (up > PRESOLVE_INF)
120        up = PRESOLVE_INF;
121      else {
122        up /= acoeff;
123        if (up > PRESOLVE_INF)
124          up = PRESOLVE_INF;
125      }
126     
127      if (clo[jcol] < lo)
128        clo[jcol] = lo;
129       
130      if (cup[jcol] > up) 
131        cup[jcol] = up;
132
133      if (fabs(cup[jcol] - clo[jcol]) < ZTOLDP) {
134        fixed_cols[nfixed_cols++] = jcol;
135      }
136
137      if (lo > up) {
138        if (lo <= up + prob->feasibilityTolerance_) {
139          // If close to integer then go there
140          double nearest = floor(lo+0.5);
141          if (fabs(nearest-lo)<2.0*prob->feasibilityTolerance_) {
142            lo = nearest;
143            up = nearest;
144          } else {
145            lo = up;
146          }
147          clo[jcol] = lo;
148          cup[jcol] = up;
149        } else {
150          prob->status_ |= 1;
151          prob->messageHandler()->message(CLP_PRESOLVE_COLINFEAS,
152                                             prob->messages())
153                                               <<jcol
154                                               <<lo
155                                               <<up
156                                               <<CoinMessageEol;
157          break;
158        }
159      }
160
161#if     DEBUG_PRESOLVE
162      printf("SINGLETON R-%d C-%d\n", irow, jcol);
163#endif
164
165      // eliminate the row entirely from the row rep
166      // rlinks don't seem to be used PRESOLVE_REMOVE_LINK(rlink, irow);
167      hinrow[irow] = 0;
168
169      // just to make things squeeky
170      rlo[irow] = 0.0;
171      rup[irow] = 0.0;
172
173      if (rowstat) {
174        // update solution and basis
175        int basisChoice=0;
176        int numberBasic=0;
177        if (prob->columnIsBasic(jcol))
178          numberBasic++;
179        if (prob->rowIsBasic(irow))
180          numberBasic++;
181        if (sol[jcol]<=clo[jcol]+ztolzb) {
182          sol[jcol] =clo[jcol];
183          prob->setColumnStatus(jcol,PrePostsolveMatrix::atLowerBound);
184        } else if (sol[jcol]>=cup[jcol]-ztolzb) {
185          sol[jcol] =cup[jcol];
186          prob->setColumnStatus(jcol,PrePostsolveMatrix::atUpperBound);
187        } else {
188          basisChoice=1;
189        }
190        if (numberBasic>1||basisChoice)
191          prob->setColumnStatus(jcol,PrePostsolveMatrix::basic);
192      }
193
194      // remove the row from this col in the col rep
195      presolve_delete_from_row(jcol, irow, mcstrt, hincol, hrow, colels);
196
197      if (nactions >= MAX_SLACK_DOUBLETONS) {
198        notFinished=true;
199        break;
200      }
201    }
202  }
203
204  if (nactions) {
205#if     PRESOLVE_SUMMARY
206    printf("SINGLETON ROWS:  %d\n", nactions);
207#endif
208    action *save_actions = new action[nactions];
209    ClpDisjointCopyN(actions, nactions, save_actions);
210    next = new slack_doubleton_action(nactions, save_actions, next);
211
212    if (nfixed_cols)
213      next = make_fixed_action::presolve(prob, fixed_cols, nfixed_cols,
214                                         true, // arbitrary
215                                         next);
216  }
217  delete [] fixed_cols;
218  return (next);
219}
220
221void slack_doubleton_action::postsolve(PostsolveMatrix *prob) const
222{
223  const action *const actions = actions_;
224  const int nactions = nactions_;
225
226  double *colels        = prob->colels_;
227  int *hrow             = prob->hrow_;
228  CoinBigIndex *mcstrt          = prob->mcstrt_;
229  int *hincol           = prob->hincol_;
230  int *link             = prob->link_;
231  //  int ncols         = prob->ncols_;
232
233  double *clo           = prob->clo_;
234  double *cup           = prob->cup_;
235
236  double *rlo           = prob->rlo_;
237  double *rup           = prob->rup_;
238
239  double *sol           = prob->sol_;
240  double *rcosts        = prob->rcosts_;
241
242  double *acts          = prob->acts_;
243  double *rowduals      = prob->rowduals_;
244
245  unsigned char *colstat                = prob->colstat_;
246  //  unsigned char *rowstat            = prob->rowstat_;
247
248  //  char *cdone               = prob->cdone_;
249  char *rdone           = prob->rdone_;
250  CoinBigIndex free_list                = prob->free_list_;
251
252  const double ztolzb   = prob->ztolzb_;
253
254  for (const action *f = &actions[nactions-1]; actions<=f; f--) {
255    int irow = f->row;
256    double lo0 = f->clo;
257    double up0 = f->cup;
258    double coeff = f->coeff;
259    int jcol = f->col;
260
261    /* the column was in the reduced problem */
262    PRESOLVEASSERT(cdone[jcol] && rdone[irow]==DROP_ROW);
263
264    rlo[irow] = f->rlo;
265    rup[irow] = f->rup;
266
267    clo[jcol] = lo0;
268    cup[jcol] = up0;
269
270    acts[irow] = coeff * sol[jcol];
271
272    // copy col to end to make room for new element
273    {
274      CoinBigIndex k = free_list;
275      free_list = link[free_list];
276
277      check_free_list(free_list);
278
279      hrow[k] = irow;
280      colels[k] = coeff;
281      link[k] = mcstrt[jcol];
282      mcstrt[jcol] = k;
283    }
284    hincol[jcol]++;     // right?
285
286    /*
287     * Have to compute rowstat and rowduals, since we are adding the row.
288     * that satisfy complentarity slackness.
289     *
290     * We may also have to modify colstat and rcosts if bounds
291     * have been relaxed.
292     */
293    if (!colstat) {
294      // ????
295      rowduals[irow] = 0.0;
296    } else {
297      if (prob->columnIsBasic(jcol)) {
298        /* variable is already basic, so slack in this row must be basic */
299        prob->setRowStatus(irow,PrePostsolveMatrix::basic);
300       
301        rowduals[irow] = 0.0;
302        /* hence no reduced costs change */
303        /* since the column was basic, it doesn't matter if it is now
304           away from its bounds. */
305        /* the slack is basic and its reduced cost is 0 */
306      } else if ((fabs(sol[jcol] - lo0) <= ztolzb &&
307                  rcosts[jcol] >= 0) ||
308                 
309                 (fabs(sol[jcol] - up0) <= ztolzb &&
310                  rcosts[jcol] <= 0)) {
311        /* up against its bound but the rcost is still ok */
312        prob->setRowStatus(irow,PrePostsolveMatrix::basic);
313       
314        rowduals[irow] = 0.0;
315        /* hence no reduced costs change */
316      } else if (! (fabs(sol[jcol] - lo0) <= ztolzb) &&
317                 ! (fabs(sol[jcol] - up0) <= ztolzb)) {
318        /* variable not marked basic,
319         * so was at a bound in the reduced problem,
320         * but its bounds were tighter in the reduced problem,
321         * so now it has to be made basic.
322         */
323        prob->setColumnStatus(jcol,PrePostsolveMatrix::basic);
324        prob->setRowStatusUsingValue(irow);
325
326        /* choose a value for rowduals[irow] that forces rcosts[jcol] to 0.0 */
327        /* new reduced cost = 0 = old reduced cost - new dual * coeff */
328        rowduals[irow] = rcosts[jcol] / coeff;
329        rcosts[jcol] = 0.0;
330
331        /* this value is provably of the right sign for the slack */
332        /* SHOULD SHOW THIS */
333      } else {
334        /* the solution is at a bound, but the rcost is wrong */
335
336        prob->setColumnStatus(jcol,PrePostsolveMatrix::basic);
337        prob->setRowStatusUsingValue(irow);
338
339        /* choose a value for rowduals[irow] that forcesd rcosts[jcol] to 0.0 */
340        rowduals[irow] = rcosts[jcol] / coeff;
341        rcosts[jcol] = 0.0;
342
343        /* this value is provably of the right sign for the slack */
344        /*rcosts[irow] = -rowduals[irow];*/
345      }
346    }
347
348    rdone[irow] = SLACK_DOUBLETON;
349  }
350  prob->free_list_ = free_list;
351}
352
Note: See TracBrowser for help on using the repository browser.