Changeset 989 for branches


Ignore:
Timestamp:
Mar 10, 2007 7:05:03 PM (13 years ago)
Author:
forrest
Message:

this may be a mistake

Location:
branches/devel/Clp/src
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • branches/devel/Clp/src/ClpMain.cpp

    r861 r989  
    950950                bool deleteModel2=false;
    951951                ClpSimplex * model2 = models+iModel;
     952                if (dualize) {
     953                  model2 = ((ClpSimplexOther *) model2)->dualOfModel();
     954                  printf("Dual of model has %d rows and %d columns\n",
     955                         model2->numberRows(),model2->numberColumns());
     956                  model2->setOptimizationDirection(1.0);
     957                  preSolve=0; // as picks up from model
     958                }
    952959                if (preSolve) {
    953960                  ClpPresolve pinfo;
  • branches/devel/Clp/src/ClpModel.cpp

    r923 r989  
    7676  lengthNames_(0),
    7777  numberThreads_(0),
     78  specialOptions_(0),
    7879#ifndef CLP_NO_STD
    7980  defaultHandler_(true),
     
    157158  eventHandler_=NULL;
    158159  whatsChanged_=0;
     160  specialOptions_ = 0;
    159161}
    160162//#############################################################################
     
    674676  userPointer_ = rhs.userPointer_;
    675677  scalingFlag_ = rhs.scalingFlag_;
     678  specialOptions_ = rhs.specialOptions_;
    676679  if (trueCopy) {
    677680#ifndef CLP_NO_STD
     
    706709    rowUpper_ = ClpCopyOfArray ( rhs.rowUpper_, numberRows_ );
    707710    columnLower_ = ClpCopyOfArray ( rhs.columnLower_, numberColumns_ );
     711    int scaleLength = ((specialOptions_&131072)==0) ? 1 : 2;
    708712    columnUpper_ = ClpCopyOfArray ( rhs.columnUpper_, numberColumns_ );
    709     rowScale_ = ClpCopyOfArray(rhs.rowScale_,numberRows_);
    710     columnScale_ = ClpCopyOfArray(rhs.columnScale_,numberColumns_);
     713    rowScale_ = ClpCopyOfArray(rhs.rowScale_,numberRows_*scaleLength);
     714    columnScale_ = ClpCopyOfArray(rhs.columnScale_,numberColumns_*scaleLength);
    711715    if (rhs.objective_)
    712716      objective_  = rhs.objective_->clone();
     
    27712775  strParam_[ClpProbName] = rhs->strParam_[ClpProbName];
    27722776#endif
    2773 
     2777  specialOptions_ = rhs->specialOptions_;
    27742778  optimizationDirection_ = rhs->optimizationDirection_;
    27752779  objectiveValue_=rhs->objectiveValue_;
     
    34363440  columnScale_ = NULL;
    34373441}
     3442void
     3443ClpModel::setSpecialOptions(unsigned int value)
     3444{
     3445  specialOptions_=value;
     3446}
     3447/* This creates a coinModel object
     3448 */
     3449CoinModel *
     3450ClpModel::createCoinModel() const
     3451{
     3452  CoinModel * coinModel = new CoinModel();
     3453  CoinPackedMatrix matrixByRow;
     3454  matrixByRow.reverseOrderedCopyOf(*matrix());
     3455  coinModel->setObjectiveOffset(objectiveOffset());
     3456  coinModel->setProblemName(problemName().c_str());
     3457
     3458  // Build by row from scratch
     3459  const double * element = matrixByRow.getElements();
     3460  const int * column = matrixByRow.getIndices();
     3461  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
     3462  const int * rowLength = matrixByRow.getVectorLengths();
     3463  int i;
     3464  for (i=0;i<numberRows_;i++) {
     3465    coinModel->addRow(rowLength[i],column+rowStart[i],
     3466                      element+rowStart[i],rowLower_[i],rowUpper_[i]);
     3467  }
     3468  // Now do column part
     3469  const double * objective = this->objective();
     3470  for (i=0;i<numberColumns_;i++) {
     3471    coinModel->setColumnBounds(i,columnLower_[i],columnUpper_[i]);
     3472    coinModel->setColumnObjective(i,objective[i]);
     3473  }
     3474  for ( i=0;i<numberColumns_;i++) {
     3475    if (isInteger(i))
     3476      coinModel->setColumnIsInteger(i,true);
     3477  }
     3478  // do names
     3479  for (i=0;i<numberRows_;i++) {
     3480    char temp[30];
     3481    strcpy(temp,rowName(i).c_str());
     3482    int length = strlen(temp);
     3483    for (int j=0;j<length;j++) {
     3484      if (temp[j]=='-')
     3485        temp[j]='_';
     3486    }
     3487    coinModel->setRowName(i,temp);
     3488  }
     3489  for (i=0;i<numberColumns_;i++) {
     3490    char temp[30];
     3491    strcpy(temp,columnName(i).c_str());
     3492    int length = strlen(temp);
     3493    for (int j=0;j<length;j++) {
     3494      if (temp[j]=='-')
     3495        temp[j]='_';
     3496    }
     3497    coinModel->setColumnName(i,temp);
     3498  }
     3499  ClpQuadraticObjective * obj = (dynamic_cast< ClpQuadraticObjective*>(objective_));
     3500  if (obj) {
     3501    const CoinPackedMatrix * quadObj = obj->quadraticObjective();
     3502    // add in quadratic
     3503    const double * element = quadObj->getElements();
     3504    const int * row = quadObj->getIndices();
     3505    const CoinBigIndex * columnStart = quadObj->getVectorStarts();
     3506    const int * columnLength = quadObj->getVectorLengths();
     3507    for (i=0;i<numberColumns_;i++) {
     3508      int nels = columnLength[i];
     3509      if (nels) {
     3510        CoinBigIndex start = columnStart[i];
     3511        double constant = coinModel->getColumnObjective(i);
     3512        char temp[100000];
     3513        char temp2[30];
     3514        sprintf(temp,"%g",constant);
     3515        for (CoinBigIndex k=start;k<start+nels;k++) {
     3516          int kColumn = row[k];
     3517          double value = element[k];
     3518#if 1
     3519          // ampl gives twice with assumed 0.5
     3520          if (kColumn<i)
     3521            continue;
     3522          else if (kColumn==i)
     3523            value *= 0.5;
     3524#endif
     3525          if (value==1.0)
     3526            sprintf(temp2,"+%s",coinModel->getColumnName(kColumn));
     3527          else if (value==-1.0)
     3528            sprintf(temp2,"-%s",coinModel->getColumnName(kColumn));
     3529          else if (value>0.0)
     3530            sprintf(temp2,"+%g*%s",value,coinModel->getColumnName(kColumn));
     3531          else
     3532            sprintf(temp2,"%g*%s",value,coinModel->getColumnName(kColumn));
     3533          strcat(temp,temp2);
     3534          assert (strlen(temp)<100000);
     3535        }
     3536        coinModel->setObjective(i,temp);
     3537        if (logLevel()>2)
     3538          printf("el for objective column %s is %s\n",coinModel->getColumnName(i),temp);
     3539      }
     3540    }
     3541  }
     3542  return coinModel;
     3543}
    34383544//#############################################################################
    34393545// Constructors / Destructor / Assignment
  • branches/devel/Clp/src/ClpModel.hpp

    r923 r989  
    268268  void setColumnName(int colIndex, std::string & name) ;
    269269#endif
     270  /** This creates a coinModel object
     271  */
     272  CoinModel * createCoinModel() const;
    270273
    271274    /** Write the problem in MPS format to the specified file.
     
    361364       4 - scaled problem optimal - unscaled problem has primal and dual infeasibilities
    362365       5 - giving up in primal with flagged variables
    363        6 - failed due to empty problem check
     366       6 - failed due to empty problem check 
    364367       7 - postSolve says not optimal
    365368       8 - failed due to bad element check
     
    775778    /// Create C++ lines to get to current state
    776779    void generateCpp( FILE * fp);
     780  /** For advanced options
     781      1 - Don't keep changing infeasibility weight
     782      2 - Keep nonLinearCost round solves
     783      4 - Force outgoing variables to exact bound (primal)
     784      8 - Safe to use dense initial factorization
     785      16 -Just use basic variables for operation if column generation
     786      32 -Clean up with primal before strong branching
     787      64 -Treat problem as feasible until last minute (i.e. minimize infeasibilities)
     788      128 - Switch off all matrix sanity checks
     789      256 - No row copy
     790      512 - If not in values pass, solution guaranteed, skip as much as possible
     791      1024 - In branch and bound
     792      2048 - Don't bother to re-factorize if < 20 iterations
     793      4096 - Skip some optimality checks
     794      8192 - Do Primal when cleaning up primal
     795      16384 - In fast dual (so we can switch off things)
     796      32678 - called from Osi
     797      65356 - keep arrays around as much as possible
     798      131072 - scale factor arrays have inverse values at end
     799      NOTE - many applications can call Clp but there may be some short cuts
     800             which are taken which are not guaranteed safe from all applications.
     801             Vetted applications will have a bit set and the code may test this
     802             At present I expect a few such applications - if too many I will
     803             have to re-think.  It is up to application owner to change the code
     804             if she/he needs these short cuts.  I will not debug unless in Coin
     805             repository.  See COIN_CLP_VETTED comments.
     806      0x01000000 is Cbc (and in branch and bound)
     807      0x02000000 is in a different branch and bound
     808  */
     809#define COIN_CBC_USING_CLP 0x01000000
     810  inline unsigned int specialOptions() const
     811  { return specialOptions_;};
     812  void setSpecialOptions(unsigned int value);
    777813  //@}
    778814
     
    903939  /// Number of threads (not very operational)
    904940  int numberThreads_;
     941  /** For advanced options
     942      See get and set for meaning
     943  */
     944  unsigned int specialOptions_;
    905945  /// Message handler
    906946  CoinMessageHandler * handler_;
  • branches/devel/Clp/src/ClpPackedMatrix.cpp

    r928 r989  
    21052105  const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
    21062106  const double * element = rowCopy->getElements();
    2107   double * rowScale = new double [numberRows];
    2108   double * columnScale = new double [numberColumns];
     2107  int scaleLength = ((model->specialOptions()&131072)==0) ? 1 : 2;
     2108  double * rowScale = new double [numberRows*scaleLength];
     2109  double * columnScale = new double [numberColumns*scaleLength];
    21092110  // we are going to mark bits we are interested in
    21102111  char * usefulRow = new char [numberRows];
     
    24902491    }
    24912492#endif
     2493    if (scaleLength>1) {
     2494      // make copy
     2495      double * inverseScale = rowScale + numberRows;
     2496      for (iRow=0;iRow<numberRows;iRow++)
     2497        inverseScale[iRow] = 1.0/rowScale[iRow] ;
     2498      inverseScale = columnScale + numberColumns;
     2499      for (iColumn=0;iColumn<numberColumns;iColumn++)
     2500        inverseScale[iColumn] = 1.0/columnScale[iColumn] ;
     2501    }
    24922502    model->setRowScale(rowScale);
    24932503    model->setColumnScale(columnScale);
  • branches/devel/Clp/src/ClpSimplex.cpp

    r976 r989  
    5151  largestPrimalError_(0.0),
    5252  largestDualError_(0.0),
    53   largestSolutionError_(0.0),
     53  alphaAccuracy_(-1.0),
    5454  dualBound_(1.0e10),
    5555  alpha_(0.0),
     
    107107  perturbation_(100),
    108108  nonLinearCost_(NULL),
    109   specialOptions_(0),
    110109  lastBadIteration_(-999999),
    111110  lastFlaggedIteration_(-999999),
     
    163162  largestPrimalError_(0.0),
    164163  largestDualError_(0.0),
    165   largestSolutionError_(0.0),
     164  alphaAccuracy_(-1.0),
    166165  dualBound_(1.0e10),
    167166  alpha_(0.0),
     
    219218  perturbation_(100),
    220219  nonLinearCost_(NULL),
    221   specialOptions_(0),
    222220  lastBadIteration_(-999999),
    223221  lastFlaggedIteration_(-999999),
     
    14621460    progressFlag_ |= 1; // making real progress
    14631461  if (sequenceIn_!=sequenceOut_) {
     1462    if (alphaAccuracy_>0.0) {
     1463      double value = fabs(alpha_);
     1464      if (value>1.0)
     1465        alphaAccuracy_ *= value;
     1466      else
     1467        alphaAccuracy_ /= value;
     1468    }
    14641469    //assert( getStatus(sequenceOut_)== basic);
    14651470    setStatus(sequenceIn_,basic);
     
    16061611  largestPrimalError_(0.0),
    16071612  largestDualError_(0.0),
    1608   largestSolutionError_(0.0),
     1613  alphaAccuracy_(-1.0),
    16091614  dualBound_(1.0e10),
    16101615  alpha_(0.0),
     
    16621667  perturbation_(100),
    16631668  nonLinearCost_(NULL),
    1664   specialOptions_(0),
    16651669  lastBadIteration_(-999999),
    16661670  lastFlaggedIteration_(-999999),
     
    16901694  primalColumnPivot_ = NULL;
    16911695  gutsOfDelete(0);
    1692   specialOptions_ =0;
    16931696  delete nonLinearCost_;
    16941697  nonLinearCost_ = NULL;
     
    17121715  largestPrimalError_(0.0),
    17131716  largestDualError_(0.0),
    1714   largestSolutionError_(0.0),
     1717  alphaAccuracy_(-1.0),
    17151718  dualBound_(1.0e10),
    17161719  alpha_(0.0),
     
    17681771  perturbation_(100),
    17691772  nonLinearCost_(NULL),
    1770   specialOptions_(0),
    17711773  lastBadIteration_(-999999),
    17721774  lastFlaggedIteration_(-999999),
     
    18071809  if (this != &rhs) {
    18081810    gutsOfDelete(0);
    1809     specialOptions_=0;
    18101811    delete nonLinearCost_;
    18111812    nonLinearCost_ = NULL;
     
    19081909  largestPrimalError_ = rhs.largestPrimalError_;
    19091910  largestDualError_ = rhs.largestDualError_;
    1910   largestSolutionError_ = rhs.largestSolutionError_;
     1911  alphaAccuracy_ = rhs.alphaAccuracy_;
    19111912  dualBound_ = rhs.dualBound_;
    19121913  alpha_ = rhs.alpha_;
     
    19441945  perturbation_ = rhs.perturbation_;
    19451946  infeasibilityCost_ = rhs.infeasibilityCost_;
    1946   specialOptions_ = rhs.specialOptions_;
    19471947  lastBadIteration_ = rhs.lastBadIteration_;
    19481948  lastFlaggedIteration_ = rhs.lastFlaggedIteration_;
     
    20482048  objectiveValue_ = 0.0;
    20492049  // now look at primal solution
    2050   columnPrimalInfeasibility_=0.0;
    2051   columnPrimalSequence_=-1;
    2052   rowPrimalInfeasibility_=0.0;
    2053   rowPrimalSequence_=-1;
    2054   largestSolutionError_=0.0;
    20552050  solution = rowActivityWork_;
    20562051  sumPrimalInfeasibilities_=0.0;
     
    20782073      numberPrimalInfeasibilities_ ++;
    20792074    }
    2080     if (infeasibility>rowPrimalInfeasibility_) {
    2081       rowPrimalInfeasibility_=infeasibility;
    2082       rowPrimalSequence_=iRow;
    2083     }
    20842075    infeasibility = fabs(rowActivities[iRow]-solution[iRow]);
    2085     if (infeasibility>largestSolutionError_)
    2086       largestSolutionError_=infeasibility;
    20872076  }
    20882077  // Check any infeasibilities from dynamic rows
     
    20992088        infeasibility=columnLowerWork_[iColumn]-solution[iColumn];
    21002089      }
    2101       if (infeasibility>columnPrimalInfeasibility_) {
    2102         columnPrimalInfeasibility_=infeasibility;
    2103         columnPrimalSequence_=iColumn;
    2104       }
    21052090      if (infeasibility>primalTolerance) {
    21062091        sumPrimalInfeasibilities_ += infeasibility-primalTolerance_;
     
    21102095      }
    21112096      infeasibility = fabs(columnActivities[iColumn]-solution[iColumn]);
    2112       if (infeasibility>largestSolutionError_)
    2113         largestSolutionError_=infeasibility;
    21142097    }
    21152098  } else {
     
    21262109        infeasibility=columnLowerWork_[iColumn]-solution[iColumn];
    21272110      }
    2128       if (infeasibility>columnPrimalInfeasibility_) {
    2129         columnPrimalInfeasibility_=infeasibility;
    2130         columnPrimalSequence_=iColumn;
    2131       }
    21322111      if (infeasibility>primalTolerance) {
    21332112        sumPrimalInfeasibilities_ += infeasibility-primalTolerance_;
     
    21372116      }
    21382117      infeasibility = fabs(columnActivities[iColumn]-solution[iColumn]);
    2139       if (infeasibility>largestSolutionError_)
    2140         largestSolutionError_=infeasibility;
    21412118    }
    21422119  }
     
    21522129  numberDualInfeasibilities_=0;
    21532130  numberDualInfeasibilitiesWithoutFree_=0;
    2154   columnDualInfeasibility_=0.0;
    2155   columnDualSequence_=-1;
    21562131  if (matrix_->skipDualCheck()&&algorithm_>0&&problemStatus_==-2) {
    21572132    // pretend we found dual infeasibilities
     
    21612136    return;
    21622137  }
    2163   rowDualInfeasibility_=0.0;
    2164   rowDualSequence_=-1;
    21652138  int firstFreePrimal = -1;
    21662139  int firstFreeDual = -1;
    21672140  int numberSuperBasicWithDj=0;
    2168   primalToleranceToGetOptimal_=CoinMax(rowPrimalInfeasibility_,
    2169                                    columnPrimalInfeasibility_);
    2170   remainingDualInfeasibility_=0.0;
    21712141  double relaxedTolerance=dualTolerance_;
    21722142  // we can't really trust infeasibilities if there is dual error
     
    22012171        if (value<0.0) {
    22022172          value = - value;
    2203           if (value>columnDualInfeasibility_) {
    2204             columnDualInfeasibility_=value;
    2205             columnDualSequence_=iColumn;
    2206           }
    22072173          if (value>dualTolerance_) {
    22082174            if (getColumnStatus(iColumn) != isFree) {
     
    22222188              }
    22232189            }
    2224             // maybe we can make feasible by increasing tolerance
    2225             if (distanceUp<largeValue_) {
    2226               if (distanceUp>primalToleranceToGetOptimal_)
    2227                 primalToleranceToGetOptimal_=distanceUp;
    2228             } else {
    2229               //gap too big for any tolerance
    2230               remainingDualInfeasibility_=
    2231                 CoinMax(remainingDualInfeasibility_,value);
    2232             }
    22332190          }
    22342191        }
     
    22382195        // should not be positive
    22392196        if (value>0.0) {
    2240           if (value>columnDualInfeasibility_) {
    2241             columnDualInfeasibility_=value;
    2242             columnDualSequence_=iColumn;
    2243           }
    22442197          if (value>dualTolerance_) {
    22452198            sumDualInfeasibilities_ += value-dualTolerance_;
     
    22502203              numberDualInfeasibilitiesWithoutFree_ ++;
    22512204            // maybe we can make feasible by increasing tolerance
    2252             if (distanceDown<largeValue_&&
    2253                 distanceDown>primalToleranceToGetOptimal_)
    2254               primalToleranceToGetOptimal_=distanceDown;
    22552205          }
    22562206        }
     
    22782228        if (value<0.0) {
    22792229          value = - value;
    2280           if (value>rowDualInfeasibility_) {
    2281             rowDualInfeasibility_=value;
    2282             rowDualSequence_=iRow;
     2230          if (value>dualTolerance_) {
     2231            sumDualInfeasibilities_ += value-dualTolerance_;
     2232            if (value>relaxedTolerance)
     2233              sumOfRelaxedDualInfeasibilities_ += value-relaxedTolerance;
     2234            numberDualInfeasibilities_ ++;
     2235            if (getRowStatus(iRow) != isFree)
     2236              numberDualInfeasibilitiesWithoutFree_ ++;
    22832237          }
     2238        }
     2239      }
     2240      if (distanceDown>primalTolerance_) {
     2241        double value = rowReducedCost_[iRow];
     2242        // should not be positive
     2243        if (value>0.0) {
    22842244          if (value>dualTolerance_) {
    22852245            sumDualInfeasibilities_ += value-dualTolerance_;
     
    22902250              numberDualInfeasibilitiesWithoutFree_ ++;
    22912251            // maybe we can make feasible by increasing tolerance
    2292             if (distanceUp<largeValue_) {
    2293               if (distanceUp>primalToleranceToGetOptimal_)
    2294                 primalToleranceToGetOptimal_=distanceUp;
    2295             } else {
    2296               //gap too big for any tolerance
    2297               remainingDualInfeasibility_=
    2298                 CoinMax(remainingDualInfeasibility_,value);
    2299             }
    2300           }
    2301         }
    2302       }
    2303       if (distanceDown>primalTolerance_) {
    2304         double value = rowReducedCost_[iRow];
    2305         // should not be positive
    2306         if (value>0.0) {
    2307           if (value>rowDualInfeasibility_) {
    2308             rowDualInfeasibility_=value;
    2309             rowDualSequence_=iRow;
    2310           }
    2311           if (value>dualTolerance_) {
    2312             sumDualInfeasibilities_ += value-dualTolerance_;
    2313             if (value>relaxedTolerance)
    2314               sumOfRelaxedDualInfeasibilities_ += value-relaxedTolerance;
    2315             numberDualInfeasibilities_ ++;
    2316             if (getRowStatus(iRow) != isFree)
    2317               numberDualInfeasibilitiesWithoutFree_ ++;
    2318             // maybe we can make feasible by increasing tolerance
    2319             if (distanceDown<largeValue_&&
    2320                 distanceDown>primalToleranceToGetOptimal_)
    2321               primalToleranceToGetOptimal_=distanceDown;
    23222252          }
    23232253        }
     
    25272457  }
    25282458}
     2459//static int scale_times[]={0,0,0,0};
    25292460bool
    25302461ClpSimplex::createRim(int what,bool makeRowCopy, int startFinishOptions)
     
    26272558      rowCopy_ = matrix_->reverseOrderedCopy();
    26282559    }
     2560#if 0
     2561    if (what==63&&(specialOptions_&131072)!=0) {
     2562      int k=rowScale_ ? 1 : 0;
     2563      if (oldMatrix)
     2564        k+=2;
     2565      scale_times[k]++;
     2566      if ((scale_times[0]+scale_times[1]+scale_times[2]+scale_times[3])%1000==0)
     2567        printf("scale counts %d %d %d %d\n",
     2568               scale_times[0],scale_times[1],scale_times[2],scale_times[3]);
     2569    }
     2570#endif
    26292571    // do scaling if needed
    26302572    if (!oldMatrix&&scalingFlag_<0) {
    26312573      if (scalingFlag_<0&&rowScale_) {
    2632         if (handler_->logLevel()>0)
     2574        //if (handler_->logLevel()>0)
    26332575          printf("How did we get scalingFlag_ %d and non NULL rowScale_? - switching off scaling\n",
    26342576                 scalingFlag_);
     
    28612803        // If scaled then do all columns later in one loop
    28622804        if (what!=63) {
    2863           for (i=0;i<numberColumns_;i++) {
    2864             double multiplier = rhsScale_/columnScale[i];
    2865             double lowerValue = columnLower_[i];
    2866             double upperValue = columnUpper_[i];
    2867             if (lowerValue>-1.0e20) {
    2868               columnLowerWork_[i]=lowerValue*multiplier;
    2869               if (upperValue>=1.0e20) {
    2870                 columnUpperWork_[i]=COIN_DBL_MAX;
    2871               } else {
    2872                 columnUpperWork_[i]=upperValue*multiplier;
    2873                 if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    2874                   if (columnLowerWork_[i]>=0.0) {
    2875                     columnUpperWork_[i] = columnLowerWork_[i];
    2876                   } else if (columnUpperWork_[i]<=0.0) {
    2877                     columnLowerWork_[i] = columnUpperWork_[i];
    2878                   } else {
    2879                     columnUpperWork_[i] = 0.0;
    2880                     columnLowerWork_[i] = 0.0;
    2881                   }
    2882                 }
    2883               }
    2884             } else if (upperValue<1.0e20) {
    2885               columnLowerWork_[i]=-COIN_DBL_MAX;
    2886               columnUpperWork_[i]=upperValue*multiplier;
    2887             } else {
    2888               // free
    2889               columnLowerWork_[i]=-COIN_DBL_MAX;
    2890               columnUpperWork_[i]=COIN_DBL_MAX;
    2891             }
    2892           }
    2893         }
     2805          if ((specialOptions_&131072)==0) {
     2806            for (i=0;i<numberColumns_;i++) {
     2807              double multiplier = rhsScale_/columnScale[i];
     2808              double lowerValue = columnLower_[i];
     2809              double upperValue = columnUpper_[i];
     2810              if (lowerValue>-1.0e20) {
     2811                columnLowerWork_[i]=lowerValue*multiplier;
     2812                if (upperValue>=1.0e20) {
     2813                  columnUpperWork_[i]=COIN_DBL_MAX;
     2814                } else {
     2815                  columnUpperWork_[i]=upperValue*multiplier;
     2816                  if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     2817                    if (columnLowerWork_[i]>=0.0) {
     2818                      columnUpperWork_[i] = columnLowerWork_[i];
     2819                    } else if (columnUpperWork_[i]<=0.0) {
     2820                      columnLowerWork_[i] = columnUpperWork_[i];
     2821                    } else {
     2822                      columnUpperWork_[i] = 0.0;
     2823                      columnLowerWork_[i] = 0.0;
     2824                    }
     2825                  }
     2826                }
     2827              } else if (upperValue<1.0e20) {
     2828                columnLowerWork_[i]=-COIN_DBL_MAX;
     2829                columnUpperWork_[i]=upperValue*multiplier;
     2830              } else {
     2831                // free
     2832                columnLowerWork_[i]=-COIN_DBL_MAX;
     2833                columnUpperWork_[i]=COIN_DBL_MAX;
     2834              }
     2835            }
     2836          } else {
     2837            const double * inverseScale = columnScale_+numberColumns_;
     2838            for (i=0;i<numberColumns_;i++) {
     2839              double multiplier = rhsScale_*inverseScale[i];
     2840              double lowerValue = columnLower_[i];
     2841              double upperValue = columnUpper_[i];
     2842              if (lowerValue>-1.0e20) {
     2843                columnLowerWork_[i]=lowerValue*multiplier;
     2844                if (upperValue>=1.0e20) {
     2845                  columnUpperWork_[i]=COIN_DBL_MAX;
     2846                } else {
     2847                  columnUpperWork_[i]=upperValue*multiplier;
     2848                  if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     2849                    if (columnLowerWork_[i]>=0.0) {
     2850                      columnUpperWork_[i] = columnLowerWork_[i];
     2851                    } else if (columnUpperWork_[i]<=0.0) {
     2852                      columnLowerWork_[i] = columnUpperWork_[i];
     2853                    } else {
     2854                      columnUpperWork_[i] = 0.0;
     2855                      columnLowerWork_[i] = 0.0;
     2856                    }
     2857                  }
     2858                }
     2859              } else if (upperValue<1.0e20) {
     2860                columnLowerWork_[i]=-COIN_DBL_MAX;
     2861                columnUpperWork_[i]=upperValue*multiplier;
     2862              } else {
     2863                // free
     2864                columnLowerWork_[i]=-COIN_DBL_MAX;
     2865                columnUpperWork_[i]=COIN_DBL_MAX;
     2866              }
     2867            }
     2868          }
     2869        }
    28942870        for (i=0;i<numberRows_;i++) {
    28952871          double multiplier = rhsScale_*rowScale[i];
     
    31783154        double primalTolerance=dblParam_[ClpPrimalTolerance];
    31793155        // on entry
    3180         for (i=0;i<numberColumns_;i++) {
    3181           CoinAssert(fabs(obj[i])<1.0e25);
    3182           double scaleFactor = columnScale_[i];
    3183           double multiplier = rhsScale_/scaleFactor;
    3184           scaleFactor *= direction;
    3185           objectiveWork_[i] = obj[i]*scaleFactor;
    3186           reducedCostWork_[i] = reducedCost_[i]*scaleFactor;
    3187           double lowerValue = columnLower_[i];
    3188           double upperValue = columnUpper_[i];
    3189           if (lowerValue>-1.0e20) {
    3190             columnLowerWork_[i]=lowerValue*multiplier;
    3191             if (upperValue>=1.0e20) {
    3192               columnUpperWork_[i]=COIN_DBL_MAX;
    3193             } else {
    3194               columnUpperWork_[i]=upperValue*multiplier;
    3195               if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
    3196                 if (columnLowerWork_[i]>=0.0) {
    3197                   columnUpperWork_[i] = columnLowerWork_[i];
    3198                 } else if (columnUpperWork_[i]<=0.0) {
    3199                   columnLowerWork_[i] = columnUpperWork_[i];
    3200                 } else {
    3201                   columnUpperWork_[i] = 0.0;
    3202                   columnLowerWork_[i] = 0.0;
    3203                 }
    3204               }
    3205             }
    3206           } else if (upperValue<1.0e20) {
    3207             columnLowerWork_[i]=-COIN_DBL_MAX;
    3208             columnUpperWork_[i]=upperValue*multiplier;
    3209           } else {
    3210             // free
    3211             columnLowerWork_[i]=-COIN_DBL_MAX;
    3212             columnUpperWork_[i]=COIN_DBL_MAX;
    3213           }
    3214           double value = columnActivity_[i] * multiplier;
    3215           if (fabs(value)>1.0e20) {
    3216             //printf("bad value of %g for column %d\n",value,i);
    3217             setColumnStatus(i,superBasic);
    3218             if (columnUpperWork_[i]<0.0) {
    3219               value=columnUpperWork_[i];
    3220             } else if (columnLowerWork_[i]>0.0) {
    3221               value=columnLowerWork_[i];
    3222             } else {
    3223               value=0.0;
    3224             }
    3225           }
    3226           columnActivityWork_[i] = value;
    3227         }
    3228         for (i=0;i<numberRows_;i++) {
    3229           dual_[i] /= rowScale_[i];
    3230           dual_[i] *= objectiveScale_;
    3231           rowReducedCost_[i] = dual_[i];
    3232           double multiplier = rhsScale_*rowScale_[i];
    3233           double value = rowActivity_[i] * multiplier;
    3234           if (fabs(value)>1.0e20) {
    3235             //printf("bad value of %g for row %d\n",value,i);
    3236             setRowStatus(i,superBasic);
    3237             if (rowUpperWork_[i]<0.0) {
    3238               value=rowUpperWork_[i];
    3239             } else if (rowLowerWork_[i]>0.0) {
    3240               value=rowLowerWork_[i];
    3241             } else {
    3242               value=0.0;
    3243             }
    3244           }
    3245           rowActivityWork_[i] = value;
    3246         }
     3156        if ((specialOptions_&131072)==0) {
     3157          for (i=0;i<numberColumns_;i++) {
     3158            CoinAssert(fabs(obj[i])<1.0e25);
     3159            double scaleFactor = columnScale_[i];
     3160            double multiplier = rhsScale_/scaleFactor;
     3161            scaleFactor *= direction;
     3162            objectiveWork_[i] = obj[i]*scaleFactor;
     3163            reducedCostWork_[i] = reducedCost_[i]*scaleFactor;
     3164            double lowerValue = columnLower_[i];
     3165            double upperValue = columnUpper_[i];
     3166            if (lowerValue>-1.0e20) {
     3167              columnLowerWork_[i]=lowerValue*multiplier;
     3168              if (upperValue>=1.0e20) {
     3169                columnUpperWork_[i]=COIN_DBL_MAX;
     3170              } else {
     3171                columnUpperWork_[i]=upperValue*multiplier;
     3172                if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3173                  if (columnLowerWork_[i]>=0.0) {
     3174                    columnUpperWork_[i] = columnLowerWork_[i];
     3175                  } else if (columnUpperWork_[i]<=0.0) {
     3176                    columnLowerWork_[i] = columnUpperWork_[i];
     3177                  } else {
     3178                    columnUpperWork_[i] = 0.0;
     3179                    columnLowerWork_[i] = 0.0;
     3180                  }
     3181                }
     3182              }
     3183            } else if (upperValue<1.0e20) {
     3184              columnLowerWork_[i]=-COIN_DBL_MAX;
     3185              columnUpperWork_[i]=upperValue*multiplier;
     3186            } else {
     3187              // free
     3188              columnLowerWork_[i]=-COIN_DBL_MAX;
     3189              columnUpperWork_[i]=COIN_DBL_MAX;
     3190            }
     3191            double value = columnActivity_[i] * multiplier;
     3192            if (fabs(value)>1.0e20) {
     3193              //printf("bad value of %g for column %d\n",value,i);
     3194              setColumnStatus(i,superBasic);
     3195              if (columnUpperWork_[i]<0.0) {
     3196                value=columnUpperWork_[i];
     3197              } else if (columnLowerWork_[i]>0.0) {
     3198                value=columnLowerWork_[i];
     3199              } else {
     3200                value=0.0;
     3201              }
     3202            }
     3203            columnActivityWork_[i] = value;
     3204          }
     3205          for (i=0;i<numberRows_;i++) {
     3206            dual_[i] /= rowScale_[i];
     3207            dual_[i] *= objectiveScale_;
     3208            rowReducedCost_[i] = dual_[i];
     3209            double multiplier = rhsScale_*rowScale_[i];
     3210            double value = rowActivity_[i] * multiplier;
     3211            if (fabs(value)>1.0e20) {
     3212              //printf("bad value of %g for row %d\n",value,i);
     3213              setRowStatus(i,superBasic);
     3214              if (rowUpperWork_[i]<0.0) {
     3215                value=rowUpperWork_[i];
     3216              } else if (rowLowerWork_[i]>0.0) {
     3217                value=rowLowerWork_[i];
     3218              } else {
     3219                value=0.0;
     3220              }
     3221            }
     3222            rowActivityWork_[i] = value;
     3223          }
     3224        } else {
     3225          const double * inverseScale = columnScale_+numberColumns_;
     3226          for (i=0;i<numberColumns_;i++) {
     3227            CoinAssert(fabs(obj[i])<1.0e25);
     3228            double scaleFactor = columnScale_[i];
     3229            double multiplier = rhsScale_*inverseScale[i];
     3230            scaleFactor *= direction;
     3231            objectiveWork_[i] = obj[i]*scaleFactor;
     3232            reducedCostWork_[i] = reducedCost_[i]*scaleFactor;
     3233            double lowerValue = columnLower_[i];
     3234            double upperValue = columnUpper_[i];
     3235            if (lowerValue>-1.0e20) {
     3236              columnLowerWork_[i]=lowerValue*multiplier;
     3237              if (upperValue>=1.0e20) {
     3238                columnUpperWork_[i]=COIN_DBL_MAX;
     3239              } else {
     3240                columnUpperWork_[i]=upperValue*multiplier;
     3241                if (fabs(columnUpperWork_[i]-columnLowerWork_[i])<=primalTolerance) {
     3242                  if (columnLowerWork_[i]>=0.0) {
     3243                    columnUpperWork_[i] = columnLowerWork_[i];
     3244                  } else if (columnUpperWork_[i]<=0.0) {
     3245                    columnLowerWork_[i] = columnUpperWork_[i];
     3246                  } else {
     3247                    columnUpperWork_[i] = 0.0;
     3248                    columnLowerWork_[i] = 0.0;
     3249                  }
     3250                }
     3251              }
     3252            } else if (upperValue<1.0e20) {
     3253              columnLowerWork_[i]=-COIN_DBL_MAX;
     3254              columnUpperWork_[i]=upperValue*multiplier;
     3255            } else {
     3256              // free
     3257              columnLowerWork_[i]=-COIN_DBL_MAX;
     3258              columnUpperWork_[i]=COIN_DBL_MAX;
     3259            }
     3260            double value = columnActivity_[i] * multiplier;
     3261            if (fabs(value)>1.0e20) {
     3262              //printf("bad value of %g for column %d\n",value,i);
     3263              setColumnStatus(i,superBasic);
     3264              if (columnUpperWork_[i]<0.0) {
     3265                value=columnUpperWork_[i];
     3266              } else if (columnLowerWork_[i]>0.0) {
     3267                value=columnLowerWork_[i];
     3268              } else {
     3269                value=0.0;
     3270              }
     3271            }
     3272            columnActivityWork_[i] = value;
     3273          }
     3274          inverseScale = rowScale_+numberRows_;
     3275          for (i=0;i<numberRows_;i++) {
     3276            dual_[i] *= inverseScale[i];
     3277            dual_[i] *= objectiveScale_;
     3278            rowReducedCost_[i] = dual_[i];
     3279            double multiplier = rhsScale_*rowScale_[i];
     3280            double value = rowActivity_[i] * multiplier;
     3281            if (fabs(value)>1.0e20) {
     3282              //printf("bad value of %g for row %d\n",value,i);
     3283              setRowStatus(i,superBasic);
     3284              if (rowUpperWork_[i]<0.0) {
     3285                value=rowUpperWork_[i];
     3286              } else if (rowLowerWork_[i]>0.0) {
     3287                value=rowLowerWork_[i];
     3288              } else {
     3289                value=0.0;
     3290              }
     3291            }
     3292            rowActivityWork_[i] = value;
     3293          }
     3294        }
    32473295      } else if (objectiveScale_!=1.0||rhsScale_!=1.0) {
    32483296        // on entry
     
    34723520ClpSimplex::deleteRim(int getRidOfFactorizationData)
    34733521{
     3522  // Just possible empty problem
     3523  int numberRows=numberRows_;
     3524  int numberColumns=numberColumns_;
     3525  if (!numberRows||!numberColumns) {
     3526    numberRows=0;
     3527    numberColumns=0;
     3528  }
    34743529  if (!auxiliaryModel_) {
    34753530    int i;
     
    34833538    {
    34843539      int nBad=0;
    3485       for (i=0;i<numberColumns_;i++) {
     3540      for (i=0;i<numberColumns;i++) {
    34863541        if (lower_[i]==upper_[i]&&getColumnStatus(i)==basic)
    34873542          nBad++;
     
    35003555      double scaleC = 1.0/objectiveScale_;
    35013556      double scaleR = 1.0/rhsScale_;
    3502       for (i=0;i<numberColumns_;i++) {
    3503         double scaleFactor = columnScale_[i];
    3504         double valueScaled = columnActivityWork_[i];
    3505         double lowerScaled = columnLowerWork_[i];
    3506         double upperScaled = columnUpperWork_[i];
    3507         if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    3508           if (valueScaled<lowerScaled-primalTolerance_||
    3509               valueScaled>upperScaled+primalTolerance_)
    3510             numberPrimalScaled++;
    3511           else
    3512             upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    3513         }
    3514         columnActivity_[i] = valueScaled*scaleFactor*scaleR;
    3515         double value = columnActivity_[i];
    3516         if (value<columnLower_[i]-primalTolerance_)
    3517           numberPrimalUnscaled++;
    3518         else if (value>columnUpper_[i]+primalTolerance_)
    3519           numberPrimalUnscaled++;
    3520         double valueScaledDual = reducedCostWork_[i];
    3521         if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    3522           numberDualScaled++;
    3523         if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    3524           numberDualScaled++;
    3525         reducedCost_[i] = (valueScaledDual*scaleC)/scaleFactor;
    3526         double valueDual = reducedCost_[i];
    3527         if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    3528           numberDualUnscaled++;
    3529         if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    3530           numberDualUnscaled++;
    3531       }
    3532       for (i=0;i<numberRows_;i++) {
    3533         double scaleFactor = rowScale_[i];
    3534         double valueScaled = rowActivityWork_[i];
    3535         double lowerScaled = rowLowerWork_[i];
    3536         double upperScaled = rowUpperWork_[i];
    3537         if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
    3538           if (valueScaled<lowerScaled-primalTolerance_||
    3539               valueScaled>upperScaled+primalTolerance_)
    3540             numberPrimalScaled++;
    3541           else
    3542             upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
    3543         }
    3544         rowActivity_[i] = (valueScaled*scaleR)/scaleFactor;
    3545         double value = rowActivity_[i];
    3546         if (value<rowLower_[i]-primalTolerance_)
    3547           numberPrimalUnscaled++;
    3548         else if (value>rowUpper_[i]+primalTolerance_)
    3549           numberPrimalUnscaled++;
    3550         double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
    3551         if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
    3552           numberDualScaled++;
    3553         if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
    3554           numberDualScaled++;
    3555         dual_[i] *= scaleFactor*scaleC;
    3556         double valueDual = dual_[i];
    3557         if (rowObjective_)
    3558           valueDual += rowObjective_[i];
    3559         if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
    3560           numberDualUnscaled++;
    3561         if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
    3562           numberDualUnscaled++;
     3557      if ((specialOptions_&131072)==0) {
     3558        for (i=0;i<numberColumns;i++) {
     3559          double scaleFactor = columnScale_[i];
     3560          double valueScaled = columnActivityWork_[i];
     3561          double lowerScaled = columnLowerWork_[i];
     3562          double upperScaled = columnUpperWork_[i];
     3563          if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     3564            if (valueScaled<lowerScaled-primalTolerance_||
     3565                valueScaled>upperScaled+primalTolerance_)
     3566              numberPrimalScaled++;
     3567            else
     3568              upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     3569          }
     3570          columnActivity_[i] = valueScaled*scaleFactor*scaleR;
     3571          double value = columnActivity_[i];
     3572          if (value<columnLower_[i]-primalTolerance_)
     3573            numberPrimalUnscaled++;
     3574          else if (value>columnUpper_[i]+primalTolerance_)
     3575            numberPrimalUnscaled++;
     3576          double valueScaledDual = reducedCostWork_[i];
     3577          if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     3578            numberDualScaled++;
     3579          if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     3580            numberDualScaled++;
     3581          reducedCost_[i] = (valueScaledDual*scaleC)/scaleFactor;
     3582          double valueDual = reducedCost_[i];
     3583          if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     3584            numberDualUnscaled++;
     3585          if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     3586            numberDualUnscaled++;
     3587        }
     3588        for (i=0;i<numberRows;i++) {
     3589          double scaleFactor = rowScale_[i];
     3590          double valueScaled = rowActivityWork_[i];
     3591          double lowerScaled = rowLowerWork_[i];
     3592          double upperScaled = rowUpperWork_[i];
     3593          if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     3594            if (valueScaled<lowerScaled-primalTolerance_||
     3595                valueScaled>upperScaled+primalTolerance_)
     3596              numberPrimalScaled++;
     3597            else
     3598              upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     3599          }
     3600          rowActivity_[i] = (valueScaled*scaleR)/scaleFactor;
     3601          double value = rowActivity_[i];
     3602          if (value<rowLower_[i]-primalTolerance_)
     3603            numberPrimalUnscaled++;
     3604          else if (value>rowUpper_[i]+primalTolerance_)
     3605            numberPrimalUnscaled++;
     3606          double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
     3607          if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     3608            numberDualScaled++;
     3609          if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     3610            numberDualScaled++;
     3611          dual_[i] *= scaleFactor*scaleC;
     3612          double valueDual = dual_[i];
     3613          if (rowObjective_)
     3614            valueDual += rowObjective_[i];
     3615          if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     3616            numberDualUnscaled++;
     3617          if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     3618            numberDualUnscaled++;
     3619        }
     3620      } else {
     3621        const double * inverseScale = columnScale_+numberColumns;
     3622        for (i=0;i<numberColumns;i++) {
     3623          double scaleFactor = columnScale_[i];
     3624          double valueScaled = columnActivityWork_[i];
     3625          double lowerScaled = columnLowerWork_[i];
     3626          double upperScaled = columnUpperWork_[i];
     3627          if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     3628            if (valueScaled<lowerScaled-primalTolerance_||
     3629                valueScaled>upperScaled+primalTolerance_)
     3630              numberPrimalScaled++;
     3631            else
     3632              upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     3633          }
     3634          columnActivity_[i] = valueScaled*scaleFactor*scaleR;
     3635          double value = columnActivity_[i];
     3636          if (value<columnLower_[i]-primalTolerance_)
     3637            numberPrimalUnscaled++;
     3638          else if (value>columnUpper_[i]+primalTolerance_)
     3639            numberPrimalUnscaled++;
     3640          double valueScaledDual = reducedCostWork_[i];
     3641          if (valueScaled>columnLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     3642            numberDualScaled++;
     3643          if (valueScaled<columnUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     3644            numberDualScaled++;
     3645          reducedCost_[i] = (valueScaledDual*scaleC)*inverseScale[i];
     3646          double valueDual = reducedCost_[i];
     3647          if (value>columnLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     3648            numberDualUnscaled++;
     3649          if (value<columnUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     3650            numberDualUnscaled++;
     3651        }
     3652        inverseScale = rowScale_+numberRows;
     3653        for (i=0;i<numberRows;i++) {
     3654          double scaleFactor = rowScale_[i];
     3655          double valueScaled = rowActivityWork_[i];
     3656          double lowerScaled = rowLowerWork_[i];
     3657          double upperScaled = rowUpperWork_[i];
     3658          if (lowerScaled>-1.0e20||upperScaled<1.0e20) {
     3659            if (valueScaled<lowerScaled-primalTolerance_||
     3660                valueScaled>upperScaled+primalTolerance_)
     3661              numberPrimalScaled++;
     3662            else
     3663              upperOut_ = CoinMax(upperOut_,CoinMin(valueScaled-lowerScaled,upperScaled-valueScaled));
     3664          }
     3665          rowActivity_[i] = (valueScaled*scaleR)*inverseScale[i];
     3666          double value = rowActivity_[i];
     3667          if (value<rowLower_[i]-primalTolerance_)
     3668            numberPrimalUnscaled++;
     3669          else if (value>rowUpper_[i]+primalTolerance_)
     3670            numberPrimalUnscaled++;
     3671          double valueScaledDual = dual_[i]+rowObjectiveWork_[i];;
     3672          if (valueScaled>rowLowerWork_[i]+primalTolerance_&&valueScaledDual>dualTolerance_)
     3673            numberDualScaled++;
     3674          if (valueScaled<rowUpperWork_[i]-primalTolerance_&&valueScaledDual<-dualTolerance_)
     3675            numberDualScaled++;
     3676          dual_[i] *= scaleFactor*scaleC;
     3677          double valueDual = dual_[i];
     3678          if (rowObjective_)
     3679            valueDual += rowObjective_[i];
     3680          if (value>rowLower_[i]+primalTolerance_&&valueDual>dualTolerance_)
     3681            numberDualUnscaled++;
     3682          if (value<rowUpper_[i]-primalTolerance_&&valueDual<-dualTolerance_)
     3683            numberDualUnscaled++;
     3684        }
    35633685      }
    35643686      if (!problemStatus_&&!secondaryStatus_) {
     
    35753697      }
    35763698      if (problemStatus_==2) {
    3577         for (i=0;i<numberColumns_;i++) {
     3699        for (i=0;i<numberColumns;i++) {
    35783700          ray_[i] *= columnScale_[i];
    35793701        }
    35803702      } else if (problemStatus_==1&&ray_) {
    3581         for (i=0;i<numberRows_;i++) {
     3703        for (i=0;i<numberRows;i++) {
    35823704          ray_[i] *= rowScale_[i];
    35833705        }
     
    35913713      double scaleC = 1.0/objectiveScale_;
    35923714      double scaleR = 1.0/rhsScale_;
    3593       for (i=0;i<numberColumns_;i++) {
     3715      for (i=0;i<numberColumns;i++) {
    35943716        double valueScaled = columnActivityWork_[i];
    35953717        double lowerScaled = columnLowerWork_[i];
     
    36203742          numberDualUnscaled++;
    36213743      }
    3622       for (i=0;i<numberRows_;i++) {
     3744      for (i=0;i<numberRows;i++) {
    36233745        double valueScaled = rowActivityWork_[i];
    36243746        double lowerScaled = rowLowerWork_[i];
     
    36653787    } else {
    36663788      if (columnActivityWork_) {
    3667         for (i=0;i<numberColumns_;i++) {
     3789        for (i=0;i<numberColumns;i++) {
    36683790          double value = columnActivityWork_[i];
    36693791          double lower = columnLowerWork_[i];
     
    36763798          reducedCost_[i] = reducedCostWork_[i];
    36773799        }
    3678         for (i=0;i<numberRows_;i++) {
     3800        for (i=0;i<numberRows;i++) {
    36793801          double value = rowActivityWork_[i];
    36803802          double lower = rowLowerWork_[i];
     
    36953817    if (optimizationDirection_!=1.0) {
    36963818      // and modify all dual signs
    3697       for (i=0;i<numberColumns_;i++)
     3819      for (i=0;i<numberColumns;i++)
    36983820        reducedCost_[i] *= optimizationDirection_;
    3699       for (i=0;i<numberRows_;i++)
     3821      for (i=0;i<numberRows;i++)
    37003822        dual_[i] *= optimizationDirection_;
    37013823    }
     
    37313853      if (auxiliaryModel_->rowScale_) {
    37323854        const double * scale = auxiliaryModel_->columnScale_;
    3733         const double * inverseScale = scale + numberColumns_;
    3734         for (i=0;i<numberColumns_;i++) {
     3855        const double * inverseScale = scale + numberColumns;
     3856        for (i=0;i<numberColumns;i++) {
    37353857          columnActivity_[i] = auxiliaryModel_->columnActivityWork_[i]*scale[i];
    37363858          reducedCost_[i] = auxiliaryModel_->reducedCostWork_[i]*inverseScale[i];
    37373859        }
    37383860        scale = auxiliaryModel_->rowScale_;
    3739         inverseScale = scale + numberRows_;
    3740         for (i=0;i<numberRows_;i++) {
     3861        inverseScale = scale + numberRows;
     3862        for (i=0;i<numberRows;i++) {
    37413863          rowActivity_[i] = auxiliaryModel_->rowActivityWork_[i]*inverseScale[i];
    37423864        }
    37433865      } else {
    3744         for (i=0;i<numberColumns_;i++) {
     3866        for (i=0;i<numberColumns;i++) {
    37453867          columnActivity_[i] = auxiliaryModel_->columnActivityWork_[i];
    37463868          reducedCost_[i] = auxiliaryModel_->reducedCostWork_[i];
    37473869        }
    3748         for (i=0;i<numberRows_;i++) {
     3870        for (i=0;i<numberRows;i++) {
    37493871          rowActivity_[i] = auxiliaryModel_->rowActivityWork_[i];
    37503872        }
     
    37523874      if (optimizationDirection_!=1.0) {
    37533875        // and modify reduced costs
    3754         for (i=0;i<numberColumns_;i++)
     3876        for (i=0;i<numberColumns;i++)
    37553877          reducedCost_[i] *= optimizationDirection_;
    37563878      }
     
    37593881      if (auxiliaryModel_->rowScale_) {
    37603882        const double * scale = auxiliaryModel_->columnScale_;
    3761         const double * inverseScale = scale + numberColumns_;
    3762         for (i=0;i<numberColumns_;i++) {
     3883        const double * inverseScale = scale + numberColumns;
     3884        for (i=0;i<numberColumns;i++) {
    37633885          double lower = auxiliaryModel_->columnLowerWork_[i];
    37643886          double upper = auxiliaryModel_->columnUpperWork_[i];
     
    37703892        }
    37713893        scale = auxiliaryModel_->rowScale_;
    3772         inverseScale = scale + numberRows_;
    3773         for (i=0;i<numberRows_;i++) {
     3894        inverseScale = scale + numberRows;
     3895        for (i=0;i<numberRows;i++) {
    37743896          double lower = auxiliaryModel_->rowLowerWork_[i];
    37753897          double upper = auxiliaryModel_->rowUpperWork_[i];
     
    37813903        }
    37823904      } else {
    3783         for (i=0;i<numberColumns_;i++) {
     3905        for (i=0;i<numberColumns;i++) {
    37843906          double lower = auxiliaryModel_->columnLowerWork_[i];
    37853907          double upper = auxiliaryModel_->columnUpperWork_[i];
     
    37903912          columnActivity_[i] = value;
    37913913        }
    3792         for (i=0;i<numberRows_;i++) {
     3914        for (i=0;i<numberRows;i++) {
    37933915          double lower = auxiliaryModel_->rowLowerWork_[i];
    37943916          double upper = auxiliaryModel_->rowUpperWork_[i];
     
    65006622                    case atLowerBound:
    65016623                      if (djValue<-dualTolerance)
    6502                         { thisWay =1;}
     6624                        { thisWay =1; }
    65036625                      else if (djValue>dualTolerance)
    65046626                        thisWay = 3;
     
    66366758    if (sequenceOut_<0) {
    66376759      if (directionIn_<0) {
    6638         assert (fabs(solution_[sequenceIn_]-upperIn_)<1.0e-7);
     6760        assert (fabs(solution_[sequenceIn_]-upperIn_)<5.0e-7);
    66396761        solution_[sequenceIn_]=upperIn_;
    66406762      } else {
    6641         assert (fabs(solution_[sequenceIn_]-lowerIn_)<1.0e-7);
     6763        assert (fabs(solution_[sequenceIn_]-lowerIn_)<5.0e-7);
    66426764        solution_[sequenceIn_]=lowerIn_;
    66436765      }
    66446766    } else {
    66456767      if (directionOut_<0) {
    6646         assert (fabs(solution_[sequenceOut_]-upperOut_)<1.0e-7);
     6768        assert (fabs(solution_[sequenceOut_]-upperOut_)<5.0e-7);
    66476769        solution_[sequenceOut_]=upperOut_;
    66486770      } else {
    6649         assert (fabs(solution_[sequenceOut_]-lowerOut_)<1.0e-7);
     6771        assert (fabs(solution_[sequenceOut_]-lowerOut_)<5.0e-7);
    66506772        solution_[sequenceOut_]=lowerOut_;
    66516773      }
     
    68346956  // sanity check
    68356957  // bad if empty (trap here to avoid using bad matrix_)
     6958#if 0
     6959  // but also check bounds
     6960  {
     6961    int badProblem = 0;
     6962    int i;
     6963    for (i=0;i<numberColumns_;i++) {
     6964      if (columnLower_[i]>columnUpper_[i])
     6965        badProblem++;
     6966    }
     6967    for (i=0;i<numberRows_;i++) {
     6968      if (rowLower_[i]>rowUpper_[i])
     6969        badProblem++;
     6970    }
     6971    if (badProblem) {
     6972      numberDualInfeasibilities_=0;
     6973      sumDualInfeasibilities_=0.0;
     6974      numberPrimalInfeasibilities_=badProblem;
     6975      sumPrimalInfeasibilities_=badProblem;
     6976      secondaryStatus_=6; // so user can see something odd
     6977      problemStatus_=1;
     6978      bool printIt = (specialOptions_&32768)==0 ? true : false; // no message if from Osi
     6979      if (printIt)
     6980        handler_->message(CLP_INFEASIBLE,messages_)
     6981          <<CoinMessageEol;
     6982      return 2;
     6983    }
     6984  }
     6985#endif
    68366986  if (!matrix_||(!matrix_->getNumElements()&&objective_->type()<2)) {
    68376987    int infeasNumber[2];
     
    70957245  otherModel.largestPrimalError_ = largestPrimalError_;
    70967246  otherModel.largestDualError_ = largestDualError_;
    7097   otherModel.largestSolutionError_ = largestSolutionError_;
     7247  otherModel.alphaAccuracy_ = alphaAccuracy_;
    70987248  otherModel.alpha_ = alpha_;
    70997249  otherModel.theta_ = theta_;
     
    89399089  return objectiveValue;
    89409090}
    8941 void
    8942 ClpSimplex::setSpecialOptions(unsigned int value)
    8943 {
    8944   specialOptions_=value;
    8945 }
    89469091// If user left factorization frequency then compute
    89479092void
  • branches/devel/Clp/src/ClpSimplex.hpp

    r918 r989  
    640640  //@}
    641641  /**@name most useful gets and sets */
    642   //@{
    643   // On reflection I doubt whether anyone uses so test
    644 private:
    645   /// Worst column primal infeasibility
    646   inline double columnPrimalInfeasibility() const
    647           { return columnPrimalInfeasibility_;} ;
    648   /// Sequence of worst (-1 if feasible)
    649   inline int columnPrimalSequence() const
    650           { return columnPrimalSequence_;} ;
    651   /// Worst row primal infeasibility
    652   inline double rowPrimalInfeasibility() const
    653           { return rowPrimalInfeasibility_;} ;
    654   /// Sequence of worst (-1 if feasible)
    655   inline int rowPrimalSequence() const
    656           { return rowPrimalSequence_;} ;
    657   /** Worst column dual infeasibility (note - these may not be as meaningful
    658       if the problem is primal infeasible */
    659   inline double columnDualInfeasibility() const
    660           { return columnDualInfeasibility_;} ;
    661   /// Sequence of worst (-1 if feasible)
    662   inline int columnDualSequence() const
    663           { return columnDualSequence_;} ;
    664   /// Worst row dual infeasibility
    665   inline double rowDualInfeasibility() const
    666           { return rowDualInfeasibility_;} ;
    667   /// Sequence of worst (-1 if feasible)
    668   inline int rowDualSequence() const
    669           { return rowDualSequence_;} ;
    670   /// Primal tolerance needed to make dual feasible (<largeTolerance)
    671   inline double primalToleranceToGetOptimal() const
    672           { return primalToleranceToGetOptimal_;} ;
    673   /// Remaining largest dual infeasibility
    674   inline double remainingDualInfeasibility() const
    675           { return remainingDualInfeasibility_;} ;
    676   /// Largest difference between input primal solution and computed
    677   inline double largestSolutionError() const
    678           { return largestSolutionError_;} ;
     642  //@{
     643public:
     644  /// Initial value for alpha accuracy calculation (-1.0 off)
     645  inline double alphaAccuracy() const
     646          { return alphaAccuracy_;} ;
     647  inline void setAlphaAccuracy(double value)
     648          { alphaAccuracy_ = value;} ;
    679649public:
    680650  /// Disaster handler
     
    958928  /// Create C++ lines to get to current state
    959929  void generateCpp( FILE * fp,bool defaultFactor=false);
    960   /** For advanced options
    961       1 - Don't keep changing infeasibility weight
    962       2 - Keep nonLinearCost round solves
    963       4 - Force outgoing variables to exact bound (primal)
    964       8 - Safe to use dense initial factorization
    965       16 -Just use basic variables for operation if column generation
    966       32 -Clean up with primal before strong branching
    967       64 -Treat problem as feasible until last minute (i.e. minimize infeasibilities)
    968       128 - Switch off all matrix sanity checks
    969       256 - No row copy
    970       512 - If not in values pass, solution guaranteed, skip as much as possible
    971       1024 - In branch and bound
    972       2048 - Don't bother to re-factorize if < 20 iterations
    973       4096 - Skip some optimality checks
    974       8192 - Do Primal when cleaning up primal
    975       16384 - In fast dual (so we can switch off things)
    976       32678 - called from Osi
    977       65356 - keep arrays around as much as possible
    978       NOTE - many applications can call Clp but there may be some short cuts
    979              which are taken which are not guaranteed safe from all applications.
    980              Vetted applications will have a bit set and the code may test this
    981              At present I expect a few such applications - if too many I will
    982              have to re-think.  It is up to application owner to change the code
    983              if she/he needs these short cuts.  I will not debug unless in Coin
    984              repository.  See COIN_CLP_VETTED comments.
    985       0x01000000 is Cbc (and in branch and bound)
    986       0x02000000 is in a different branch and bound
    987   */
    988 #define COIN_CBC_USING_CLP 0x01000000
    989   inline unsigned int specialOptions() const
    990   { return specialOptions_;};
    991   void setSpecialOptions(unsigned int value);
    992930  /// Gets clean and emptyish factorization
    993931  ClpFactorization * getEmptyFactorization();
     
    11361074  /// Largest error on basic duals
    11371075  double largestDualError_;
    1138   /// Largest difference between input primal solution and computed
    1139   double largestSolutionError_;
     1076  /// For computing whether to re-factorize
     1077  double alphaAccuracy_;
    11401078  /// Dual bound
    11411079  double dualBound_;
     
    12681206  */
    12691207  ClpNonLinearCost * nonLinearCost_;
    1270   /** For advanced options
    1271       See get and set for meaning
    1272   */
    1273   unsigned int specialOptions_;
    12741208  /// So we know when to be cautious
    12751209  int lastBadIteration_;
  • branches/devel/Clp/src/ClpSimplexDual.cpp

    r920 r989  
    484484    CoinMemcpyN(dual_,numberRows_,saveDuals);
    485485  }
     486  if (alphaAccuracy_!=-1.0)
     487    alphaAccuracy_ = 1.0;
    486488  int returnCode = startupSolve(ifValuesPass,saveDuals,startFinishOptions);
    487489  // Save so can see if doing after primal
     
    33003302  return status;
    33013303}
     3304//static int count_alpha=0;
    33023305/* Checks if finished.  Updates status */
    33033306void
     
    33233326  double realDualInfeasibilities=0.0;
    33243327  if (type==2) {
     3328    if (alphaAccuracy_!=-1.0)
     3329      alphaAccuracy_=-2.0;
    33253330    // trouble - restore solution
    33263331    CoinMemcpyN(saveStatus_,numberColumns_+numberRows_,status_);
     
    33403345  double changeCost;
    33413346  bool unflagVariables = true;
    3342   if (problemStatus_>-3||factorization_->pivots()) {
    3343     // factorize
    3344     // later on we will need to recover from singularities
    3345     // also we could skip if first time
    3346     // save dual weights
    3347     dualRowPivot_->saveWeights(this,1);
    3348     if (type) {
    3349       // is factorization okay?
    3350       if (internalFactorize(1)) {
    3351         // no - restore previous basis
    3352         unflagVariables = false;
    3353         assert (type==1);
    3354         changeMade_++; // say something changed
    3355         // Keep any flagged variables
    3356         int i;
    3357         for (i=0;i<numberRows_+numberColumns_;i++) {
    3358           if (flagged(i))
    3359             saveStatus_[i] |= 64; //say flagged
    3360         }
    3361         CoinMemcpyN(saveStatus_,numberColumns_+numberRows_,status_);
    3362         CoinMemcpyN(savedSolution_+numberColumns_ ,
    3363                numberRows_,rowActivityWork_);
    3364         CoinMemcpyN(savedSolution_ ,
    3365                numberColumns_,columnActivityWork_);
    3366         // restore extra stuff
    3367         int dummy;
    3368         matrix_->generalExpanded(this,6,dummy);
    3369         // get correct bounds on all variables
    3370         resetFakeBounds();
    3371         // need to reject something
    3372         char x = isColumn(sequenceOut_) ? 'C' :'R';
    3373         handler_->message(CLP_SIMPLEX_FLAG,messages_)
    3374           <<x<<sequenceWithin(sequenceOut_)
    3375           <<CoinMessageEol;
    3376         setFlagged(sequenceOut_);
    3377         progress_->clearBadTimes();
    3378        
    3379         // Go to safe
    3380         factorization_->pivotTolerance(0.99);
    3381         forceFactorization_=1; // a bit drastic but ..
    3382         type = 2;
    3383         //assert (internalFactorize(1)==0);
     3347  bool weightsSaved=false;
     3348  if (alphaAccuracy_<0.0||!numberPivots||alphaAccuracy_>1.0e4||factorization_->pivots()>20) {
     3349    if (problemStatus_>-3||numberPivots) {
     3350      // factorize
     3351      // later on we will need to recover from singularities
     3352      // also we could skip if first time
     3353      // save dual weights
     3354      dualRowPivot_->saveWeights(this,1);
     3355      weightsSaved=true;
     3356      if (type) {
     3357        // is factorization okay?
    33843358        if (internalFactorize(1)) {
     3359          // no - restore previous basis
     3360          unflagVariables = false;
     3361          assert (type==1);
     3362          changeMade_++; // say something changed
     3363          // Keep any flagged variables
     3364          int i;
     3365          for (i=0;i<numberRows_+numberColumns_;i++) {
     3366            if (flagged(i))
     3367              saveStatus_[i] |= 64; //say flagged
     3368          }
    33853369          CoinMemcpyN(saveStatus_,numberColumns_+numberRows_,status_);
    33863370          CoinMemcpyN(savedSolution_+numberColumns_ ,
    3387                  numberRows_,rowActivityWork_);
     3371                      numberRows_,rowActivityWork_);
    33883372          CoinMemcpyN(savedSolution_ ,
    3389                  numberColumns_,columnActivityWork_);
     3373                      numberColumns_,columnActivityWork_);
    33903374          // restore extra stuff
    33913375          int dummy;
    33923376          matrix_->generalExpanded(this,6,dummy);
    3393           // debug
    3394           int returnCode = internalFactorize(1);
    3395           while (returnCode) {
    3396             // ouch
    3397             // switch off dense
    3398             int saveDense = factorization_->denseThreshold();
    3399             factorization_->setDenseThreshold(0);
    3400             // Go to safe
    3401             factorization_->pivotTolerance(0.99);
    3402             // make sure will do safe factorization
    3403             pivotVariable_[0]=-1;
    3404             returnCode=internalFactorize(2);
    3405             factorization_->setDenseThreshold(saveDense);
    3406           }
    34073377          // get correct bounds on all variables
    34083378          resetFakeBounds();
    3409         }
    3410       }
    3411     }
    3412     if (problemStatus_!=-4||factorization_->pivots()>10)
    3413       problemStatus_=-3;
     3379          // need to reject something
     3380          char x = isColumn(sequenceOut_) ? 'C' :'R';
     3381          handler_->message(CLP_SIMPLEX_FLAG,messages_)
     3382            <<x<<sequenceWithin(sequenceOut_)
     3383            <<CoinMessageEol;
     3384          setFlagged(sequenceOut_);
     3385          progress_->clearBadTimes();
     3386         
     3387          // Go to safe
     3388          factorization_->pivotTolerance(0.99);
     3389          forceFactorization_=1; // a bit drastic but ..
     3390          type = 2;
     3391          //assert (internalFactorize(1)==0);
     3392          if (internalFactorize(1)) {
     3393            CoinMemcpyN(saveStatus_,numberColumns_+numberRows_,status_);
     3394            CoinMemcpyN(savedSolution_+numberColumns_ ,
     3395                        numberRows_,rowActivityWork_);
     3396            CoinMemcpyN(savedSolution_ ,
     3397                        numberColumns_,columnActivityWork_);
     3398            // restore extra stuff
     3399            int dummy;
     3400            matrix_->generalExpanded(this,6,dummy);
     3401            // debug
     3402            int returnCode = internalFactorize(1);
     3403            while (returnCode) {
     3404              // ouch
     3405              // switch off dense
     3406              int saveDense = factorization_->denseThreshold();
     3407              factorization_->setDenseThreshold(0);
     3408              // Go to safe
     3409              factorization_->pivotTolerance(0.99);
     3410              // make sure will do safe factorization
     3411              pivotVariable_[0]=-1;
     3412              returnCode=internalFactorize(2);
     3413              factorization_->setDenseThreshold(saveDense);
     3414            }
     3415            // get correct bounds on all variables
     3416            resetFakeBounds();
     3417          }
     3418        }
     3419      }
     3420      if (problemStatus_!=-4||numberPivots>10)
     3421        problemStatus_=-3;
     3422    }
     3423  } else {
     3424    //printf("testing with accuracy of %g and status of %d\n",alphaAccuracy_,problemStatus_);
     3425    //count_alpha++;
     3426    //if ((count_alpha%5000)==0)
     3427    //printf("count alpha %d\n",count_alpha);
    34143428  }
    34153429  // at this stage status is -3 or -4 if looks infeasible
     
    34493463    factorization_->pivotTolerance(newTolerance);
    34503464    forceFactorization_=1; // a bit drastic but ..
     3465    if (alphaAccuracy_!=-1.0)
     3466      alphaAccuracy_=-2.0;
    34513467    type = 2;
    34523468    //assert (internalFactorize(1)==0);
     
    35183534        if(factorization_->pivotTolerance()<0.2)
    35193535          factorization_->pivotTolerance(0.2);
     3536        if (alphaAccuracy_!=-1.0)
     3537          alphaAccuracy_=-2.0;
    35203538        if (internalFactorize(1)) {
    35213539          CoinMemcpyN(saveStatus_,numberColumns_+numberRows_,status_);
     
    39944012    matrix_->generalExpanded(this,5,dummy);
    39954013  }
    3996 
    3997   // restore weights (if saved) - also recompute infeasibility list
    3998   if (tentativeStatus>-3)
    3999     dualRowPivot_->saveWeights(this,(type <2) ? 2 : 4);
    4000   else
    4001     dualRowPivot_->saveWeights(this,3);
     4014  if (weightsSaved) {
     4015    // restore weights (if saved) - also recompute infeasibility list
     4016    if (tentativeStatus>-3)
     4017      dualRowPivot_->saveWeights(this,(type <2) ? 2 : 4);
     4018    else
     4019      dualRowPivot_->saveWeights(this,3);
     4020  }
    40024021  // unflag all variables (we may want to wait a bit?)
    40034022  if ((tentativeStatus!=-2&&tentativeStatus!=-1)&&unflagVariables) {
     
    41424161  int fake=-999; // signal sort
    41434162  matrix_->correctSequence(this,fake,fake);
     4163  if (alphaAccuracy_>0.0)
     4164      alphaAccuracy_=1.0;
    41444165}
    41454166/* While updateDualsInDual sees what effect is of flip
     
    49084929  double saveDualBound = dualBound_;
    49094930
     4931  if (alphaAccuracy_!=-1.0)
     4932    alphaAccuracy_ = 1.0;
    49104933  double objectiveChange;
    49114934  // for dual we will change bounds using dualBound_
  • branches/devel/Clp/src/ClpSimplexOther.cpp

    r911 r989  
    888888  if (numberExtraRows) {
    889889    CoinPackedMatrix newCopy;
     890    newCopy.setExtraGap(0.0);
     891    newCopy.setExtraMajor(0.0);
    890892    newCopy.submatrixOfWithDuplicates(rowCopy,kRow,which);
    891893    rowCopy=newCopy;
     
    12261228    small->setPerturbation(perturbation_);
    12271229    small->defaultFactorizationFrequency();
     1230    small->setAlphaAccuracy(alphaAccuracy_);
    12281231    // If no rows left then no tightening!
    12291232    if (!numberRows2||!numberColumns2)
  • branches/devel/Clp/src/ClpSimplexPrimal.cpp

    r889 r989  
    191191  int initialNegDjs=-1;
    192192  // initialize - maybe values pass and algorithm_ is +1
     193#if 0
     194  // if so - put in any superbasic costed slacks
     195  if (ifValuesPass&&specialOptions_<0x01000000) {
     196    // Get column copy
     197    const CoinPackedMatrix * columnCopy = matrix();
     198    const int * row = columnCopy->getIndices();
     199    const CoinBigIndex * columnStart = columnCopy->getVectorStarts();
     200    const int * columnLength = columnCopy->getVectorLengths();
     201    //const double * element = columnCopy->getElements();
     202    int n=0;
     203    for (int iColumn = 0;iColumn<numberColumns_;iColumn++) {
     204      if (columnLength[iColumn]==1) {
     205        Status status = getColumnStatus(iColumn);
     206        if (status!=basic&&status!=isFree) {
     207          double value = columnActivity_[iColumn];
     208          if (fabs(value-columnLower_[iColumn])>primalTolerance_&&
     209              fabs(value-columnUpper_[iColumn])>primalTolerance_) {
     210            int iRow = row[columnStart[iColumn]];
     211            if (getRowStatus(iRow)==basic) {
     212              setRowStatus(iRow,superBasic);
     213              setColumnStatus(iColumn,basic);
     214              n++;
     215            }
     216          }
     217        }   
     218      }
     219    }
     220    printf("%d costed slacks put in basis\n",n);
     221  }
     222#endif
    193223  if (!startup(ifValuesPass,startFinishOptions)) {
    194224   
     
    435465  //printf("XXXXY final cost %g\n",infeasibilityCost_);
    436466  progress_->initialWeight_=0.0;
    437   if (problemStatus_==1) {
     467  if (problemStatus_==1&&secondaryStatus_!=6) {
    438468    infeasibilityCost_=0.0;
    439469    createRim(1+4);
  • branches/devel/Clp/src/ClpSolve.cpp

    r905 r989  
    15991599      if (interrupt)
    16001600        currentModel = &small;
     1601      small.defaultFactorizationFrequency();
    16011602      if (dynamic_cast< ClpPackedMatrix*>(matrix_)) {
    16021603        // See if original wanted vector
     
    21912192      saveLower=NULL;
    21922193      saveUpper=NULL;
    2193       model2->primal(1);
     2194      if (method!=ClpSolve::useBarrierNoCross)
     2195        model2->primal(1);
    21942196    }
    21952197    model2->setPerturbation(savePerturbation);
  • branches/devel/Clp/src/Idiot.cpp

    r987 r989  
    461461  dj=new double[ncols];
    462462  delete [] whenUsed_;
     463  bool oddSlacks=false;
    463464  // See if any costed slacks
    464465  int numberSlacks=0;
     
    467468      numberSlacks++;
    468469  }
    469   if (!numberSlacks||true) {
     470  if (!numberSlacks) {
    470471    whenUsed_=new int[ncols];
    471472  } else {
    472473    printf("%d slacks\n",numberSlacks);
    473     strategy_ |= 16384;
     474    oddSlacks=true;
    474475    int extra = (int) (nrows*sizeof(double)/sizeof(int));
    475476    whenUsed_=new int[2*ncols+2*nrows+extra];
     
    10441045  }
    10451046  muAtExit_=mu;
     1047  // For last iteration make as feasible as possible
     1048  if (oddSlacks);
     1049    strategy_ |= 16384;
    10461050  // not scaled
    10471051  n = cleanIteration(iteration, ordStart,ordEnd,
     
    10491053                     model_->rowLower(), model_->rowUpper(),
    10501054                     cost, element, fixTolerance,lastResult.objval,lastResult.infeas);
    1051   if ((strategy_&16384)!=0) {
    1052     int * posSlack = whenUsed_+ncols;
    1053     int * negSlack = posSlack+nrows;
    1054     int * nextSlack = negSlack + nrows;
    1055     double * rowsol2 = (double *) (nextSlack+ncols);
    1056     for (i=0;i<nrows;i++)
    1057       rowsol[i] += rowsol2[i];
    1058   }
    1059   if ((logLevel&1)==0) {
     1055  if ((logLevel&1)==0||(strategy_&16384)!=0) {
    10601056    printf(
    10611057            "%d - mu %g, infeasibility %g, objective %g, %d interior\n",
     
    12841280    }
    12851281    model_->setColumnStatus(i,ClpSimplex::superBasic);
     1282  }
     1283  if ((strategy_&16384)!=0) {
     1284    // put in basis
     1285    int * posSlack = whenUsed_+ncols;
     1286    int * negSlack = posSlack+nrows;
     1287    int * nextSlack = negSlack + nrows;
     1288    for (i=0;i<nrows;i++) {
     1289      int n=0;
     1290      int iCol;
     1291      iCol =posSlack[i];
     1292      if (iCol>=0) {
     1293        if (colsol[iCol]>lower[iCol]+1.0e-8&&
     1294            colsol[iCol]<upper[iCol]+1.0e-8) {
     1295          model_->setColumnStatus(iCol,ClpSimplex::basic);
     1296          n++;
     1297        }
     1298        while (nextSlack[iCol]>=0) {
     1299          iCol = nextSlack[iCol];
     1300          if (colsol[iCol]>lower[iCol]+1.0e-8&&
     1301              colsol[iCol]<upper[iCol]+1.0e-8) {
     1302            model_->setColumnStatus(iCol,ClpSimplex::basic);
     1303            n++;
     1304          }
     1305        }
     1306      }
     1307      iCol =negSlack[i];
     1308      if (iCol>=0) {
     1309        if (colsol[iCol]>lower[iCol]+1.0e-8&&
     1310            colsol[iCol]<upper[iCol]+1.0e-8) {
     1311          model_->setColumnStatus(iCol,ClpSimplex::basic);
     1312          n++;
     1313        }
     1314        while (nextSlack[iCol]>=0) {
     1315          iCol = nextSlack[iCol];
     1316          if (colsol[iCol]>lower[iCol]+1.0e-8&&
     1317              colsol[iCol]<upper[iCol]+1.0e-8) {
     1318            model_->setColumnStatus(iCol,ClpSimplex::basic);
     1319            n++;
     1320          }
     1321        }
     1322      }
     1323      if (n) {
     1324        if (fabs(rowsol[i]-rowlower[i])<fabs(rowsol[i]-rowupper[i]))
     1325          model_->setRowStatus(i,ClpSimplex::atLowerBound);
     1326        else
     1327          model_->setRowStatus(i,ClpSimplex::atUpperBound);
     1328        if (n>1)
     1329          printf("%d basic on row %d!\n",n,i);
     1330      }
     1331    }
    12861332  }
    12871333  double maxmin;
Note: See TracChangeset for help on using the changeset viewer.