Changeset 640 for trunk/Cbc/examples


Ignore:
Timestamp:
Jun 26, 2007 5:17:15 AM (13 years ago)
Author:
forrest
Message:

trunk from branches/devel

Location:
trunk
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • trunk

    • Property svn:externals
      •  

        old new  
        1 MSVisualStudio   https://projects.coin-or.org/svn/MSVisualStudio/trunk/ExternalsDirs/Cbc
        2 BuildTools    https://projects.coin-or.org/svn/BuildTools/stable/0.5
         1MSVisualStudio   https://projects.coin-or.org/svn/MSVisualStudio/branches/devel/ExternalsDirs/Cbc
         2BuildTools    https://projects.coin-or.org/svn/BuildTools/trunk
        33ThirdParty/ASL https://projects.coin-or.org/svn/BuildTools/ThirdParty/ASL/stable/1.0
         4ThirdParty/Blas https://projects.coin-or.org/svn/BuildTools/ThirdParty/Blas/stable/1.0
         5ThirdParty/Lapack https://projects.coin-or.org/svn/BuildTools/ThirdParty/Lapack/stable/1.0
        46Data/Netlib   https://projects.coin-or.org/svn/Data/stable/1.0/Netlib
        57Data/Sample   https://projects.coin-or.org/svn/Data/stable/1.0/Sample
  • trunk/Cbc/examples/CbcBranchLink.cpp

    r134 r640  
    1515#include "CbcBranchLink.hpp"
    1616#include "CoinError.hpp"
     17#include "CoinPackedMatrix.hpp"
    1718
    1819// Default Constructor
     
    2223    numberMembers_(0),
    2324    numberLinks_(0),
    24     first_(-1)
     25    which_(NULL),
     26    sosType_(1)
    2527{
    2628}
     
    3234    numberMembers_(numberMembers),
    3335    numberLinks_(numberLinks),
    34     first_(first)
     36    which_(NULL),
     37    sosType_(1)
    3538{
    3639  id_=identifier;
    3740  if (numberMembers_) {
    3841    weights_ = new double[numberMembers_];
     42    which_ = new int[numberMembers_*numberLinks_];
    3943    if (weights) {
    4044      memcpy(weights_,weights,numberMembers_*sizeof(double));
     
    5054      last=weights_[i];
    5155    }
     56    for (i=0;i<numberMembers_*numberLinks_;i++) {
     57      which_[i]=first+i;
     58    }
     59  } else {
     60    weights_ = NULL;
     61  }
     62}
     63
     64// Useful constructor (which are indices)
     65CbcLink::CbcLink (CbcModel * model,  int numberMembers,
     66           int numberLinks, int sosType, const int * which , const double * weights, int identifier)
     67  : CbcObject(model),
     68    numberMembers_(numberMembers),
     69    numberLinks_(numberLinks),
     70    which_(NULL),
     71    sosType_(sosType)
     72{
     73  id_=identifier;
     74  if (numberMembers_) {
     75    weights_ = new double[numberMembers_];
     76    which_ = new int[numberMembers_*numberLinks_];
     77    if (weights) {
     78      memcpy(weights_,weights,numberMembers_*sizeof(double));
     79    } else {
     80      for (int i=0;i<numberMembers_;i++)
     81        weights_[i]=i;
     82    }
     83    // weights must be increasing
     84    int i;
     85    double last=-COIN_DBL_MAX;
     86    for (i=0;i<numberMembers_;i++) {
     87      assert (weights_[i]>last+1.0e-12);
     88      last=weights_[i];
     89    }
     90    for (i=0;i<numberMembers_*numberLinks_;i++) {
     91      which_[i]= which[i];
     92    }
    5293  } else {
    5394    weights_ = NULL;
     
    61102  numberMembers_ = rhs.numberMembers_;
    62103  numberLinks_ = rhs.numberLinks_;
    63   first_ = rhs.first_;
     104  sosType_ = rhs.sosType_;
    64105  if (numberMembers_) {
    65     weights_ = new double[numberMembers_];
    66     memcpy(weights_,rhs.weights_,numberMembers_*sizeof(double));
     106    weights_ = CoinCopyOfArray(rhs.weights_,numberMembers_);
     107    which_ = CoinCopyOfArray(rhs.which_,numberMembers_*numberLinks_);
    67108  } else {
    68109    weights_ = NULL;
     110    which_ = NULL;
    69111  }
    70112}
     
    84126    CbcObject::operator=(rhs);
    85127    delete [] weights_;
     128    delete [] which_;
    86129    numberMembers_ = rhs.numberMembers_;
    87130    numberLinks_ = rhs.numberLinks_;
    88     first_ = rhs.first_;
     131    sosType_ = rhs.sosType_;
    89132    if (numberMembers_) {
    90       weights_ = new double[numberMembers_];
    91       memcpy(weights_,rhs.weights_,numberMembers_*sizeof(double));
     133      weights_ = CoinCopyOfArray(rhs.weights_,numberMembers_);
     134      which_ = CoinCopyOfArray(rhs.which_,numberMembers_*numberLinks_);
    92135    } else {
    93136      weights_ = NULL;
     137      which_ = NULL;
    94138    }
    95139  }
     
    101145{
    102146  delete [] weights_;
     147  delete [] which_;
    103148}
    104149
     
    112157  OsiSolverInterface * solver = model_->solver();
    113158  const double * solution = model_->testSolution();
    114   const double * lower = solver->getColLower();
     159  //const double * lower = solver->getColLower();
    115160  const double * upper = solver->getColUpper();
    116161  double integerTolerance =
     
    121166  // check bounds etc
    122167  double lastWeight=-1.0e100;
    123   int base=first_;
     168  int base=0;
    124169  for (j=0;j<numberMembers_;j++) {
    125170    for (int k=0;k<numberLinks_;k++) {
    126       int iColumn = base+k;
    127       if (lower[iColumn])
    128         throw CoinError("Non zero lower bound in CBCLink","infeasibility","CbcLink");
     171      int iColumn = which_[base+k];
     172      //if (lower[iColumn])
     173      //throw CoinError("Non zero lower bound in CBCLink","infeasibility","CbcLink");
    129174      if (lastWeight>=weights_[j]-1.0e-7)
    130175        throw CoinError("Weights too close together in CBCLink","infeasibility","CbcLink");
     
    141186#endif
    142187        }
    143         weight += weights_[j]*CoinMin(value,upper[iColumn]);
     188        value = CoinMin(value,upper[iColumn]);
     189        weight += weights_[j]*value;
    144190        if (firstNonZero<0)
    145191          firstNonZero=j;
     
    149195    base += numberLinks_;
    150196  }
     197  double valueInfeasibility;
    151198  preferredWay=1;
    152   if (lastNonZero-firstNonZero>=1) {
     199  if (lastNonZero-firstNonZero>=sosType_) {
    153200    // find where to branch
    154201    assert (sum>0.0);
    155202    weight /= sum;
    156     double value = lastNonZero-firstNonZero+1;
    157     value *= 0.5/((double) numberMembers_);
    158     return value;
     203    valueInfeasibility = lastNonZero-firstNonZero+1;
     204    valueInfeasibility *= 0.5/((double) numberMembers_);
     205    //#define DISTANCE
     206#ifdef DISTANCE
     207    assert (sosType_==1); // code up
     208    /* may still be satisfied.
     209       For LOS type 2 we might wish to move coding around
     210       and keep initial info in model_ for speed
     211    */
     212    int iWhere;
     213    bool possible=false;
     214    for (iWhere=firstNonZero;iWhere<=lastNonZero;iWhere++) {
     215      if (fabs(weight-weights_[iWhere])<1.0e-8) {
     216        possible=true;
     217        break;
     218      }
     219    }
     220    if (possible) {
     221      // One could move some of this (+ arrays) into model_
     222      const CoinPackedMatrix * matrix = solver->getMatrixByCol();
     223      const double * element = matrix->getMutableElements();
     224      const int * row = matrix->getIndices();
     225      const CoinBigIndex * columnStart = matrix->getVectorStarts();
     226      const int * columnLength = matrix->getVectorLengths();
     227      const double * rowSolution = solver->getRowActivity();
     228      const double * rowLower = solver->getRowLower();
     229      const double * rowUpper = solver->getRowUpper();
     230      int numberRows = matrix->getNumRows();
     231      double * array = new double [numberRows];
     232      CoinZeroN(array,numberRows);
     233      int * which = new int [numberRows];
     234      int n=0;
     235      int base=numberLinks_*firstNonZero;
     236      for (j=firstNonZero;j<=lastNonZero;j++) {
     237        for (int k=0;k<numberLinks_;k++) {
     238          int iColumn = which_[base+k];
     239          double value = CoinMax(0.0,solution[iColumn]);
     240          if (value>integerTolerance&&upper[iColumn]) {
     241            value = CoinMin(value,upper[iColumn]);
     242            for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
     243              int iRow = row[j];
     244              double a = array[iRow];
     245              if (a) {
     246                a += value*element[j];
     247                if (!a)
     248                  a = 1.0e-100;
     249              } else {
     250                which[n++]=iRow;
     251                a=value*element[j];
     252                assert (a);
     253              }
     254              array[iRow]=a;
     255            }
     256          }
     257        }
     258        base += numberLinks_;
     259      }
     260      base=numberLinks_*iWhere;
     261      for (int k=0;k<numberLinks_;k++) {
     262        int iColumn = which_[base+k];
     263        const double value = 1.0;
     264        for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
     265          int iRow = row[j];
     266          double a = array[iRow];
     267          if (a) {
     268            a -= value*element[j];
     269            if (!a)
     270              a = 1.0e-100;
     271          } else {
     272            which[n++]=iRow;
     273            a=-value*element[j];
     274            assert (a);
     275          }
     276          array[iRow]=a;
     277        }
     278      }
     279      for (j=0;j<n;j++) {
     280        int iRow = which[j];
     281        // moving to point will increase row solution by this
     282        double distance = array[iRow];
     283        if (distance>1.0e-8) {
     284          if (distance+rowSolution[iRow]>rowUpper[iRow]+1.0e-8) {
     285            possible=false;
     286            break;
     287          }
     288        } else if (distance<-1.0e-8) {
     289          if (distance+rowSolution[iRow]<rowLower[iRow]-1.0e-8) {
     290            possible=false;
     291            break;
     292          }
     293        }
     294      }
     295      for (j=0;j<n;j++)
     296        array[which[j]]=0.0;
     297      delete [] array;
     298      delete [] which;
     299      if (possible) {
     300        valueInfeasibility=0.0;
     301        printf("possible %d %d %d\n",firstNonZero,lastNonZero,iWhere);
     302      }
     303    }
     304#endif
    159305  } else {
    160     return 0.0; // satisfied
    161   }
     306    valueInfeasibility = 0.0; // satisfied
     307  }
     308  return valueInfeasibility;
    162309}
    163310
     
    177324  double sum =0.0;
    178325
    179   int base=first_;
     326  int base=0;
    180327  for (j=0;j<numberMembers_;j++) {
    181328    for (int k=0;k<numberLinks_;k++) {
    182       int iColumn = base+k;
     329      int iColumn = which_[base+k];
    183330      double value = CoinMax(0.0,solution[iColumn]);
    184331      sum += value;
     
    192339    base += numberLinks_;
    193340  }
    194   assert (lastNonZero-firstNonZero==0) ;
    195   base=first_;
     341#ifdef DISTANCE
     342  if (lastNonZero-firstNonZero>sosType_-1) {
     343    /* may still be satisfied.
     344       For LOS type 2 we might wish to move coding around
     345       and keep initial info in model_ for speed
     346    */
     347    int iWhere;
     348    bool possible=false;
     349    for (iWhere=firstNonZero;iWhere<=lastNonZero;iWhere++) {
     350      if (fabs(weight-weights_[iWhere])<1.0e-8) {
     351        possible=true;
     352        break;
     353      }
     354    }
     355    if (possible) {
     356      // One could move some of this (+ arrays) into model_
     357      const CoinPackedMatrix * matrix = solver->getMatrixByCol();
     358      const double * element = matrix->getMutableElements();
     359      const int * row = matrix->getIndices();
     360      const CoinBigIndex * columnStart = matrix->getVectorStarts();
     361      const int * columnLength = matrix->getVectorLengths();
     362      const double * rowSolution = solver->getRowActivity();
     363      const double * rowLower = solver->getRowLower();
     364      const double * rowUpper = solver->getRowUpper();
     365      int numberRows = matrix->getNumRows();
     366      double * array = new double [numberRows];
     367      CoinZeroN(array,numberRows);
     368      int * which = new int [numberRows];
     369      int n=0;
     370      int base=numberLinks_*firstNonZero;
     371      for (j=firstNonZero;j<=lastNonZero;j++) {
     372        for (int k=0;k<numberLinks_;k++) {
     373          int iColumn = which_[base+k];
     374          double value = CoinMax(0.0,solution[iColumn]);
     375          if (value>integerTolerance&&upper[iColumn]) {
     376            value = CoinMin(value,upper[iColumn]);
     377            for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
     378              int iRow = row[j];
     379              double a = array[iRow];
     380              if (a) {
     381                a += value*element[j];
     382                if (!a)
     383                  a = 1.0e-100;
     384              } else {
     385                which[n++]=iRow;
     386                a=value*element[j];
     387                assert (a);
     388              }
     389              array[iRow]=a;
     390            }
     391          }
     392        }
     393        base += numberLinks_;
     394      }
     395      base=numberLinks_*iWhere;
     396      for (int k=0;k<numberLinks_;k++) {
     397        int iColumn = which_[base+k];
     398        const double value = 1.0;
     399        for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
     400          int iRow = row[j];
     401          double a = array[iRow];
     402          if (a) {
     403            a -= value*element[j];
     404            if (!a)
     405              a = 1.0e-100;
     406          } else {
     407            which[n++]=iRow;
     408            a=-value*element[j];
     409            assert (a);
     410          }
     411          array[iRow]=a;
     412        }
     413      }
     414      for (j=0;j<n;j++) {
     415        int iRow = which[j];
     416        // moving to point will increase row solution by this
     417        double distance = array[iRow];
     418        if (distance>1.0e-8) {
     419          if (distance+rowSolution[iRow]>rowUpper[iRow]+1.0e-8) {
     420            possible=false;
     421            break;
     422          }
     423        } else if (distance<-1.0e-8) {
     424          if (distance+rowSolution[iRow]<rowLower[iRow]-1.0e-8) {
     425            possible=false;
     426            break;
     427          }
     428        }
     429      }
     430      for (j=0;j<n;j++)
     431        array[which[j]]=0.0;
     432      delete [] array;
     433      delete [] which;
     434      if (possible) {
     435        printf("possible feas region %d %d %d\n",firstNonZero,lastNonZero,iWhere);
     436        firstNonZero=iWhere;
     437        lastNonZero=iWhere;
     438      }
     439    }
     440  }
     441#else
     442  assert (lastNonZero-firstNonZero<sosType_) ;
     443#endif
     444  base=0;
    196445  for (j=0;j<firstNonZero;j++) {
    197446    for (int k=0;k<numberLinks_;k++) {
    198       int iColumn = base+k;
     447      int iColumn = which_[base+k];
    199448      solver->setColUpper(iColumn,0.0);
    200449    }
     
    205454  for (j=lastNonZero+1;j<numberMembers_;j++) {
    206455    for (int k=0;k<numberLinks_;k++) {
    207       int iColumn = base+k;
     456      int iColumn = which_[base+k];
    208457      solver->setColUpper(iColumn,0.0);
    209458    }
     
    229478  double weight = 0.0;
    230479  double sum =0.0;
    231   int base=first_;
     480  int base=0;
    232481  for (j=0;j<numberMembers_;j++) {
    233482    for (int k=0;k<numberLinks_;k++) {
    234       int iColumn = base+k;
     483      int iColumn = which_[base+k];
    235484      if (upper[iColumn]) {
    236485        double value = CoinMax(0.0,solution[iColumn]);
     
    249498    base += numberLinks_;
    250499  }
    251   assert (lastNonZero-firstNonZero>=1) ;
     500  assert (lastNonZero-firstNonZero>=sosType_) ;
    252501  // find where to branch
    253502  assert (sum>0.0);
     
    258507    if (weight<weights_[iWhere+1])
    259508      break;
    260   separator = 0.5 *(weights_[iWhere]+weights_[iWhere+1]);
     509  if (sosType_==1) {
     510    // SOS 1
     511    separator = 0.5 *(weights_[iWhere]+weights_[iWhere+1]);
     512  } else {
     513    // SOS 2
     514    if (iWhere==firstNonFixed)
     515      iWhere++;;
     516    if (iWhere==lastNonFixed-1)
     517      iWhere = lastNonFixed-2;
     518    separator = weights_[iWhere+1];
     519  }
    261520  // create object
    262521  CbcBranchingObject * branch;
     
    305564}
    306565double
    307 CbcLinkBranchingObject::branch(bool normalBranch)
    308 {
    309   if (model_->messageHandler()->logLevel()>2&&normalBranch)
    310     print(normalBranch);
    311   numberBranchesLeft_--;
     566CbcLinkBranchingObject::branch()
     567{
     568  decrementNumberBranchesLeft();
    312569  int numberMembers = set_->numberMembers();
    313570  int numberLinks = set_->numberLinks();
    314571  const double * weights = set_->weights();
     572  const int * which = set_->which();
    315573  OsiSolverInterface * solver = model_->solver();
    316574  //const double * lower = solver->getColLower();
     
    324582    }
    325583    assert (i<numberMembers);
    326     int base=set_->first()+i*numberLinks;;
     584    int base=i*numberLinks;;
    327585    for (;i<numberMembers;i++) {
    328586      for (int k=0;k<numberLinks;k++) {
    329         int iColumn = base+k;
     587        int iColumn = which[base+k];
    330588        solver->setColUpper(iColumn,0.0);
    331589      }
     
    335593  } else {
    336594    int i;
    337     int base=set_->first();
     595    int base=0;
    338596    for ( i=0;i<numberMembers;i++) {
    339597      if (weights[i] >= separator_) {
     
    341599      } else {
    342600        for (int k=0;k<numberLinks;k++) {
    343           int iColumn = base+k;
     601          int iColumn = which[base+k];
    344602          solver->setColUpper(iColumn,0.0);
    345603        }
     
    354612// Print what would happen 
    355613void
    356 CbcLinkBranchingObject::print(bool normalBranch)
     614CbcLinkBranchingObject::print()
    357615{
    358616  int numberMembers = set_->numberMembers();
    359617  int numberLinks = set_->numberLinks();
    360618  const double * weights = set_->weights();
     619  const int * which = set_->which();
    361620  OsiSolverInterface * solver = model_->solver();
    362621  const double * upper = solver->getColUpper();
     
    366625  int numberOther=0;
    367626  int i;
    368   int base=set_->first();
     627  int base=0;
    369628  for ( i=0;i<numberMembers;i++) {
    370629    for (int k=0;k<numberLinks;k++) {
    371       int iColumn = base+k;
     630      int iColumn = which[base+k];
    372631      double bound = upper[iColumn];
    373632      if (bound) {
     
    379638  }
    380639  // *** for way - up means fix all those in down section
    381   base=set_->first();
     640  base=0;
    382641  if (way_<0) {
    383642    printf("SOS Down");
    384643    for ( i=0;i<numberMembers;i++) {
     644      if (weights[i] > separator_)
     645        break;
    385646      for (int k=0;k<numberLinks;k++) {
    386         int iColumn = base+k;
     647        int iColumn = which[base+k];
    387648        double bound = upper[iColumn];
    388         if (weights[i] > separator_)
    389           break;
    390         else if (bound)
     649        if (bound)
    391650          numberOther++;
    392651      }
     
    396655    for (;i<numberMembers;i++) {
    397656      for (int k=0;k<numberLinks;k++) {
    398         int iColumn = base+k;
     657        int iColumn = which[base+k];
    399658        double bound = upper[iColumn];
    400659        if (bound)
     
    406665    printf("SOS Up");
    407666    for ( i=0;i<numberMembers;i++) {
     667      if (weights[i] >= separator_)
     668        break;
    408669      for (int k=0;k<numberLinks;k++) {
    409         int iColumn = base+k;
     670        int iColumn = which[base+k];
    410671        double bound = upper[iColumn];
    411         if (weights[i] >= separator_)
    412           break;
    413         else if (bound)
     672        if (bound)
    414673          numberFixed++;
    415674      }
     
    419678    for (;i<numberMembers;i++) {
    420679      for (int k=0;k<numberLinks;k++) {
    421         int iColumn = base+k;
     680        int iColumn = which[base+k];
    422681        double bound = upper[iColumn];
    423682        if (bound)
     
    428687  }
    429688  assert ((numberFixed%numberLinks)==0);
    430   assert ((numberFixed%numberOther)==0);
     689  assert ((numberOther%numberLinks)==0);
    431690  printf(" - at %g, free range %d (%g) => %d (%g), %d would be fixed, %d other way\n",
    432691         separator_,first,weights[first],last,weights[last],numberFixed/numberLinks,
  • trunk/Cbc/examples/CbcBranchLink.hpp

    r129 r640  
    2121      apart from k*numberLink to (k+1)*numberLink-1 where k is 0 through
    2222      numberInSet-1.  The length of weights array is numberInSet.
    23       For this simple version the variables in matrix are the numberInSet*numberLink
     23      For this constructor the variables in matrix are the numberInSet*numberLink
    2424      starting at first. If weights null then 0,1,2..
    2525  */
    2626  CbcLink (CbcModel * model, int numberMembers,
    2727           int numberLinks, int first,
     28           const double * weights, int setNumber);
     29  /** Useful constructor - A valid solution is if all variables are zero
     30      apart from k*numberLink to (k+1)*numberLink-1 where k is 0 through
     31      numberInSet-1.  The length of weights array is numberInSet.
     32      For this constructor the variables are given by list - grouped.
     33      If weights null then 0,1,2..
     34  */
     35  CbcLink (CbcModel * model, int numberMembers,
     36           int numberLinks, int typeSOS, const int * which,
    2837           const double * weights, int setNumber);
    2938 
     
    5665  {return numberLinks_;};
    5766
    58   /// First variable in matrix
    59   inline int first() const
    60   {return first_;};
     67  /// Which variables
     68  inline const int * which() const
     69  {return which_;};
    6170
    6271  /** Array of weights */
     
    7483  /// Number of links
    7584   int numberLinks_;
    76   /// First member
    77   int first_;
     85  /// Members
     86  int * which_;
     87  /// Type 1 or 2
     88  int sosType_;
    7889};
    7990/** Branching object for Special ordered sets
     
    107118 
    108119  /// Does next branch and updates state
    109   virtual double branch(bool normalBranch=false);
     120  virtual double branch();
    110121
    111122  /** \brief Print something about branch - only if log level high
    112123  */
    113   virtual void print(bool normalBranch);
     124  virtual void print();
    114125private:
    115126  /// data
  • trunk/Cbc/examples/CbcBranchUser.cpp

    r2 r640  
    323323  return whichObject;
    324324}
     325/** Default Constructor
     326
     327  Equivalent to an unspecified binary variable.
     328*/
     329CbcSimpleIntegerFixed::CbcSimpleIntegerFixed ()
     330  : CbcSimpleInteger()
     331{
     332}
     333
     334/** Useful constructor
     335
     336  Loads actual upper & lower bounds for the specified variable.
     337*/
     338CbcSimpleIntegerFixed::CbcSimpleIntegerFixed (CbcModel * model,
     339                                    int iColumn, double breakEven)
     340  : CbcSimpleInteger(model,iColumn,breakEven)
     341{
     342}
     343// Constructor from simple
     344CbcSimpleIntegerFixed::CbcSimpleIntegerFixed (const CbcSimpleInteger & rhs)
     345  : CbcSimpleInteger(rhs)
     346{
     347}
     348
     349// Copy constructor
     350CbcSimpleIntegerFixed::CbcSimpleIntegerFixed ( const CbcSimpleIntegerFixed & rhs)
     351  :CbcSimpleInteger(rhs)
     352
     353{
     354}
     355
     356// Clone
     357CbcObject *
     358CbcSimpleIntegerFixed::clone() const
     359{
     360  return new CbcSimpleIntegerFixed(*this);
     361}
     362
     363// Assignment operator
     364CbcSimpleIntegerFixed &
     365CbcSimpleIntegerFixed::operator=( const CbcSimpleIntegerFixed& rhs)
     366{
     367  if (this!=&rhs) {
     368    CbcSimpleInteger::operator=(rhs);
     369  }
     370  return *this;
     371}
     372
     373// Destructor
     374CbcSimpleIntegerFixed::~CbcSimpleIntegerFixed ()
     375{
     376}
     377
     378// Infeasibility - large is 0.5
     379double
     380CbcSimpleIntegerFixed::infeasibility(int & preferredWay) const
     381{
     382  OsiSolverInterface * solver = model_->solver();
     383  const double * solution = model_->testSolution();
     384  const double * lower = solver->getColLower();
     385  const double * upper = solver->getColUpper();
     386  double value = solution[columnNumber_];
     387  value = CoinMax(value, lower[columnNumber_]);
     388  value = CoinMin(value, upper[columnNumber_]);
     389  /*printf("%d %g %g %g %g\n",columnNumber_,value,lower[columnNumber_],
     390    solution[columnNumber_],upper[columnNumber_]);*/
     391  double nearest = floor(value+(1.0-breakEven_));
     392  assert (breakEven_>0.0&&breakEven_<1.0);
     393  double integerTolerance =
     394    model_->getDblParam(CbcModel::CbcIntegerTolerance);
     395  if (nearest>value)
     396    preferredWay=1;
     397  else
     398    preferredWay=-1;
     399  if (preferredWay_)
     400    preferredWay=preferredWay_;
     401  double weight = fabs(value-nearest);
     402  // normalize so weight is 0.5 at break even
     403  if (nearest<value)
     404    weight = (0.5/breakEven_)*weight;
     405  else
     406    weight = (0.5/(1.0-breakEven_))*weight;
     407  if (fabs(value-nearest)<=integerTolerance) {
     408    if (upper[columnNumber_]==lower[columnNumber_])
     409      return 0.0;
     410    else
     411      return 1.0e-5;
     412  } else {
     413    return weight;
     414  }
     415}
     416// Creates a branching object
     417CbcBranchingObject *
     418CbcSimpleIntegerFixed::createBranch(OsiSolverInterface * solver,
     419                                            const OsiBranchingInformation * info, int way) 
     420{
     421  const double * solution = model_->testSolution();
     422  const double * lower = solver->getColLower();
     423  const double * upper = solver->getColUpper();
     424  double value = solution[columnNumber_];
     425  value = CoinMax(value, lower[columnNumber_]);
     426  value = CoinMin(value, upper[columnNumber_]);
     427  assert (upper[columnNumber_]>lower[columnNumber_]);
     428  if (!model_->hotstartSolution()) {
     429    double nearest = floor(value+0.5);
     430    double integerTolerance =
     431    model_->getDblParam(CbcModel::CbcIntegerTolerance);
     432    if (fabs(value-nearest)<integerTolerance) {
     433      // adjust value
     434      if (nearest!=upper[columnNumber_])
     435        value = nearest+2.0*integerTolerance;
     436      else
     437        value = nearest-2.0*integerTolerance;
     438    }
     439  } else {
     440    const double * hotstartSolution = model_->hotstartSolution();
     441    double targetValue = hotstartSolution[columnNumber_];
     442    if (way>0)
     443      value = targetValue-0.1;
     444    else
     445      value = targetValue+0.1;
     446  }
     447  CbcBranchingObject * branch = new CbcIntegerBranchingObject(model_,columnNumber_,way,
     448                                             value);
     449  branch->setOriginalObject(this);
     450  return branch;
     451}
  • trunk/Cbc/examples/CbcBranchUser.hpp

    r2 r640  
    55
    66#include "CbcBranchBase.hpp"
     7#include "CbcBranchActual.hpp"
    78
    89/** Branching decision user class */
     
    5556};
    5657
     58/// Define a single integer class where branching is forced until fixed
     59
     60
     61class CbcSimpleIntegerFixed : public CbcSimpleInteger {
     62
     63public:
     64
     65  // Default Constructor
     66  CbcSimpleIntegerFixed ();
     67
     68  // Useful constructor - passed integer index and model index
     69  CbcSimpleIntegerFixed (CbcModel * model, int iColumn, double breakEven=0.5);
     70 
     71  // Constructor from simple
     72  CbcSimpleIntegerFixed (const CbcSimpleInteger & simple);
     73 
     74  // Copy constructor
     75  CbcSimpleIntegerFixed ( const CbcSimpleIntegerFixed &);
     76   
     77  /// Clone
     78  virtual CbcObject * clone() const;
     79
     80  // Assignment operator
     81  CbcSimpleIntegerFixed & operator=( const CbcSimpleIntegerFixed& rhs);
     82
     83  // Destructor
     84  ~CbcSimpleIntegerFixed ();
     85 
     86  /// Infeasibility - large is 0.5
     87  virtual double infeasibility(int & preferredWay) const;
     88
     89  /** Creates a branching object
     90
     91    The preferred direction is set by \p way, -1 for down, +1 for up.
     92  */
     93  //virtual CbcBranchingObject * createBranch(int way) ;
     94  /** Create a branching object and indicate which way to branch first.
     95     
     96      The branching object has to know how to create branches (fix
     97      variables, etc.)
     98  */
     99  virtual CbcBranchingObject * createBranch(OsiSolverInterface * solver,
     100                                            const OsiBranchingInformation * info, int way) ;
     101
     102protected:
     103  /// data
     104};
     105
    57106#endif
  • trunk/Cbc/examples/CbcCompareUser.cpp

    r180 r640  
    8282{
    8383}
    84 
     84// For moment go to default
     85#if 0
    8586// Returns true if y better than x
    8687bool
     
    179180  }
    180181}
     182#else
     183
     184// Returns true if y better than x
     185bool
     186CbcCompareUser::test (CbcNode * x, CbcNode * y)
     187{
     188  if (weight_==-1.0&&(y->depth()>7||x->depth()>7)) {
     189    // before solution
     190    /* printf("x %d %d %g, y %d %d %g\n",
     191       x->numberUnsatisfied(),x->depth(),x->objectiveValue(),
     192       y->numberUnsatisfied(),y->depth(),y->objectiveValue()); */
     193    if (x->numberUnsatisfied() > y->numberUnsatisfied()) {
     194      return true;
     195    } else if (x->numberUnsatisfied() < y->numberUnsatisfied()) {
     196      return false;
     197    } else {
     198      int testX = x->depth();
     199      int testY = y->depth();
     200      if (testX!=testY)
     201        return testX < testY;
     202      else
     203        return equalityTest(x,y); // so ties will be broken in consistent manner
     204    }
     205  } else {
     206    // after solution
     207    double weight = CoinMax(weight_,0.0);
     208    double testX =  x->objectiveValue()+ weight*x->numberUnsatisfied();
     209    double testY = y->objectiveValue() + weight*y->numberUnsatisfied();
     210    if (testX!=testY)
     211      return testX > testY;
     212    else
     213      return equalityTest(x,y); // so ties will be broken in consistent manner
     214  }
     215}
     216// This allows method to change behavior as it is called
     217// after each solution
     218void
     219CbcCompareUser::newSolution(CbcModel * model,
     220                               double objectiveAtContinuous,
     221                               int numberInfeasibilitiesAtContinuous)
     222{
     223  if (model->getSolutionCount()==model->getNumberHeuristicSolutions()&&
     224      model->getSolutionCount()<5&&model->getNodeCount()<500)
     225    return; // solution was got by rounding
     226  // set to get close to this solution
     227  double costPerInteger =
     228    (model->getObjValue()-objectiveAtContinuous)/
     229    ((double) numberInfeasibilitiesAtContinuous);
     230  weight_ = 0.95*costPerInteger;
     231  saveWeight_ = 0.95*weight_;
     232  numberSolutions_++;
     233  if (numberSolutions_>5)
     234    weight_ =0.0; // this searches on objective
     235}
     236// This allows method to change behavior
     237bool
     238CbcCompareUser::every1000Nodes(CbcModel * model, int numberNodes)
     239{
     240  double saveWeight=weight_;
     241  int numberNodes1000 = numberNodes/1000;
     242  if (numberNodes>10000) {
     243    weight_ =0.0; // this searches on objective
     244    // but try a bit of other stuff
     245    if ((numberNodes1000%4)==1)
     246      weight_=saveWeight_;
     247  } else if (numberNodes==1000&&weight_==-2.0) {
     248    weight_=-1.0; // Go to depth first
     249  }
     250  // get size of tree
     251  treeSize_ = model->tree()->size();
     252  if (treeSize_>10000) {
     253    int n1 = model->solver()->getNumRows()+model->solver()->getNumCols();
     254    int n2 = model->numberObjects();
     255    double size = n1*0.1 + n2*2.0;
     256    // set weight to reduce size most of time
     257    if (treeSize_*size>5.0e7)
     258      weight_=-1.0;
     259    else if ((numberNodes1000%4)==0&&treeSize_*size>1.0e6)
     260      weight_=-1.0;
     261    else if ((numberNodes1000%4)==1)
     262      weight_=0.0;
     263    else
     264      weight_=saveWeight_;
     265  }
     266  return (weight_!=saveWeight);
     267}
     268// Returns true if wants code to do scan with alternate criterion
     269bool
     270CbcCompareUser::fullScan() const
     271{
     272  return false;
     273}
     274// This is alternate test function
     275bool
     276CbcCompareUser::alternateTest (CbcNode * x, CbcNode * y)
     277{
     278  // not used
     279  abort();
     280  return false;
     281}
     282#endif
  • trunk/Cbc/examples/Makefile.in

    r461 r640  
    2323        CbcSolver2.@OBJEXT@ \
    2424        CbcSolver3.@OBJEXT@ \
     25        CbcSolverLink.@OBJEXT@ \
     26        ClpDynamicInterface.@OBJEXT@ \
    2527        ClpQuadInterface.@OBJEXT@ \
     28        CoinWarmStartBasisDynamic.@OBJEXT@ \
    2629        CbcBranchFollow2.@OBJEXT@ \
     30        CbcBranchLotsizeSimple.@OBJEXT@ \
    2731        CbcBranchUser.@OBJEXT@ \
    2832        CbcBranchLink.@OBJEXT@
     
    8185        driver2.@OBJEXT@ driver2@EXEEXT@ \
    8286        driver.@OBJEXT@ driver@EXEEXT@ \
     87        dynamic.@OBJEXT@ dynamic@EXEEXT@ \
    8388        fast0507b.@OBJEXT@ fast0507b@EXEEXT@ \
    8489        fast0507.@OBJEXT@ fast0507@EXEEXT@ \
     
    111116
    112117clean:
    113         rm -rf $(CLEANFILES)
     118        rm -rf $(CLEANFILES) $(OBJS)
    114119
    115120.cpp.o:
  • trunk/Cbc/examples/qmip2.cpp

    r394 r640  
    188188  // This clones solver
    189189  CbcModel model(solver1);
    190   // Add stored cuts
    191   model.addCutGenerator(&stored,-1,"Stored");
     190  // Add stored cuts (making sure at all depths)
     191  model.addCutGenerator(&stored,1,"Stored",true,false,false,-100,1,-1);
    192192  /*  You need the next few lines -
    193193      a) so that cut generator will always be called again if it generated cuts
    194194      b) it is known that matrix is not enough to define problem so do cuts even
    195195         if it looks integer feasible at continuous optimum.
    196       c) a solution found ny strong branching will be ignored.
     196      c) a solution found by strong branching will be ignored.
     197      d) don't recompute a solution once found
    197198  */
    198199  // Make sure cut generator called correctly (a)
     
    205206  CbcFeasibilityNoStrong noStrong;
    206207  model.setProblemFeasibility(noStrong);
     208  // Say don't recompute solution d)
     209  model.setSpecialOptions(4);
    207210 
    208211  double time1 = CoinCpuTime();
Note: See TracChangeset for help on using the changeset viewer.