source: branches/devel-1/PresolveEmpty.cpp @ 2351

Last change on this file since 2351 was 49, checked in by ladanyi, 17 years ago

everything compiles

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 9.4 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 "PresolveMatrix.hpp"
9
10#include "PresolveEmpty.hpp"
11#include "ClpMessage.hpp"
12
13
14const PresolveAction *drop_empty_cols_action::presolve(PresolveMatrix *prob,
15                                                        int *ecols,
16                                                        int necols,
17                                                        const PresolveAction *next)
18{
19  int ncols             = prob->ncols_;
20  CoinBigIndex *mcstrt          = prob->mcstrt_;
21  int *hincol           = prob->hincol_;
22  //  int *hcol         = prob->hcol_;
23
24  // We know we are not going to need row copy again
25  //presolvehlink *clink        = prob->clink_;
26
27  double *clo           = prob->clo_;
28  double *cup           = prob->cup_;
29  double *dcost         = prob->cost_;
30
31  //int *mrstrt         = prob->mrstrt_;
32  //int *hinrow         = prob->hinrow_;
33  //int *hrow           = prob->hrow_;
34
35  char * integerType     = prob->integerType_;
36  int * originalColumn  = prob->originalColumn_;
37
38  const double maxmin   = prob->maxmin_;
39
40  double * sol = prob->sol_;
41  unsigned char * colstat = prob->colstat_;
42
43  action *actions       = new action[necols];
44  int * colmapping = new int [ncols];
45
46  memset(colmapping,0,ncols*sizeof(int));
47  int i;
48  for (i=necols-1; i>=0; i--) {
49    int jcol = ecols[i];
50    colmapping[jcol]=-1;
51    action &e = actions[i];
52
53    e.jcol      = jcol;
54    e.clo       = clo[jcol];
55    e.cup       = cup[jcol];
56    e.cost      = dcost[jcol];
57
58    // there are no more constraints on this variable,
59    // so we had better be able to compute the answer now
60
61    if (dcost[jcol] * maxmin == 0.0) {
62      // hopefully, we can make this non-basic
63      // what does OSL currently do in this case???
64      e.sol = (-PRESOLVE_INF < clo[jcol]
65               ? clo[jcol]
66               : cup[jcol] < PRESOLVE_INF
67               ? cup[jcol]
68               : 0.0);
69    } else if (dcost[jcol] * maxmin > 0.0) {
70      if (-PRESOLVE_INF < clo[jcol])
71        e.sol = clo[jcol];
72      else {
73          prob->messageHandler()->message(CLP_PRESOLVE_COLUMNBOUNDB,
74                                             prob->messages())
75                                               <<jcol
76                                               <<CoinMessageEol;
77        prob->status_ |= 2;
78        break;
79      }
80    } else {
81      if (cup[jcol] < PRESOLVE_INF)
82        e.sol = cup[jcol];
83      else {
84          prob->messageHandler()->message(CLP_PRESOLVE_COLUMNBOUNDA,
85                                             prob->messages())
86                                               <<jcol
87                                               <<CoinMessageEol;
88        prob->status_ |= 2;
89        break;
90      }
91    }
92       
93#if     DEBUG_PRESOLVE
94    if (e.sol * dcost[jcol]) {
95      //printf("\a>>>NON-ZERO COST FOR EMPTY COL %d:  %g\n", jcol, dcost[jcol]);
96    }
97#endif
98    prob->change_bias(e.sol * dcost[jcol]);
99
100
101  }
102  int ncols2=0;
103
104  // now move remaining ones down
105  for (i=0;i<ncols;i++) {
106    if (!colmapping[i]) {
107      mcstrt[ncols2] = mcstrt[i];
108      hincol[ncols2] = hincol[i];
109   
110      clo[ncols2]   = clo[i];
111      cup[ncols2]   = cup[i];
112
113      dcost[ncols2] = dcost[i];
114      if (sol) {
115        sol[ncols2] = sol[i];
116        colstat[ncols2] = colstat[i];
117      }
118
119      integerType[ncols2] = integerType[i];
120      originalColumn[ncols2] = originalColumn[i];
121      ncols2++;
122    }
123  }
124 
125  delete [] colmapping;
126  prob->ncols_ = ncols2;
127
128  return (new drop_empty_cols_action(necols, actions, next));
129}
130
131
132const PresolveAction *drop_empty_cols_action::presolve(PresolveMatrix *prob,
133                                                          const PresolveAction *next)
134{
135  const int *hincol     = prob->hincol_;
136  //  const double *clo = prob->clo_;
137  //  const double *cup = prob->cup_;
138  int ncols             = prob->ncols_;
139  int i;
140  int nempty            = 0;
141  int * empty = new int [ncols];
142
143  // count empty cols
144  for (i=0; i<ncols; i++)
145    if (hincol[i] == 0) {
146#if     DEBUG_PRESOLVE
147      if (nempty==0)
148        printf("UNUSED COLS:  ");
149      printf("%d ", i);
150#endif
151      empty[nempty++] = i;
152    }
153
154  if (nempty) {
155#if     DEBUG_PRESOLVE
156      printf("\ndropped %d cols\n", nempty);
157#endif
158    next = drop_empty_cols_action::presolve(prob, empty, nempty, next);
159  }
160  delete [] empty;
161  return (next);
162}
163
164
165void drop_empty_cols_action::postsolve(PostsolveMatrix *prob) const
166{
167  const int nactions    = nactions_;
168  const action *const actions = actions_;
169
170  int ncols             = prob->ncols_;
171  //  CoinBigIndex nelems               = prob->nelems_;
172
173  CoinBigIndex *mcstrt  = prob->mcstrt_;
174  int *hincol   = prob->hincol_;
175  //  int *hrow = prob->hrow_;
176
177  double *clo   = prob->clo_;
178  double *cup   = prob->cup_;
179
180  double *sol   = prob->sol_;
181  double *cost  = prob->cost_;
182  double *rcosts        = prob->rcosts_;
183  unsigned char *colstat        = prob->colstat_;
184  const double maxmin = prob->maxmin_;
185
186  int ncols2 = ncols+nactions;
187  int * colmapping = new int [ncols2];
188
189  memset(colmapping,0,ncols2*sizeof(int));
190  char *cdone   = prob->cdone_;
191
192  for (int action_i = 0; action_i < nactions; action_i++) {
193    const action *e = &actions[action_i];
194    int jcol = e->jcol;
195    colmapping[jcol]=-1;
196  }
197
198  int i;
199
200  // now move remaining ones up
201  for (i=ncols2-1;i>=0;i--) {
202    if (!colmapping[i]) {
203      ncols--;
204      mcstrt[i] = mcstrt[ncols];
205      hincol[i] = hincol[ncols];
206   
207      clo[i]   = clo[ncols];
208      cup[i]   = cup[ncols];
209
210      cost[i] = cost[ncols];
211
212      sol[i] = sol[ncols];
213     
214      rcosts[i] = rcosts[ncols];
215     
216      if (colstat) 
217        colstat[i] = colstat[ncols];
218      cdone[i]   = cdone[ncols];
219    }
220  }
221  assert (!ncols);
222 
223  delete [] colmapping;
224
225  for (int action_i = 0; action_i < nactions; action_i++) {
226    const action *e = &actions[action_i];
227    int jcol = e->jcol;
228   
229    // now recreate jcol
230    clo[jcol] = e->clo;
231    cup[jcol] = e->cup;
232    sol[jcol] = e->sol;
233    cost[jcol] = e->cost;
234
235    rcosts[jcol] = maxmin*cost[jcol];
236
237    hincol[jcol] = 0;
238
239    cdone[jcol] = DROP_COL;
240    if (colstat) 
241      prob->setColumnStatusUsingValue(jcol);
242  }
243
244  prob->ncols_ += nactions;
245}
246
247
248
249const PresolveAction *drop_empty_rows_action::presolve(PresolveMatrix *prob,
250                                       const PresolveAction *next)
251{
252  int ncols     = prob->ncols_;
253  CoinBigIndex *mcstrt  = prob->mcstrt_;
254  int *hincol   = prob->hincol_;
255  int *hrow     = prob->hrow_;
256
257  int nrows     = prob->nrows_;
258  // This is done after row copy needed
259  //int *mrstrt = prob->mrstrt_;
260  int *hinrow   = prob->hinrow_;
261  //int *hcol   = prob->hcol_;
262 
263  double *rlo   = prob->rlo_;
264  double *rup   = prob->rup_;
265
266  unsigned char *rowstat        = prob->rowstat_;
267  double *acts  = prob->acts_;
268
269  //presolvehlink *rlink = prob->rlink_;
270 
271
272  int i;
273  int nactions = 0;
274  for (i=0; i<nrows; i++)
275    if (hinrow[i] == 0)
276      nactions++;
277
278  if (nactions == 0)
279    return next;
280  else {
281    action *actions     = new action[nactions];
282    int * rowmapping = new int [nrows];
283
284    nactions = 0;
285    int nrows2=0;
286    for (i=0; i<nrows; i++) {
287      if (hinrow[i] == 0) {
288        action &e = actions[nactions];
289        nactions++;
290
291#if     DEBUG_PRESOLVE
292        if (nactions==1)
293          printf("unused rows:  ");
294
295        printf("%d ", i);
296#endif
297        if (rlo[i] > 0.0 || rup[i] < 0.0) {
298          if (fabs(rlo[i])<=prob->feasibilityTolerance_ &&
299              fabs(rup[i])<=prob->feasibilityTolerance_) {
300            rlo[i]=0.0;
301            rup[i]=0.0;
302          } else {
303            prob->status_|= 1;
304          prob->messageHandler()->message(CLP_PRESOLVE_ROWINFEAS,
305                                             prob->messages())
306                                               <<i
307                                               <<rlo[i]
308                                               <<rup[i]
309                                               <<CoinMessageEol;
310            break;
311          }
312        }
313        e.row   = i;
314        e.rlo   = rlo[i];
315        e.rup   = rup[i];
316        rowmapping[i]=-1;
317
318      } else {
319        // move down - we want to preserve order
320        rlo[nrows2]=rlo[i];
321        rup[nrows2]=rup[i];
322        if (acts) {
323          acts[nrows2]=acts[i];
324          rowstat[nrows2]=rowstat[i];
325        }
326        rowmapping[i]=nrows2++;
327      }
328    }
329
330    // remap matrix
331    for (i=0;i<ncols;i++) {
332      int j;
333      for (j=mcstrt[i];j<mcstrt[i]+hincol[i];j++) 
334        hrow[j] = rowmapping[hrow[j]];
335    }
336    delete [] rowmapping;
337
338    prob->nrows_ = nrows2;
339
340#if     DEBUG_PRESOLVE
341    if (nactions)
342      printf("\ndropped %d rows\n", nactions);
343#endif
344    return (new drop_empty_rows_action(nactions, actions, next));
345  }
346}
347
348void drop_empty_rows_action::postsolve(PostsolveMatrix *prob) const
349{
350  const int nactions    = nactions_;
351  const action *const actions = actions_;
352
353  int ncols     = prob->ncols_;
354  CoinBigIndex *mcstrt  = prob->mcstrt_;
355  int *hincol   = prob->hincol_;
356
357  int *hrow     = prob->hrow_;
358
359  double *rlo   = prob->rlo_;
360  double *rup   = prob->rup_;
361  unsigned char *rowstat        = prob->rowstat_;
362  double *rowduals = prob->rowduals_;
363  double *acts  = prob->acts_;
364  char *rdone   = prob->rdone_;
365
366  int nrows0    = prob->nrows0_;
367  int nrows     = prob->nrows_;
368
369  int * rowmapping = new int [nrows0];
370
371  int i, action_i;
372  for (i=0; i<nrows0; i++)
373    rowmapping[i] = 0;
374 
375  for (action_i = 0; action_i<nactions; action_i++) {
376    const action *e = &actions[action_i];
377    int hole = e->row;
378    rowmapping[hole]=-1;
379  }
380
381  // move data
382  for (i=nrows0-1; i>=0; i--) {
383    if (!rowmapping[i]) {
384      // not a hole
385      nrows--;
386      rlo[i]=rlo[nrows];
387      rup[i]=rup[nrows];
388      acts[i]=acts[nrows];
389      rowduals[i]=rowduals[nrows];
390      if (rowstat)
391        rowstat[i] = rowstat[nrows];
392    }
393  }
394  assert (!nrows);
395  // set up mapping for matrix
396  for (i=0;i<nrows0;i++) {
397    if (!rowmapping[i])
398      rowmapping[nrows++]=i;
399  }
400
401  for (int j=0; j<ncols; j++) {
402    CoinBigIndex start = mcstrt[j];
403    CoinBigIndex end   = start + hincol[j];
404   
405    for (CoinBigIndex k=start; k<end; ++k) {
406      hrow[k] = rowmapping[hrow[k]];
407    }
408  }
409
410
411  delete [] rowmapping;
412
413  for (action_i = 0; action_i < nactions; action_i++) {
414    const action *e = &actions[action_i];
415    int irow = e->row;
416
417    // Now recreate irow
418    rlo[irow] = e->rlo;
419    rup[irow] = e->rup;
420
421    if (rowstat)
422      prob->setRowStatus(irow,PrePostsolveMatrix::basic);
423    rowduals[irow] = 0.0;       // ???
424    acts[irow] = 0.0;
425
426    // hinrow is not maintained during postsolve
427    //hinrow[irow] = 0;
428
429    rdone[irow] = DROP_ROW;
430  }
431
432  prob->nrows_ = prob->nrows_+nactions;
433}
434
Note: See TracBrowser for help on using the repository browser.