source: trunk/ClpPresolve.cpp @ 225

Last change on this file since 225 was 225, checked in by forrest, 16 years ago

This should break everything

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 36.2 KB
Line 
1// Copyright (C) 2002, International Business Machines
2// Corporation and others.  All Rights Reserved.
3
4//#define       CHECK_CONSISTENCY       1
5
6#include <stdio.h>
7
8#include <assert.h>
9#include <iostream>
10
11#include "CoinHelperFunctions.hpp"
12
13#include "CoinPackedMatrix.hpp"
14#include "ClpSimplex.hpp"
15
16#include "ClpPresolve.hpp"
17#include "CoinPresolveMatrix.hpp"
18
19#include "CoinPresolveEmpty.hpp"
20#include "CoinPresolveFixed.hpp"
21#include "CoinPresolvePsdebug.hpp"
22#include "CoinPresolveSingleton.hpp"
23#include "CoinPresolveDoubleton.hpp"
24#include "CoinPresolveTripleton.hpp"
25#include "CoinPresolveZeros.hpp"
26#include "CoinPresolveSubst.hpp"
27#include "CoinPresolveForcing.hpp"
28#include "CoinPresolveDual.hpp"
29#include "CoinPresolveTighten.hpp"
30#include "CoinPresolveUseless.hpp"
31#include "CoinPresolveDupcol.hpp"
32#include "CoinPresolveImpliedFree.hpp"
33#include "CoinPresolveIsolated.hpp"
34#include "CoinMessage.hpp"
35
36
37
38ClpPresolve::ClpPresolve() :
39  originalModel_(NULL),
40  presolvedModel_(NULL),
41  nonLinearValue_(0.0),
42  originalColumn_(NULL),
43  originalRow_(NULL),
44  paction_(0),
45  ncols_(0),
46  nrows_(0),
47  nelems_(0),
48  numberPasses_(5)
49{
50}
51
52ClpPresolve::~ClpPresolve()
53{
54  gutsOfDestroy();
55}
56// Gets rid of presolve actions (e.g.when infeasible)
57void 
58ClpPresolve::gutsOfDestroy()
59{
60 const CoinPresolveAction *paction = paction_;
61  while (paction) {
62    const CoinPresolveAction *next = paction->next;
63    delete paction;
64    paction = next;
65  }
66  delete [] originalColumn_;
67  delete [] originalRow_;
68  paction_=NULL;
69  originalColumn_=NULL;
70  originalRow_=NULL;
71}
72
73/* This version of presolve returns a pointer to a new presolved
74   model.  NULL if infeasible
75*/
76ClpSimplex * 
77ClpPresolve::presolvedModel(ClpSimplex & si,
78                         double feasibilityTolerance,
79                         bool keepIntegers,
80                         int numberPasses,
81                         bool dropNames)
82{
83  ncols_ = si.getNumCols();
84  nrows_ = si.getNumRows();
85  nelems_ = si.getNumElements();
86  numberPasses_ = numberPasses;
87
88  double maxmin = si.getObjSense();
89  originalModel_ = &si;
90  delete [] originalColumn_;
91  originalColumn_ = new int[ncols_];
92  delete [] originalRow_;
93  originalRow_ = new int[nrows_];
94
95  // Check matrix
96  if (!si.clpMatrix()->allElementsInRange(&si,si.getSmallElementValue(),
97                                          1.0e20))
98    return NULL;
99
100  // result is 0 - okay, 1 infeasible, -1 go round again
101  int result = -1;
102 
103  // User may have deleted - its their responsibility
104  presolvedModel_=NULL;
105  // Messages
106  CoinMessages messages = CoinMessage(si.messages().language());
107  while (result==-1) {
108
109    // make new copy
110    delete presolvedModel_;
111    presolvedModel_ = new ClpSimplex(si);
112    presolvedModel_->dropNames();
113
114    // drop integer information if wanted
115    if (!keepIntegers)
116      presolvedModel_->deleteIntegerInformation();
117
118   
119    CoinPresolveMatrix prob(ncols_,
120                        maxmin,
121                        presolvedModel_,
122                        nrows_, nelems_,true,nonLinearValue_);
123    // make sure row solution correct
124    {
125      double *colels    = prob.colels_;
126      int *hrow         = prob.hrow_;
127      CoinBigIndex *mcstrt              = prob.mcstrt_;
128      int *hincol               = prob.hincol_;
129      int ncols         = prob.ncols_;
130     
131     
132      double * csol = prob.sol_;
133      double * acts = prob.acts_;
134      int nrows = prob.nrows_;
135
136      int colx;
137
138      memset(acts,0,nrows*sizeof(double));
139     
140      for (colx = 0; colx < ncols; ++colx) {
141        double solutionValue = csol[colx];
142        for (int i=mcstrt[colx]; i<mcstrt[colx]+hincol[colx]; ++i) {
143          int row = hrow[i];
144          double coeff = colels[i];
145          acts[row] += solutionValue*coeff;
146        }
147      }
148    }
149    prob.handler_ = presolvedModel_->messageHandler();
150    prob.messages_ = CoinMessage(presolvedModel_->messages().language());
151
152    // move across feasibility tolerance
153    prob.feasibilityTolerance_ = feasibilityTolerance;
154
155    // Do presolve
156
157    paction_ = presolve(&prob);
158
159    result =0; 
160
161    if (prob.status_==0&&paction_) {
162      // Looks feasible but double check to see if anything slipped through
163      int n             = prob.ncols_;
164      double * lo = prob.clo_;
165      double * up = prob.cup_;
166      int i;
167     
168      for (i=0;i<n;i++) {
169        if (up[i]<lo[i]) {
170          if (up[i]<lo[i]-1.0e-8) {
171            // infeasible
172            prob.status_=1;
173          } else {
174            up[i]=lo[i];
175          }
176        }
177      }
178     
179      n = prob.nrows_;
180      lo = prob.rlo_;
181      up = prob.rup_;
182
183      for (i=0;i<n;i++) {
184        if (up[i]<lo[i]) {
185          if (up[i]<lo[i]-1.0e-8) {
186            // infeasible
187            prob.status_=1;
188          } else {
189            up[i]=lo[i];
190          }
191        }
192      }
193    }
194    if (prob.status_==0&&paction_) {
195      // feasible
196   
197      prob.update_model(presolvedModel_, nrows_, ncols_, nelems_);
198      // copy status and solution
199      memcpy(presolvedModel_->primalColumnSolution(),
200             prob.sol_,prob.ncols_*sizeof(double));
201      memcpy(presolvedModel_->primalRowSolution(),
202             prob.acts_,prob.nrows_*sizeof(double));
203      memcpy(presolvedModel_->statusArray(),
204             prob.colstat_,prob.ncols_*sizeof(unsigned char));
205      memcpy(presolvedModel_->statusArray()+prob.ncols_,
206             prob.rowstat_,prob.nrows_*sizeof(unsigned char));
207      delete [] prob.sol_;
208      delete [] prob.acts_;
209      delete [] prob.colstat_;
210      prob.sol_=NULL;
211      prob.acts_=NULL;
212      prob.colstat_=NULL;
213     
214      int ncolsNow = presolvedModel_->getNumCols();
215      memcpy(originalColumn_,prob.originalColumn_,ncolsNow*sizeof(int));
216      delete [] prob.originalColumn_;
217      prob.originalColumn_=NULL;
218      int nrowsNow = presolvedModel_->getNumRows();
219      memcpy(originalRow_,prob.originalRow_,nrowsNow*sizeof(int));
220      delete [] prob.originalRow_;
221      prob.originalRow_=NULL;
222      if (!dropNames&&si.lengthNames()) {
223        // Redo names
224        int iRow;
225        std::vector<std::string> rowNames;
226        rowNames.reserve(nrowsNow);
227        for (iRow=0;iRow<nrowsNow;iRow++) {
228          int kRow = originalRow_[iRow];
229          rowNames.push_back(si.rowName(kRow));
230        }
231     
232        int iColumn;
233        std::vector<std::string> columnNames;
234        columnNames.reserve(ncolsNow);
235        for (iColumn=0;iColumn<ncolsNow;iColumn++) {
236          int kColumn = originalColumn_[iColumn];
237          columnNames.push_back(si.columnName(kColumn));
238        }
239        presolvedModel_->copyNames(rowNames,columnNames);
240      }
241      // now clean up integer variables.  This can modify original
242      int i;
243      const char * information = presolvedModel_->integerInformation();
244      if (information) {
245        int numberChanges=0;
246        double * lower0 = originalModel_->columnLower();
247        double * upper0 = originalModel_->columnUpper();
248        double * lower = presolvedModel_->columnLower();
249        double * upper = presolvedModel_->columnUpper();
250        for (i=0;i<ncolsNow;i++) {
251          if (!information[i])
252            continue;
253          int iOriginal = originalColumn_[i];
254          double lowerValue0 = lower0[iOriginal];
255          double upperValue0 = upper0[iOriginal];
256          double lowerValue = ceil(lower[i]-1.0e-5);
257          double upperValue = floor(upper[i]+1.0e-5);
258          lower[i]=lowerValue;
259          upper[i]=upperValue;
260          if (lowerValue>upperValue) {
261            numberChanges++;
262            presolvedModel_->messageHandler()->message(COIN_PRESOLVE_COLINFEAS,
263                                                       messages)
264                                                         <<iOriginal
265                                                         <<lowerValue
266                                                         <<upperValue
267                                                         <<CoinMessageEol;
268            result=1;
269          } else {
270            if (lowerValue>lowerValue0+1.0e-8) {
271              lower0[iOriginal] = lowerValue;
272              numberChanges++;
273            }
274            if (upperValue<upperValue0-1.0e-8) {
275              upper0[iOriginal] = upperValue;
276              numberChanges++;
277            }
278          }       
279        }
280        if (numberChanges) {
281          presolvedModel_->messageHandler()->message(COIN_PRESOLVE_INTEGERMODS,
282                                                     messages)
283                                                       <<numberChanges
284                                                       <<CoinMessageEol;
285          if (!result) {
286            result = -1; // round again
287          }
288        }
289      }
290    } else {
291      // infeasible
292      delete [] prob.sol_;
293      delete [] prob.acts_;
294      delete [] prob.colstat_;
295      result=1;
296    }
297  }
298  if (!result) {
299    int nrowsAfter = presolvedModel_->getNumRows();
300    int ncolsAfter = presolvedModel_->getNumCols();
301    CoinBigIndex nelsAfter = presolvedModel_->getNumElements();
302    presolvedModel_->messageHandler()->message(COIN_PRESOLVE_STATS,
303                                               messages)
304                                                 <<nrowsAfter<< -(nrows_ - nrowsAfter)
305                                                 << ncolsAfter<< -(ncols_ - ncolsAfter)
306                                                 <<nelsAfter<< -(nelems_ - nelsAfter)
307                                                 <<CoinMessageEol;
308  } else {
309    gutsOfDestroy();
310    delete presolvedModel_;
311    presolvedModel_=NULL;
312  }
313  return presolvedModel_;
314}
315
316// Return pointer to presolved model
317ClpSimplex * 
318ClpPresolve::model() const
319{
320  return presolvedModel_;
321}
322// Return pointer to original model
323ClpSimplex * 
324ClpPresolve::originalModel() const
325{
326  return originalModel_;
327}
328void 
329ClpPresolve::postsolve(bool updateStatus)
330{
331  // Messages
332  CoinMessages messages = CoinMessage(presolvedModel_->messages().language());
333  if (!presolvedModel_->isProvenOptimal()) {
334    presolvedModel_->messageHandler()->message(COIN_PRESOLVE_NONOPTIMAL,
335                                             messages)
336                                               <<CoinMessageEol;
337  }
338
339  // this is the size of the original problem
340  const int ncols0  = ncols_;
341  const int nrows0  = nrows_;
342  const CoinBigIndex nelems0 = nelems_;
343
344  // reality check
345  assert(ncols0==originalModel_->getNumCols());
346  assert(nrows0==originalModel_->getNumRows());
347
348  // this is the reduced problem
349  int ncols = presolvedModel_->getNumCols();
350  int nrows = presolvedModel_->getNumRows();
351
352  double *acts = originalModel_->primalRowSolution();
353  double *sol  = originalModel_->primalColumnSolution();
354
355  unsigned char * rowstat=NULL;
356  unsigned char * colstat = NULL;
357  if (updateStatus) {
358    unsigned char *status = originalModel_->statusArray();
359    rowstat = status + ncols0;
360    colstat = status;
361    memcpy(colstat, presolvedModel_->statusArray(), ncols);
362    memcpy(rowstat, presolvedModel_->statusArray()+ncols, nrows);
363  }
364
365
366  CoinPostsolveMatrix prob(presolvedModel_,
367                       ncols0,
368                       nrows0,
369                       nelems0,
370                       presolvedModel_->getObjSense(),
371                       // end prepost
372                       
373                       sol, acts,
374                       colstat, rowstat);
375   
376  postsolve(prob);
377
378}
379
380// return pointer to original columns
381const int * 
382ClpPresolve::originalColumns() const
383{
384  return originalColumn_;
385}
386// return pointer to original rows
387const int * 
388ClpPresolve::originalRows() const
389{
390  return originalRow_;
391}
392// Set pointer to original model
393void 
394ClpPresolve::setOriginalModel(ClpSimplex * model)
395{
396  originalModel_=model;
397}
398#if 0
399// A lazy way to restrict which transformations are applied
400// during debugging.
401static int ATOI(const char *name)
402{
403 return true;
404#if     DEBUG_PRESOLVE || PRESOLVE_SUMMARY
405  if (getenv(name)) {
406    int val = atoi(getenv(name));
407    printf("%s = %d\n", name, val);
408    return (val);
409  } else {
410    if (strcmp(name,"off"))
411      return (true);
412    else
413      return (false);
414  }
415#else
416  return (true);
417#endif
418}
419#endif
420//#define DEBUG_PRESOLVE 1
421#if DEBUG_PRESOLVE
422void check_sol(CoinPresolveMatrix *prob,double tol)
423{
424  double *colels        = prob->colels_;
425  int *hrow             = prob->hrow_;
426  int *mcstrt           = prob->mcstrt_;
427  int *hincol           = prob->hincol_;
428  int *hinrow           = prob->hinrow_;
429  int ncols             = prob->ncols_;
430
431
432  double * csol = prob->sol_;
433  double * acts = prob->acts_;
434  double * clo = prob->clo_;
435  double * cup = prob->cup_;
436  int nrows = prob->nrows_;
437  double * rlo = prob->rlo_;
438  double * rup = prob->rup_;
439
440  int colx;
441
442  double * rsol = new double[nrows];
443  memset(rsol,0,nrows*sizeof(double));
444
445  for (colx = 0; colx < ncols; ++colx) {
446    if (1) {
447      CoinBigIndex k = mcstrt[colx];
448      int nx = hincol[colx];
449      double solutionValue = csol[colx];
450      for (int i=0; i<nx; ++i) {
451        int row = hrow[k];
452        double coeff = colels[k];
453        k++;
454        rsol[row] += solutionValue*coeff;
455      }
456      if (csol[colx]<clo[colx]-tol) {
457        printf("low CSOL:  %d  - %g %g %g\n",
458                   colx, clo[colx], csol[colx], cup[colx]);
459      } else if (csol[colx]>cup[colx]+tol) {
460        printf("high CSOL:  %d  - %g %g %g\n",
461                   colx, clo[colx], csol[colx], cup[colx]);
462      } 
463    }
464  }
465  int rowx;
466  for (rowx = 0; rowx < nrows; ++rowx) {
467    if (hinrow[rowx]) {
468      if (fabs(rsol[rowx]-acts[rowx])>tol)
469        printf("inacc RSOL:  %d - %g %g (acts_ %g) %g\n",
470                   rowx,  rlo[rowx], rsol[rowx], acts[rowx], rup[rowx]);
471      if (rsol[rowx]<rlo[rowx]-tol) {
472        printf("low RSOL:  %d - %g %g %g\n",
473                   rowx,  rlo[rowx], rsol[rowx], rup[rowx]);
474      } else if (rsol[rowx]>rup[rowx]+tol ) {
475        printf("high RSOL:  %d - %g %g %g\n",
476                   rowx,  rlo[rowx], rsol[rowx], rup[rowx]);
477      } 
478    }
479  }
480  delete [] rsol;
481}
482#endif
483// This is the presolve loop.
484// It is a separate virtual function so that it can be easily
485// customized by subclassing CoinPresolve.
486const CoinPresolveAction *ClpPresolve::presolve(CoinPresolveMatrix *prob)
487{
488  paction_ = 0;
489
490  prob->status_=0; // say feasible
491
492  paction_ = make_fixed(prob, paction_);
493  // if integers then switch off dual stuff
494  // later just do individually
495  bool doDualStuff = (presolvedModel_->integerInformation()==NULL);
496
497#if     CHECK_CONSISTENCY
498  presolve_links_ok(prob->rlink_, prob->mrstrt_, prob->hinrow_, prob->nrows_);
499#endif
500
501  if (!prob->status_) {
502#if 0
503    const bool slackd = ATOI("SLACKD")!=0;
504    //const bool forcing = ATOI("FORCING")!=0;
505    const bool doubleton = ATOI("DOUBLETON")!=0;
506    const bool forcing = ATOI("off")!=0;
507    const bool ifree = ATOI("off")!=0;
508    const bool zerocost = ATOI("off")!=0;
509    const bool dupcol = ATOI("off")!=0;
510    const bool duprow = ATOI("off")!=0;
511    const bool dual = ATOI("off")!=0;
512#else
513    // normal
514#if 0
515    const bool slackd = true;
516    const bool doubleton = false;
517    const bool tripleton = true;
518    const bool forcing = false;
519    const bool ifree = true;
520    const bool zerocost = true;
521    const bool dupcol = false;
522    const bool duprow = false;
523    const bool dual = doDualStuff;
524#else
525    const bool slackd = true;
526    const bool doubleton = true;
527    const bool tripleton = true;
528    const bool forcing = true;
529    const bool ifree = true;
530    const bool zerocost = true;
531    const bool dupcol = true;
532    const bool duprow = true;
533    const bool dual = doDualStuff;
534#endif
535#endif
536   
537    // some things are expensive so just do once (normally)
538
539    int i;
540    // say look at all
541    if (!prob->anyProhibited()) {
542      for (i=0;i<nrows_;i++) 
543        prob->rowsToDo_[i]=i;
544      prob->numberRowsToDo_=nrows_;
545      for (i=0;i<ncols_;i++) 
546        prob->colsToDo_[i]=i;
547      prob->numberColsToDo_=ncols_;
548    } else {
549      // some stuff must be left alone
550      prob->numberRowsToDo_=0;
551      for (i=0;i<nrows_;i++) 
552        if (!prob->rowProhibited(i))
553            prob->rowsToDo_[prob->numberRowsToDo_++]=i;
554      prob->numberColsToDo_=0;
555      for (i=0;i<ncols_;i++) 
556        if (!prob->colProhibited(i))
557            prob->colsToDo_[prob->numberColsToDo_++]=i;
558    }
559
560
561    int iLoop;
562#if     DEBUG_PRESOLVE
563    check_sol(prob,1.0e0);
564#endif
565
566    // Check number rows dropped
567    int lastDropped=0;
568    for (iLoop=0;iLoop<numberPasses_;iLoop++) {
569#ifdef PRESOLVE_SUMMARY
570      printf("Starting major pass %d\n",iLoop+1);
571#endif
572      const CoinPresolveAction * const paction0 = paction_;
573      // look for substitutions with no fill
574      int fill_level=2;
575      //fill_level=10;
576      //printf("** fill_level == 10 !\n");
577      int whichPass=0;
578      while (1) {
579        whichPass++;
580        const CoinPresolveAction * const paction1 = paction_;
581
582        if (slackd) {
583          bool notFinished = true;
584          while (notFinished) 
585            paction_ = slack_doubleton_action::presolve(prob, paction_,
586                                                        notFinished);
587          if (prob->status_)
588            break;
589        }
590
591        if (doubleton) {
592          paction_ = doubleton_action::presolve(prob, paction_);
593          if (prob->status_)
594            break;
595        }
596
597        if (tripleton) {
598          paction_ = tripleton_action::presolve(prob, paction_);
599          if (prob->status_)
600            break;
601        }
602
603        if (zerocost) {
604          paction_ = do_tighten_action::presolve(prob, paction_);
605          if (prob->status_)
606            break;
607        }
608
609        if (forcing) {
610          paction_ = forcing_constraint_action::presolve(prob, paction_);
611          if (prob->status_)
612            break;
613        }
614
615        if (ifree) {
616          paction_ = implied_free_action::presolve(prob, paction_,fill_level);
617          if (prob->status_)
618            break;
619        }
620
621#if     DEBUG_PRESOLVE
622        check_sol(prob,1.0e0);
623#endif
624
625#if     CHECK_CONSISTENCY
626        presolve_links_ok(prob->rlink_, prob->mrstrt_, prob->hinrow_, 
627                          prob->nrows_);
628#endif
629
630#if     DEBUG_PRESOLVE
631        presolve_no_zeros(prob->mcstrt_, prob->colels_, prob->hincol_, 
632                          prob->ncols_);
633#endif
634#if     CHECK_CONSISTENCY
635        prob->consistent();
636#endif
637
638         
639        // set up for next pass
640        // later do faster if many changes i.e. memset and memcpy
641        prob->numberRowsToDo_ = prob->numberNextRowsToDo_;
642        int kcheck;
643        bool found=false;
644        kcheck=-1;
645        for (i=0;i<prob->numberNextRowsToDo_;i++) {
646          int index = prob->nextRowsToDo_[i];
647          prob->unsetRowChanged(index);
648          prob->rowsToDo_[i] = index;
649          if (index==kcheck) {
650            printf("row %d on list after pass %d\n",kcheck,
651                   whichPass);
652            found=true;
653          }
654        }
655        if (!found&&kcheck>=0)
656          prob->rowsToDo_[prob->numberRowsToDo_++]=kcheck;
657        prob->numberNextRowsToDo_=0;
658        prob->numberColsToDo_ = prob->numberNextColsToDo_;
659        kcheck=-1;
660        found=false;
661        for (i=0;i<prob->numberNextColsToDo_;i++) {
662          int index = prob->nextColsToDo_[i];
663          prob->unsetColChanged(index);
664          prob->colsToDo_[i] = index;
665          if (index==kcheck) {
666            printf("col %d on list after pass %d\n",kcheck,
667                   whichPass);
668            found=true;
669          }
670        }
671        if (!found&&kcheck>=0)
672          prob->colsToDo_[prob->numberColsToDo_++]=kcheck;
673        prob->numberNextColsToDo_=0;
674        if (paction_ == paction1&&fill_level>0)
675          break;
676      }
677      // say look at all
678      int i;
679      if (!prob->anyProhibited()) {
680        for (i=0;i<nrows_;i++) 
681          prob->rowsToDo_[i]=i;
682        prob->numberRowsToDo_=nrows_;
683        for (i=0;i<ncols_;i++) 
684          prob->colsToDo_[i]=i;
685        prob->numberColsToDo_=ncols_;
686      } else {
687        // some stuff must be left alone
688        prob->numberRowsToDo_=0;
689        for (i=0;i<nrows_;i++) 
690          if (!prob->rowProhibited(i))
691            prob->rowsToDo_[prob->numberRowsToDo_++]=i;
692        prob->numberColsToDo_=0;
693        for (i=0;i<ncols_;i++) 
694          if (!prob->colProhibited(i))
695            prob->colsToDo_[prob->numberColsToDo_++]=i;
696      }
697      // now expensive things
698      // this caused world.mps to run into numerical difficulties
699#ifdef PRESOLVE_SUMMARY
700      printf("Starting expensive\n");
701#endif
702
703      if (dual) {
704        int itry;
705        for (itry=0;itry<5;itry++) {
706          const CoinPresolveAction * const paction2 = paction_;
707          paction_ = remove_dual_action::presolve(prob, paction_);
708          if (prob->status_)
709            break;
710          if (ifree) {
711            int fill_level=0; // switches off substitution
712            paction_ = implied_free_action::presolve(prob, paction_,fill_level);
713            if (prob->status_)
714              break;
715          }
716          if (paction_ == paction2)
717            break;
718        }
719      }
720#if     DEBUG_PRESOLVE
721      check_sol(prob,1.0e0);
722#endif
723      if (dupcol) {
724        paction_ = dupcol_action::presolve(prob, paction_);
725        if (prob->status_)
726          break;
727      }
728#if     DEBUG_PRESOLVE
729        check_sol(prob,1.0e0);
730#endif
731     
732      if (duprow) {
733        paction_ = duprow_action::presolve(prob, paction_);
734        if (prob->status_)
735          break;
736      }
737#if     DEBUG_PRESOLVE
738      check_sol(prob,1.0e0);
739#endif
740      {
741        int * hinrow = prob->hinrow_;
742        int numberDropped=0;
743        for (i=0;i<nrows_;i++) 
744          if (!hinrow[i])
745            numberDropped++;
746        //printf("%d rows dropped after pass %d\n",numberDropped,
747        //     iLoop+1);
748        if (numberDropped==lastDropped)
749          break;
750        else
751          lastDropped = numberDropped;
752      }
753      if (paction_ == paction0)
754        break;
755         
756    }
757  }
758  if (!prob->status_) {
759    paction_ = drop_zero_coefficients(prob, paction_);
760#if     DEBUG_PRESOLVE
761        check_sol(prob,1.0e0);
762#endif
763
764    paction_ = drop_empty_cols_action::presolve(prob, paction_);
765    paction_ = drop_empty_rows_action::presolve(prob, paction_);
766#if     DEBUG_PRESOLVE
767        check_sol(prob,1.0e0);
768#endif
769  }
770 
771  // Messages
772  CoinMessages messages = CoinMessage(prob->messages().language());
773  if (prob->status_) {
774    if (prob->status_==1)
775          prob->messageHandler()->message(COIN_PRESOLVE_INFEAS,
776                                             messages)
777                                               <<prob->feasibilityTolerance_
778                                               <<CoinMessageEol;
779    else if (prob->status_==2)
780          prob->messageHandler()->message(COIN_PRESOLVE_UNBOUND,
781                                             messages)
782                                               <<CoinMessageEol;
783    else
784          prob->messageHandler()->message(COIN_PRESOLVE_INFEASUNBOUND,
785                                             messages)
786                                               <<CoinMessageEol;
787    // get rid of data
788    gutsOfDestroy();
789  }
790  return (paction_);
791}
792
793void check_djs(CoinPostsolveMatrix *prob);
794
795
796// We could have implemented this by having each postsolve routine
797// directly call the next one, but this may make it easier to add debugging checks.
798void ClpPresolve::postsolve(CoinPostsolveMatrix &prob)
799{
800  const CoinPresolveAction *paction = paction_;
801
802  if (prob.colstat_)
803    prob.check_nbasic();
804 
805#if     DEBUG_PRESOLVE
806  check_djs(&prob);
807#endif
808 
809 
810  while (paction) {
811#if     DEBUG_PRESOLVE
812    printf("POSTSOLVING %s\n", paction->name());
813#endif
814
815    paction->postsolve(&prob);
816   
817#if     DEBUG_PRESOLVE
818    if (prob.colstat_)
819      prob.check_nbasic();
820#endif
821    paction = paction->next;
822#if     DEBUG_PRESOLVE
823    check_djs(&prob);
824#endif
825  }   
826 
827#if     0 && DEBUG_PRESOLVE
828  for (i=0; i<ncols0; i++) {
829    if (!cdone[i]) {
830      printf("!cdone[%d]\n", i);
831      abort();
832    }
833  }
834 
835  for (i=0; i<nrows0; i++) {
836    if (!rdone[i]) {
837      printf("!rdone[%d]\n", i);
838      abort();
839    }
840  }
841 
842 
843  for (i=0; i<ncols0; i++) {
844    if (sol[i] < -1e10 || sol[i] > 1e10)
845      printf("!!!%d %g\n", i, sol[i]);
846   
847  }
848 
849 
850#endif
851 
852#if     0 && DEBUG_PRESOLVE
853  // debug check:  make sure we ended up with same original matrix
854  {
855    int identical = 1;
856   
857    for (int i=0; i<ncols0; i++) {
858      PRESOLVEASSERT(hincol[i] == &prob->mcstrt0[i+1] - &prob->mcstrt0[i]);
859      CoinBigIndex kcs0 = &prob->mcstrt0[i];
860      CoinBigIndex kcs = mcstrt[i];
861      int n = hincol[i];
862      for (int k=0; k<n; k++) {
863        CoinBigIndex k1 = presolve_find_row1(&prob->hrow0[kcs0+k], kcs, kcs+n, hrow);
864
865        if (k1 == kcs+n) {
866          printf("ROW %d NOT IN COL %d\n", &prob->hrow0[kcs0+k], i);
867          abort();
868        }
869
870        if (colels[k1] != &prob->dels0[kcs0+k])
871          printf("BAD COLEL[%d %d %d]:  %g\n",
872                 k1, i, &prob->hrow0[kcs0+k], colels[k1] - &prob->dels0[kcs0+k]);
873
874        if (kcs0+k != k1)
875          identical=0;
876      }
877    }
878    printf("identical? %d\n", identical);
879  }
880#endif
881  // put back duals
882  memcpy(originalModel_->dualRowSolution(),prob.rowduals_,
883         nrows_*sizeof(double));
884  double maxmin = originalModel_->getObjSense();
885  if (maxmin<0.0) {
886    // swap signs
887    int i;
888    double * pi = originalModel_->dualRowSolution();
889    for (i=0;i<nrows_;i++)
890      pi[i] = -pi[i];
891  }
892  // Now check solution
893  memcpy(originalModel_->dualColumnSolution(),
894         originalModel_->objective(),ncols_*sizeof(double));
895  originalModel_->transposeTimes(-1.0,
896                                 originalModel_->dualRowSolution(),
897                                 originalModel_->dualColumnSolution());
898  memset(originalModel_->primalRowSolution(),0,nrows_*sizeof(double));
899  originalModel_->times(1.0,originalModel_->primalColumnSolution(),
900                        originalModel_->primalRowSolution());
901  originalModel_->checkSolution();
902  // Messages
903  CoinMessages messages = CoinMessage(presolvedModel_->messages().language());
904  presolvedModel_->messageHandler()->message(COIN_PRESOLVE_POSTSOLVE,
905                                            messages)
906                                              <<originalModel_->objectiveValue()
907                                              <<originalModel_->sumDualInfeasibilities()
908                                              <<originalModel_->numberDualInfeasibilities()
909                                              <<originalModel_->sumPrimalInfeasibilities()
910                                              <<originalModel_->numberPrimalInfeasibilities()
911                                               <<CoinMessageEol;
912 
913  //originalModel_->objectiveValue_=objectiveValue_;
914  originalModel_->setNumberIterations(presolvedModel_->numberIterations());
915  if (!presolvedModel_->status()) {
916    if (!originalModel_->numberDualInfeasibilities()&&
917        !originalModel_->numberPrimalInfeasibilities()) {
918      originalModel_->setProblemStatus( 0);
919    } else {
920      originalModel_->setProblemStatus( -1);
921      presolvedModel_->messageHandler()->message(COIN_PRESOLVE_NEEDS_CLEANING,
922                                            messages)
923                                              <<CoinMessageEol;
924    }
925  } else {
926    originalModel_->setProblemStatus( presolvedModel_->status());
927  }
928}
929
930
931
932
933
934
935
936
937#if     DEBUG_PRESOLVE
938void check_djs(CoinPostsolveMatrix *prob)
939{
940  //return;
941  double *colels        = prob->colels_;
942  int *hrow             = prob->hrow_;
943  int *mcstrt           = prob->mcstrt_;
944  int *hincol           = prob->hincol_;
945  int *link             = prob->link_;
946  int ncols             = prob->ncols_;
947
948  double *dcost = prob->cost_;
949
950  double *rcosts        = prob->rcosts_;
951
952  double *rowduals = prob->rowduals_;
953
954  const double maxmin   = prob->maxmin_;
955
956  char *cdone   = prob->cdone_;
957
958  double * csol = prob->sol_;
959  double * clo = prob->clo_;
960  double * cup = prob->cup_;
961  int nrows = prob->nrows_;
962  double * rlo = prob->rlo_;
963  double * rup = prob->rup_;
964  char *rdone   = prob->rdone_;
965
966  int colx;
967
968  double * rsol = new double[nrows];
969  memset(rsol,0,nrows*sizeof(double));
970
971  for (colx = 0; colx < ncols; ++colx) {
972    if (cdone[colx]) {
973      CoinBigIndex k = mcstrt[colx];
974      int nx = hincol[colx];
975      double dj = maxmin * dcost[colx];
976      double solutionValue = csol[colx];
977      for (int i=0; i<nx; ++i) {
978        int row = hrow[k];
979        double coeff = colels[k];
980        k = link[k];
981        dj -= rowduals[row] * coeff;
982        rsol[row] += solutionValue*coeff;
983      }
984      if (! (fabs(rcosts[colx] - dj) < 1.0e-4))
985        printf("BAD DJ:  %d %g %g\n",
986               colx, rcosts[colx], dj);
987      if (cup[colx]-clo[colx]>1.0e-6) {
988        if (csol[colx]<clo[colx]+1.0e-6) {
989          if (dj <-1.0e-6)
990            printf("neg DJ:  %d %g  - %g %g %g\n",
991                   colx, dj, clo[colx], csol[colx], cup[colx]);
992        } else if (csol[colx]>cup[colx]-1.0e-6) {
993          if (dj > 1.0e-6)
994            printf("pos DJ:  %d %g  - %g %g %g\n",
995                   colx, dj, clo[colx], csol[colx], cup[colx]);
996        } else {
997          if (fabs(dj) >1.0e-6)
998            printf("nonzero DJ:  %d %g  - %g %g %g\n",
999                   colx, dj, clo[colx], csol[colx], cup[colx]);
1000        }
1001      } 
1002    }
1003  }
1004  int rowx;
1005  for (rowx = 0; rowx < nrows; ++rowx) {
1006    if (rdone[rowx]) {
1007      if (rup[rowx]-rlo[rowx]>1.0e-6) {
1008        double dj = rowduals[rowx];
1009        if (rsol[rowx]<rlo[rowx]+1.0e-6) {
1010          if (dj <-1.0e-5)
1011            printf("neg rDJ:  %d %g  - %g %g %g\n",
1012                   rowx, dj, rlo[rowx], rsol[rowx], rup[rowx]);
1013        } else if (rsol[rowx]>rup[rowx]-1.0e-6) {
1014          if (dj > 1.0e-5)
1015            printf("pos rDJ:  %d %g  - %g %g %g\n",
1016                   rowx, dj, rlo[rowx], rsol[rowx], rup[rowx]);
1017        } else {
1018          if (fabs(dj) >1.0e-5)
1019            printf("nonzero rDJ:  %d %g  - %g %g %g\n",
1020                   rowx, dj, rlo[rowx], rsol[rowx], rup[rowx]);
1021        }
1022      } 
1023    }
1024  }
1025  delete [] rsol;
1026}
1027#endif
1028
1029static inline double getTolerance(const ClpSimplex  *si, ClpDblParam key)
1030{
1031  double tol;
1032  if (! si->getDblParam(key, tol)) {
1033    CoinPresolveAction::throwCoinError("getDblParam failed",
1034                                      "CoinPrePostsolveMatrix::CoinPrePostsolveMatrix");
1035  }
1036  return (tol);
1037}
1038
1039
1040// Assumptions:
1041// 1. nrows>=m.getNumRows()
1042// 2. ncols>=m.getNumCols()
1043//
1044// In presolve, these values are equal.
1045// In postsolve, they may be inequal, since the reduced problem
1046// may be smaller, but we need room for the large problem.
1047// ncols may be larger than si.getNumCols() in postsolve,
1048// this at that point si will be the reduced problem,
1049// but we need to reserve enough space for the original problem.
1050CoinPrePostsolveMatrix::CoinPrePostsolveMatrix(const ClpSimplex * si,
1051                                             int ncols_in,
1052                                             int nrows_in,
1053                                             CoinBigIndex nelems_in) :
1054  ncols_(si->getNumCols()),
1055  ncols0_(ncols_in),
1056  nelems_(si->getNumElements()),
1057
1058  mcstrt_(new CoinBigIndex[ncols_in+1]),
1059  hincol_(new int[ncols_in+1]),
1060  hrow_  (new int   [2*nelems_in]),
1061  colels_(new double[2*nelems_in]),
1062
1063  cost_(new double[ncols_in]),
1064  clo_(new double[ncols_in]),
1065  cup_(new double[ncols_in]),
1066  rlo_(new double[nrows_in]),
1067  rup_(new double[nrows_in]),
1068  originalColumn_(new int[ncols_in]),
1069  originalRow_(new int[nrows_in]),
1070
1071  ztolzb_(getTolerance(si, ClpPrimalTolerance)),
1072  ztoldj_(getTolerance(si, ClpDualTolerance)),
1073
1074  maxmin_(si->getObjSense())
1075
1076{
1077  si->getDblParam(ClpObjOffset,originalOffset_);
1078  int ncols = si->getNumCols();
1079  int nrows = si->getNumRows();
1080
1081  ClpDisjointCopyN(si->getColLower(), ncols, clo_);
1082  ClpDisjointCopyN(si->getColUpper(), ncols, cup_);
1083  ClpDisjointCopyN(si->getObjCoefficients(), ncols, cost_);
1084  ClpDisjointCopyN(si->getRowLower(), nrows,  rlo_);
1085  ClpDisjointCopyN(si->getRowUpper(), nrows,  rup_);
1086  int i;
1087  for (i=0;i<ncols_in;i++) 
1088    originalColumn_[i]=i;
1089  for (i=0;i<nrows_in;i++) 
1090    originalRow_[i]=i;
1091  sol_=NULL;
1092  rowduals_=NULL;
1093  acts_=NULL;
1094
1095  rcosts_=NULL;
1096  colstat_=NULL;
1097  rowstat_=NULL;
1098}
1099
1100// I am not familiar enough with CoinPackedMatrix to be confident
1101// that I will implement a row-ordered version of toColumnOrderedGapFree
1102// properly.
1103static bool isGapFree(const CoinPackedMatrix& matrix)
1104{
1105  const CoinBigIndex * start = matrix.getVectorStarts();
1106  const int * length = matrix.getVectorLengths();
1107  int i;
1108  for (i = matrix.getSizeVectorLengths() - 1; i >= 0; --i) {
1109    if (start[i+1] - start[i] != length[i])
1110      break;
1111  }
1112  return (! (i >= 0));
1113}
1114#if     DEBUG_PRESOLVE
1115static void matrix_bounds_ok(const double *lo, const double *up, int n)
1116{
1117  int i;
1118  for (i=0; i<n; i++) {
1119    PRESOLVEASSERT(lo[i] <= up[i]);
1120    PRESOLVEASSERT(lo[i] < PRESOLVE_INF);
1121    PRESOLVEASSERT(-PRESOLVE_INF < up[i]);
1122  }
1123}
1124#endif
1125CoinPresolveMatrix::CoinPresolveMatrix(int ncols0_in,
1126                                     double maxmin_,
1127                                     // end prepost members
1128
1129                                     ClpSimplex * si,
1130
1131                                     // rowrep
1132                                     int nrows_in,
1133                                     CoinBigIndex nelems_in,
1134                               bool doStatus,
1135                               double nonLinearValue) :
1136
1137  CoinPrePostsolveMatrix(si,
1138                        ncols0_in, nrows_in, nelems_in),
1139  clink_(new presolvehlink[ncols0_in+1]),
1140  rlink_(new presolvehlink[nrows_in+1]),
1141
1142  dobias_(0.0),
1143
1144  nrows_(si->getNumRows()),
1145
1146  // temporary init
1147  mrstrt_(new CoinBigIndex[nrows_in+1]),
1148  hinrow_(new int[nrows_in+1]),
1149  rowels_(new double[2*nelems_in]),
1150  hcol_(new int[2*nelems_in]),
1151  integerType_(new char[ncols0_in]),
1152  feasibilityTolerance_(0.0),
1153  status_(-1),
1154  rowsToDo_(new int [nrows_in]),
1155  numberRowsToDo_(0),
1156  nextRowsToDo_(new int[nrows_in]),
1157  numberNextRowsToDo_(0),
1158  colsToDo_(new int [ncols0_in]),
1159  numberColsToDo_(0),
1160  nextColsToDo_(new int[ncols0_in]),
1161  numberNextColsToDo_(0)
1162
1163{
1164  const int bufsize = 2*nelems_in;
1165
1166  // Set up change bits etc
1167  rowChanged_ = new unsigned char[nrows_];
1168  memset(rowChanged_,0,nrows_);
1169  colChanged_ = new unsigned char[ncols_];
1170  memset(colChanged_,0,ncols_);
1171  CoinPackedMatrix * m = si->matrix();
1172
1173  // The coefficient matrix is a big hunk of stuff.
1174  // Do the copy here to try to avoid running out of memory.
1175
1176  const CoinBigIndex * start = m->getVectorStarts();
1177  const int * length = m->getVectorLengths();
1178  const int * row = m->getIndices();
1179  const double * element = m->getElements();
1180  int icol,nel=0;
1181  mcstrt_[0]=0;
1182  for (icol=0;icol<ncols_;icol++) {
1183    int j;
1184    for (j=start[icol];j<start[icol]+length[icol];j++) {
1185      hrow_[nel]=row[j];
1186      colels_[nel++]=element[j];
1187    }
1188    mcstrt_[icol+1]=nel;
1189  }
1190  assert(mcstrt_[ncols_] == nelems_);
1191  ClpDisjointCopyN(m->getVectorLengths(),ncols_,  hincol_);
1192
1193  // same thing for row rep
1194  m = new CoinPackedMatrix();
1195  m->reverseOrderedCopyOf(*si->matrix());
1196  m->removeGaps();
1197
1198
1199  ClpDisjointCopyN(m->getVectorStarts(),  nrows_,  mrstrt_);
1200  mrstrt_[nrows_] = nelems_;
1201  ClpDisjointCopyN(m->getVectorLengths(), nrows_,  hinrow_);
1202  ClpDisjointCopyN(m->getIndices(),       nelems_, hcol_);
1203  ClpDisjointCopyN(m->getElements(),      nelems_, rowels_);
1204
1205  delete m;
1206  if (si->integerInformation()) {
1207    memcpy(integerType_,si->integerInformation(),ncols_*sizeof(char));
1208  } else {
1209    ClpFillN<char>(integerType_, ncols_, 0);
1210  }
1211
1212  // Set up prohibited bits if needed
1213  if (nonLinearValue) {
1214    anyProhibited_ = true;
1215    for (icol=0;icol<ncols_;icol++) {
1216      int j;
1217      bool nonLinearColumn = false;
1218      if (cost_[icol]==nonLinearValue)
1219        nonLinearColumn=true;
1220      for (j=mcstrt_[icol];j<mcstrt_[icol+1];j++) {
1221        if (colels_[j]==nonLinearValue) {
1222          nonLinearColumn=true;
1223          setRowProhibited(hrow_[j]);
1224        }
1225      }
1226      if (nonLinearColumn)
1227        setColProhibited(icol);
1228    }
1229  } else {
1230    anyProhibited_ = false;
1231  }
1232
1233  if (doStatus) {
1234    // allow for status and solution
1235    sol_ = new double[ncols_];
1236    memcpy(sol_,si->primalColumnSolution(),ncols_*sizeof(double));;
1237    acts_ = new double [nrows_];
1238    memcpy(acts_,si->primalRowSolution(),nrows_*sizeof(double));
1239    if (!si->statusArray())
1240      si->createStatus();
1241    colstat_ = new unsigned char [nrows_+ncols_];
1242    memcpy(colstat_,si->statusArray(),
1243           (nrows_+ncols_)*sizeof(unsigned char));
1244    rowstat_ = colstat_+ncols_;
1245  }
1246
1247  // the original model's fields are now unneeded - free them
1248 
1249  si->resize(0,0);
1250
1251#if     DEBUG_PRESOLVE
1252  matrix_bounds_ok(rlo_, rup_, nrows_);
1253  matrix_bounds_ok(clo_, cup_, ncols_);
1254#endif
1255
1256#if 0
1257  for (i=0; i<nrows; ++i)
1258    printf("NR: %6d\n", hinrow[i]);
1259  for (int i=0; i<ncols; ++i)
1260    printf("NC: %6d\n", hincol[i]);
1261#endif
1262
1263  presolve_make_memlists(mcstrt_, hincol_, clink_, ncols_);
1264  presolve_make_memlists(mrstrt_, hinrow_, rlink_, nrows_);
1265
1266  // this allows last col/row to expand up to bufsize-1 (22);
1267  // this must come after the calls to presolve_prefix
1268  mcstrt_[ncols_] = bufsize-1;
1269  mrstrt_[nrows_] = bufsize-1;
1270
1271#if     CHECK_CONSISTENCY
1272  consistent(false);
1273#endif
1274}
1275
1276void CoinPresolveMatrix::update_model(ClpSimplex * si,
1277                                     int nrows0,
1278                                     int ncols0,
1279                                     CoinBigIndex nelems0)
1280{
1281  si->loadProblem(ncols_, nrows_, mcstrt_, hrow_, colels_, hincol_,
1282                 clo_, cup_, cost_, rlo_, rup_);
1283
1284  delete [] si->integerInformation();
1285  int numberIntegers=0;
1286  for (int i=0; i<ncols_; i++) {
1287    if (integerType_[i])
1288      numberIntegers++;
1289  }
1290  if (numberIntegers) 
1291    si->copyInIntegerInformation(integerType_);
1292  else
1293    si->copyInIntegerInformation(NULL);
1294
1295#if     PRESOLVE_SUMMARY
1296  printf("NEW NCOL/NROW/NELS:  %d(-%d) %d(-%d) %d(-%d)\n",
1297         ncols_, ncols0-ncols_,
1298         nrows_, nrows0-nrows_,
1299         si->getNumElements(), nelems0-si->getNumElements());
1300#endif
1301  si->setDblParam(ClpObjOffset,originalOffset_-dobias_);
1302
1303}
1304
1305
1306
1307
1308
1309
1310
1311
1312
1313
1314
1315////////////////  POSTSOLVE
1316
1317CoinPostsolveMatrix::CoinPostsolveMatrix(ClpSimplex*  si,
1318                                       int ncols0_in,
1319                                       int nrows0_in,
1320                                       CoinBigIndex nelems0,
1321                                   
1322                                       double maxmin_,
1323                                       // end prepost members
1324
1325                                       double *sol_in,
1326                                       double *acts_in,
1327
1328                                       unsigned char *colstat_in,
1329                                       unsigned char *rowstat_in) :
1330  CoinPrePostsolveMatrix(si,
1331                        ncols0_in, nrows0_in, nelems0),
1332
1333  free_list_(0),
1334  // link, free_list, maxlink
1335  maxlink_(2*nelems0),
1336  link_(new int[/*maxlink*/ 2*nelems0]),
1337     
1338  cdone_(new char[ncols0_]),
1339  rdone_(new char[nrows0_in]),
1340
1341  nrows_(si->getNumRows()),
1342  nrows0_(nrows0_in)
1343{
1344
1345  sol_=sol_in;
1346  rowduals_=NULL;
1347  acts_=acts_in;
1348
1349  rcosts_=NULL;
1350  colstat_=colstat_in;
1351  rowstat_=rowstat_in;
1352
1353  // this is the *reduced* model, which is probably smaller
1354  int ncols1 = si->getNumCols();
1355  int nrows1 = si->getNumRows();
1356
1357  const CoinPackedMatrix * m = si->matrix();
1358
1359  if (! isGapFree(*m)) {
1360    CoinPresolveAction::throwCoinError("Matrix not gap free",
1361                                      "CoinPostsolveMatrix");
1362  }
1363
1364  const CoinBigIndex nelemsr = m->getNumElements();
1365
1366  ClpDisjointCopyN(m->getVectorStarts(), ncols1, mcstrt_);
1367  mcstrt_[ncols_] = nelems0;    // ??
1368  ClpDisjointCopyN(m->getVectorLengths(),ncols1,  hincol_);
1369  ClpDisjointCopyN(m->getIndices(),      nelemsr, hrow_);
1370  ClpDisjointCopyN(m->getElements(),     nelemsr, colels_);
1371
1372
1373#if     0 && DEBUG_PRESOLVE
1374  presolve_check_costs(model, &colcopy);
1375#endif
1376
1377  // This determines the size of the data structure that contains
1378  // the matrix being postsolved.  Links are taken from the free_list
1379  // to recreate matrix entries that were presolved away,
1380  // and links are added to the free_list when entries created during
1381  // presolve are discarded.  There is never a need to gc this list.
1382  // Naturally, it should contain
1383  // exactly nelems0 entries "in use" when postsolving is done,
1384  // but I don't know whether the matrix could temporarily get
1385  // larger during postsolving.  Substitution into more than two
1386  // rows could do that, in principle.  I am being very conservative
1387  // here by reserving much more than the amount of space I probably need.
1388  // If this guess is wrong, check_free_list may be called.
1389  //  int bufsize = 2*nelems0;
1390
1391  memset(cdone_, -1, ncols0_);
1392  memset(rdone_, -1, nrows0_);
1393
1394  rowduals_ = new double[nrows0_];
1395  ClpDisjointCopyN(si->getRowPrice(), nrows1, rowduals_);
1396
1397  rcosts_ = new double[ncols0_];
1398  ClpDisjointCopyN(si->getReducedCost(), ncols1, rcosts_);
1399  if (maxmin_<0.0) {
1400    // change so will look as if minimize
1401    int i;
1402    for (i=0;i<nrows1;i++)
1403      rowduals_[i] = - rowduals_[i];
1404    for (i=0;i<ncols1;i++) {
1405      rcosts_[i] = - rcosts_[i];
1406    }
1407  }
1408
1409  //ClpDisjointCopyN(si->getRowUpper(), nrows1, rup_);
1410  //ClpDisjointCopyN(si->getRowLower(), nrows1, rlo_);
1411
1412  ClpDisjointCopyN(si->getColSolution(), ncols1, sol_);
1413  si->setDblParam(ClpObjOffset,originalOffset_);
1414
1415  for (int j=0; j<ncols1; j++) {
1416    CoinBigIndex kcs = mcstrt_[j];
1417    CoinBigIndex kce = kcs + hincol_[j];
1418    for (CoinBigIndex k=kcs; k<kce; ++k) {
1419      link_[k] = k+1;
1420    }
1421  }
1422  {
1423    int ml = maxlink_;
1424    for (CoinBigIndex k=nelemsr; k<ml; ++k)
1425      link_[k] = k+1;
1426    link_[ml-1] = NO_LINK;
1427  }
1428  free_list_ = nelemsr;
1429}
1430
1431
Note: See TracBrowser for help on using the repository browser.