Changeset 1366 for trunk/Clp


Ignore:
Timestamp:
May 16, 2009 4:23:18 AM (11 years ago)
Author:
forrest
Message:

try and improve quadratic and barrier

Location:
trunk/Clp/src
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • trunk/Clp/src/ClpMain.cpp

    r1344 r1366  
    2020#include "CoinSort.hpp"
    2121// History since 1.0 at end
    22 #define CLPVERSION "1.09.00"
     22#define CLPVERSION "1.10.00"
    2323
    2424#include "CoinMpsIO.hpp"
     
    26052605          //printf("possible %d - but no integers\n",numberPossible);
    26062606        }
    2607         if (numberPossible>1&&(numberIntegers||true)) {
     2607        if (numberPossible>1&&(numberIntegers||false)) {
    26082608          //
    26092609          printf("possible %d - %d integers\n",numberPossible,numberIntegers);
  • trunk/Clp/src/ClpPredictorCorrector.cpp

    r1321 r1366  
    195195  int lastGoodIteration=0;
    196196  double bestObjectiveGap=COIN_DBL_MAX;
     197  double bestObjective=COIN_DBL_MAX;
    197198  int saveIteration=-1;
    198199  bool sloppyOptimal=false;
     
    248249    historyInfeasibility_[LENGTH_HISTORY-1]=complementarityGap_;
    249250    // switch off saved if changes
    250     if (saveIteration+10<numberIterations_&&
    251         complementarityGap_*2.0<historyInfeasibility_[0])
    252       saveIteration=-1;
     251    //if (saveIteration+10<numberIterations_&&
     252    //complementarityGap_*2.0<historyInfeasibility_[0])
     253    //saveIteration=-1;
    253254    lastStep = CoinMin(actualPrimalStep_,actualDualStep_);
    254255    double goodGapChange;
     256    //#define KEEP_GOING_IF_FIXED 0
     257#ifndef KEEP_GOING_IF_FIXED
     258#define KEEP_GOING_IF_FIXED 10000
     259#endif
    255260    if (!sloppyOptimal) {
    256261      goodGapChange=0.93;
    257262    } else {
    258263      goodGapChange=0.7;
     264      if (numberFixed>KEEP_GOING_IF_FIXED)
     265        goodGapChange=0.99; // make more likely to carry on
    259266    }
    260267    double gapO;
     
    275282        <<CoinMessageEol;
    276283      //start saving best
    277       if (gapO<bestObjectiveGap) {
     284      if (gapO<bestObjectiveGap)
     285        bestObjectiveGap=gapO;
     286      if (primalObjective_<bestObjective) {
    278287        saveIteration=numberIterations_;
    279         bestObjectiveGap=gapO;
     288        bestObjective=primalObjective_;
    280289        if (!savePi) {
    281290          savePi=new double[numberRows_];
     
    333342        }
    334343      }
     344      if (complementarityGap_>0.5*checkGap&&primalObjective_>
     345          bestObjective+1.0e-9&&
     346          (numberIterations_>saveIteration+5||numberIterations_>100)) {
     347        handler_->message(CLP_BARRIER_EXIT2,messages_)
     348          <<saveIteration
     349          <<CoinMessageEol;
     350        break;
     351      }
    335352    }
    336353    if ((gapO<1.0e-6||(gapO<1.0e-4&&complementarityGap_<0.1))&&!sloppyOptimal) {
     
    365382        <<complementarityGap_<<"not decreasing"
    366383        <<CoinMessageEol;
    367       if (gapO>0.75*lastGood) {
     384      if (gapO>0.75*lastGood&&numberFixed<KEEP_GOING_IF_FIXED) {
    368385        break;
    369386      }
     
    831848  delete [] saveSU;
    832849  if (savePi) {
    833     //std::cout<<"Restoring from iteration "<<saveIteration<<std::endl;
     850    std::cout<<"Restoring from iteration "<<saveIteration<<std::endl;
    834851    CoinMemcpyN(savePi,numberRows_,dual_);
    835852    CoinMemcpyN(savePrimal,numberTotal,solution_);
     
    34323449      } else {
    34333450        numberKilled++;
     3451        if (solution_[iColumn]!=lower_[iColumn]&&
     3452            solution_[iColumn]!=upper_[iColumn]) {
     3453          printf("%d %g %g %g\n",iColumn,lower_[iColumn],
     3454                 solution_[iColumn],upper_[iColumn]);
     3455        }
    34343456        diagonal_[iColumn]=0.0;
    34353457        zVec_[iColumn]=0.0;
     
    35653587    solutionNorm_=rhsNorm_;
    35663588  }
    3567   double scaledRHSError=maximumRHSError_/solutionNorm_;
     3589  double scaledRHSError=maximumRHSError_/(solutionNorm_+10.0);
    35683590  bool dualFeasible=true;
     3591#if KEEP_GOING_IF_FIXED > 5
    35693592  if (maximumBoundInfeasibility_>primalTolerance()||
    35703593      scaledRHSError>primalTolerance())
    35713594    primalFeasible=false;
     3595#else
     3596  if (maximumBoundInfeasibility_>primalTolerance()||
     3597      scaledRHSError>CoinMax(CoinMin(100.0*primalTolerance(),1.0e-5),
     3598                             primalTolerance())
     3599    primalFeasible=false;
     3600#endif
    35723601  // relax dual test if obj big and gap smallish
    35733602  double gap=fabs(primalObjective_-dualObjective_);
  • trunk/Clp/src/ClpSimplexNonlinear.cpp

    r1321 r1366  
    655655                         array1);
    656656    if (result) {
     657      if (result==2) {
     658        // does not look good
     659        double currentObj;
     660        double thetaObj;
     661        double predictedObj;
     662        objective_->stepLength(this,solution_,solution_,0.0,
     663                               currentObj,thetaObj,predictedObj);
     664        if (currentObj==predictedObj) {
     665#ifdef CLP_INVESTIGATE
     666        printf("looks bad - no change in obj %g\n",currentObj);
     667#endif
     668          if (factorization_->pivots())
     669            result=3;
     670          else
     671            problemStatus_=0;
     672        }
     673      }
    657674      if (result==3)
    658675        break; // null vector not accurate
     
    663680        double predictedObj;
    664681        objective_->stepLength(this,solution_,solution_,0.0,
    665                                currentObj,predictedObj,thetaObj);
     682                               currentObj,thetaObj,predictedObj);
    666683        printf("obj %g after interior move\n",currentObj);
    667684      }
     
    20622079      printf("out of loop after %d passes\n",nPasses);
    20632080#endif
     2081    if (nPasses>=1000)
     2082      returnCode=2;
    20642083    bool ordinaryDj=false;
    20652084    //if(sequenceIn_>=0&&numberNonBasic==1&&theta_<1.0e-7&&theta_==basicTheta)
     
    22712290          } else {
    22722291            returnCode=2; // do single incoming
    2273             //returnCode=1;
     2292            returnCode=1;
    22742293          }
    22752294        }
  • trunk/Clp/src/ClpSimplexOther.cpp

    r1344 r1366  
    6565    //columnArray_[0]->checkClear();
    6666    int iSequence = which[i];
     67    if (iSequence<0) {
     68      costIncreased[i] = 0.0;
     69      sequenceIncreased[i] = -1;
     70      costDecreased[i] = 0.0;
     71      sequenceDecreased[i] = -1;
     72      continue;
     73    }
    6774    double costIncrease=COIN_DBL_MAX;
    6875    double costDecrease=COIN_DBL_MAX;
     
    196203    }
    197204  }
    198   //rowArray_[0]->clear();
     205  rowArray_[0]->clear();
    199206  //rowArray_[1]->clear();
    200207  //columnArray_[1]->clear();
    201   //columnArray_[0]->clear();
     208  columnArray_[0]->clear();
    202209  //rowArray_[3]->clear();
    203210  if (!optimizationDirection_)
     
    31303137      int nRow=0;
    31313138      double obj=0.0;
    3132       CoinZeroN(rowActivity,nRow);
     3139      CoinZeroN(rowActivity,numberRows_);
    31333140      for (iColumn=0;iColumn<numJ;iColumn++) {
    31343141        int iValue = stack[iColumn];
  • trunk/Clp/src/ClpSolve.cpp

    r1344 r1366  
    18751875    bool scale=false;
    18761876    bool doKKT=false;
     1877    bool forceFixing=false;
    18771878    if (barrierOptions&16) {
    18781879      barrierOptions &= ~16;
     
    18861887      barrierOptions &= ~256;
    18871888      presolveInCrossover=true;
     1889    }
     1890    if (barrierOptions&512) {
     1891      barrierOptions &= ~512;
     1892      forceFixing=true;
    18881893    }
    18891894    if (barrierOptions&8) {
     
    20532058      CoinMemcpyN(barrier.rowUpper(),numberRows,saveUpper+numberColumns);
    20542059    }
    2055     if (numberFixed*20>barrier.numberRows()&&numberFixed>5000&&
     2060    if (((numberFixed*20>barrier.numberRows()&&numberFixed>5000)||forceFixing)&&
    20562061        presolveInCrossover) {
    20572062      // may as well do presolve
    2058       barrier.fixFixed();
     2063      if (!forceFixing) {
     2064        barrier.fixFixed();
     2065      } else {
     2066        // Fix
     2067        int n=barrier.numberColumns();
     2068        double * lower = barrier.columnLower();
     2069        double * upper = barrier.columnUpper();
     2070        double * solution = barrier.primalColumnSolution();
     2071        int nFix=0;
     2072        for (int i=0;i<n;i++) {
     2073          if (barrier.fixedOrFree(i)&&lower[i]<upper[i]) {
     2074            double value = solution[i];
     2075            if (value<lower[i]+1.0e-6&&value-lower[i]<upper[i]-value) {
     2076              solution[i]=lower[i];
     2077              upper[i]=lower[i];
     2078              nFix++;
     2079            } else if (value>upper[i]-1.0e-6&&value-lower[i]>upper[i]-value) {
     2080              solution[i]=upper[i];
     2081              lower[i]=upper[i];
     2082              nFix++;
     2083            }
     2084          }
     2085        }
     2086        printf("%d columns fixed\n",nFix);
     2087        int nr=barrier.numberRows();
     2088        lower = barrier.rowLower();
     2089        upper = barrier.rowUpper();
     2090        solution = barrier.primalRowSolution();
     2091        nFix=0;
     2092        for (int i=0;i<nr;i++) {
     2093          if (barrier.fixedOrFree(i+n)&&lower[i]<upper[i]) {
     2094            double value = solution[i];
     2095            if (value<lower[i]+1.0e-6&&value-lower[i]<upper[i]-value) {
     2096              solution[i]=lower[i];
     2097              upper[i]=lower[i];
     2098              nFix++;
     2099            } else if (value>upper[i]-1.0e-6&&value-lower[i]>upper[i]-value) {
     2100              solution[i]=upper[i];
     2101              lower[i]=upper[i];
     2102              nFix++;
     2103            }
     2104          }
     2105        }
     2106        printf("%d row slacks fixed\n",nFix);
     2107      }
    20592108      saveModel2=model2;
    20602109      extraPresolve=true;
    20612110    } else if (numberFixed) {
    20622111      // Set fixed to bounds (may have restored earlier solution)
    2063       barrier.fixFixed(false);
     2112      if (!forceFixing) {
     2113        barrier.fixFixed(false);
     2114      } else {
     2115        // Fix
     2116        int n=barrier.numberColumns();
     2117        double * lower = barrier.columnLower();
     2118        double * upper = barrier.columnUpper();
     2119        double * solution = barrier.primalColumnSolution();
     2120        int nFix=0;
     2121        for (int i=0;i<n;i++) {
     2122          if (barrier.fixedOrFree(i)&&lower[i]<upper[i]) {
     2123            double value = solution[i];
     2124            if (value<lower[i]+1.0e-8&&value-lower[i]<upper[i]-value) {
     2125              solution[i]=lower[i];
     2126              upper[i]=lower[i];
     2127              nFix++;
     2128            } else if (value>upper[i]-1.0e-8&&value-lower[i]>upper[i]-value) {
     2129              solution[i]=upper[i];
     2130              lower[i]=upper[i];
     2131              nFix++;
     2132            } else {
     2133              //printf("fixcol %d %g <= %g <= %g\n",
     2134              //     i,lower[i],solution[i],upper[i]);
     2135            }
     2136          }
     2137        }
     2138        printf("%d columns fixed\n",nFix);
     2139        int nr=barrier.numberRows();
     2140        lower = barrier.rowLower();
     2141        upper = barrier.rowUpper();
     2142        solution = barrier.primalRowSolution();
     2143        nFix=0;
     2144        for (int i=0;i<nr;i++) {
     2145          if (barrier.fixedOrFree(i+n)&&lower[i]<upper[i]) {
     2146            double value = solution[i];
     2147            if (value<lower[i]+1.0e-5&&value-lower[i]<upper[i]-value) {
     2148              solution[i]=lower[i];
     2149              upper[i]=lower[i];
     2150              nFix++;
     2151            } else if (value>upper[i]-1.0e-5&&value-lower[i]>upper[i]-value) {
     2152              solution[i]=upper[i];
     2153              lower[i]=upper[i];
     2154              nFix++;
     2155            } else {
     2156              printf("fixrow %d %g <= %g <= %g\n",
     2157                     i,lower[i],solution[i],upper[i]);
     2158            }
     2159          }
     2160        }
     2161        printf("%d row slacks fixed\n",nFix);
     2162      }
    20642163    }
    20652164#ifdef BORROW   
     
    21212220        model2->createStatus();
    21222221        // solve
    2123         model2->setPerturbation(100);
     2222        if (!forceFixing)
     2223          model2->setPerturbation(100);
    21242224        if (model2->factorizationFrequency()==200) {
    21252225          // User did not touch preset
     
    21282228#if 1
    21292229        // throw some into basis
    2130         {
     2230        if(!forceFixing) {
    21312231          int numberRows = model2->numberRows();
    21322232          int numberColumns = model2->numberColumns();
     
    21692269          delete [] sort;
    21702270          delete [] dsort;
    2171         }
    2172         // model2->allSlackBasis();
    2173         if (gap<1.0e-3*static_cast<double> (numberRows+numberColumns)) {
    2174           if (saveUpper) {
    2175             int numberRows = model2->numberRows();
    2176             int numberColumns = model2->numberColumns();
    2177             CoinMemcpyN(saveLower,numberColumns,model2->columnLower());
    2178             CoinMemcpyN(saveLower+numberColumns,numberRows,model2->rowLower());
    2179             delete [] saveLower;
    2180             CoinMemcpyN(saveUpper,numberColumns,model2->columnUpper());
    2181             CoinMemcpyN(saveUpper+numberColumns,numberRows,model2->rowUpper());
    2182             delete [] saveUpper;
    2183             saveLower=NULL;
    2184             saveUpper=NULL;
    2185           }
    2186           int numberRows = model2->numberRows();
    2187           int numberColumns = model2->numberColumns();
    2188           // just primal values pass
    2189           double saveScale = model2->objectiveScale();
    2190           model2->setObjectiveScale(1.0e-3);
    2191           model2->primal(2);
    2192           model2->setObjectiveScale(saveScale);
    2193           // save primal solution and copy back dual
    2194           CoinMemcpyN(model2->primalRowSolution(),
    2195                       numberRows,rowPrimal);
    2196           CoinMemcpyN(rowDual,
    2197                       numberRows,model2->dualRowSolution());
    2198           CoinMemcpyN(model2->primalColumnSolution(),
    2199                       numberColumns,columnPrimal);
    2200           CoinMemcpyN(columnDual,
    2201                       numberColumns,model2->dualColumnSolution());
    2202           //model2->primal(1);
    2203           // clean up reduced costs and flag variables
    2204           {
    2205             double * dj = model2->dualColumnSolution();
    2206             double * cost = model2->objective();
    2207             double * saveCost = new double[numberColumns];
    2208             CoinMemcpyN(cost,numberColumns,saveCost);
    2209             double * saveLower = new double[numberColumns];
    2210             double * lower = model2->columnLower();
    2211             CoinMemcpyN(lower,numberColumns,saveLower);
    2212             double * saveUpper = new double[numberColumns];
    2213             double * upper = model2->columnUpper();
    2214             CoinMemcpyN(upper,numberColumns,saveUpper);
    2215             int i;
    2216             double tolerance = 10.0*dualTolerance_;
    2217             for ( i=0;i<numberColumns;i++) {
    2218               if (model2->getStatus(i)==basic) {
    2219                 dj[i]=0.0;
    2220               } else if (model2->getStatus(i)==atLowerBound) {
    2221                 if (optimizationDirection_*dj[i]<tolerance) {
    2222                   if (optimizationDirection_*dj[i]<0.0) {
    2223                     //if (dj[i]<-1.0e-3)
    2224                     //printf("bad dj at lb %d %g\n",i,dj[i]);
    2225                     cost[i] -= dj[i];
    2226                     dj[i]=0.0;
     2271          // model2->allSlackBasis();
     2272          if (gap<1.0e-3*static_cast<double> (numberRows+numberColumns)) {
     2273            if (saveUpper) {
     2274              int numberRows = model2->numberRows();
     2275              int numberColumns = model2->numberColumns();
     2276              CoinMemcpyN(saveLower,numberColumns,model2->columnLower());
     2277              CoinMemcpyN(saveLower+numberColumns,numberRows,model2->rowLower());
     2278              CoinMemcpyN(saveUpper,numberColumns,model2->columnUpper());
     2279              CoinMemcpyN(saveUpper+numberColumns,numberRows,model2->rowUpper());
     2280              delete [] saveLower;
     2281              delete [] saveUpper;
     2282              saveLower=NULL;
     2283              saveUpper=NULL;
     2284            }
     2285            int numberRows = model2->numberRows();
     2286            int numberColumns = model2->numberColumns();
     2287            // just primal values pass
     2288            double saveScale = model2->objectiveScale();
     2289            model2->setObjectiveScale(1.0e-3);
     2290            model2->primal(2);
     2291            model2->setObjectiveScale(saveScale);
     2292            // save primal solution and copy back dual
     2293            CoinMemcpyN(model2->primalRowSolution(),
     2294                        numberRows,rowPrimal);
     2295            CoinMemcpyN(rowDual,
     2296                        numberRows,model2->dualRowSolution());
     2297            CoinMemcpyN(model2->primalColumnSolution(),
     2298                        numberColumns,columnPrimal);
     2299            CoinMemcpyN(columnDual,
     2300                        numberColumns,model2->dualColumnSolution());
     2301            //model2->primal(1);
     2302            // clean up reduced costs and flag variables
     2303            {
     2304              double * dj = model2->dualColumnSolution();
     2305              double * cost = model2->objective();
     2306              double * saveCost = new double[numberColumns];
     2307              CoinMemcpyN(cost,numberColumns,saveCost);
     2308              double * saveLower = new double[numberColumns];
     2309              double * lower = model2->columnLower();
     2310              CoinMemcpyN(lower,numberColumns,saveLower);
     2311              double * saveUpper = new double[numberColumns];
     2312              double * upper = model2->columnUpper();
     2313              CoinMemcpyN(upper,numberColumns,saveUpper);
     2314              int i;
     2315              double tolerance = 10.0*dualTolerance_;
     2316              for ( i=0;i<numberColumns;i++) {
     2317                if (model2->getStatus(i)==basic) {
     2318                  dj[i]=0.0;
     2319                } else if (model2->getStatus(i)==atLowerBound) {
     2320                  if (optimizationDirection_*dj[i]<tolerance) {
     2321                    if (optimizationDirection_*dj[i]<0.0) {
     2322                      //if (dj[i]<-1.0e-3)
     2323                      //printf("bad dj at lb %d %g\n",i,dj[i]);
     2324                      cost[i] -= dj[i];
     2325                      dj[i]=0.0;
     2326                    }
     2327                  } else {
     2328                    upper[i]=lower[i];
    22272329                  }
    2228                 } else {
    2229                   upper[i]=lower[i];
    2230                 }
    2231               } else if (model2->getStatus(i)==atUpperBound) {
    2232                 if (optimizationDirection_*dj[i]>tolerance) {
    2233                   if (optimizationDirection_*dj[i]>0.0) {
    2234                     //if (dj[i]>1.0e-3)
    2235                     //printf("bad dj at ub %d %g\n",i,dj[i]);
    2236                     cost[i] -= dj[i];
    2237                     dj[i]=0.0;
     2330                } else if (model2->getStatus(i)==atUpperBound) {
     2331                  if (optimizationDirection_*dj[i]>tolerance) {
     2332                    if (optimizationDirection_*dj[i]>0.0) {
     2333                      //if (dj[i]>1.0e-3)
     2334                      //printf("bad dj at ub %d %g\n",i,dj[i]);
     2335                      cost[i] -= dj[i];
     2336                      dj[i]=0.0;
     2337                    }
     2338                  } else {
     2339                    lower[i]=upper[i];
    22382340                  }
    2239                 } else {
    2240                   lower[i]=upper[i];
    22412341                }
    22422342              }
     2343              // just dual values pass
     2344              //model2->setLogLevel(63);
     2345              //model2->setFactorizationFrequency(1);
     2346              model2->dual(2);
     2347              CoinMemcpyN(saveCost,numberColumns,cost);
     2348              delete [] saveCost;
     2349              CoinMemcpyN(saveLower,numberColumns,lower);
     2350              delete [] saveLower;
     2351              CoinMemcpyN(saveUpper,numberColumns,upper);
     2352              delete [] saveUpper;
    22432353            }
    2244             // just dual values pass
    2245             //model2->setLogLevel(63);
    2246             //model2->setFactorizationFrequency(1);
    2247             model2->dual(2);
    2248             CoinMemcpyN(saveCost,numberColumns,cost);
    2249             delete [] saveCost;
    2250             CoinMemcpyN(saveLower,numberColumns,lower);
    2251             delete [] saveLower;
    2252             CoinMemcpyN(saveUpper,numberColumns,upper);
    2253             delete [] saveUpper;
    22542354          }
    22552355          // and finish
     
    22952395    }
    22962396    if (saveUpper) {
    2297       int numberRows = model2->numberRows();
    2298       int numberColumns = model2->numberColumns();
    2299       CoinMemcpyN(saveLower,numberColumns,model2->columnLower());
    2300       CoinMemcpyN(saveLower+numberColumns,numberRows,model2->rowLower());
     2397      if (!forceFixing) {
     2398        int numberRows = model2->numberRows();
     2399        int numberColumns = model2->numberColumns();
     2400        CoinMemcpyN(saveLower,numberColumns,model2->columnLower());
     2401        CoinMemcpyN(saveLower+numberColumns,numberRows,model2->rowLower());
     2402        CoinMemcpyN(saveUpper,numberColumns,model2->columnUpper());
     2403        CoinMemcpyN(saveUpper+numberColumns,numberRows,model2->rowUpper());
     2404      }
    23012405      delete [] saveLower;
    2302       CoinMemcpyN(saveUpper,numberColumns,model2->columnUpper());
    2303       CoinMemcpyN(saveUpper+numberColumns,numberRows,model2->rowUpper());
    23042406      delete [] saveUpper;
    23052407      saveLower=NULL;
Note: See TracChangeset for help on using the changeset viewer.