Ignore:
Timestamp:
Jan 6, 2019 6:17:46 PM (10 months ago)
Author:
unxusr
Message:

formatting

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Cbc/examples/CbcBranchLink.cpp

    r1900 r2469  
    1717#include "CoinPackedMatrix.hpp"
    1818
    19 // Default Constructor 
    20 CbcLink::CbcLink ()
    21   : CbcObject(),
    22     weights_(NULL),
    23     numberMembers_(0),
    24     numberLinks_(0),
    25     which_(NULL),
    26     sosType_(1)
     19// Default Constructor
     20CbcLink::CbcLink()
     21  : CbcObject()
     22  , weights_(NULL)
     23  , numberMembers_(0)
     24  , numberLinks_(0)
     25  , which_(NULL)
     26  , sosType_(1)
    2727{
    2828}
    2929
    3030// Useful constructor (which are indices)
    31 CbcLink::CbcLink (CbcModel * model, int numberMembers,
    32            int numberLinks, int first , const double * weights, int identifier)
    33   : CbcObject(model),
    34     numberMembers_(numberMembers),
    35     numberLinks_(numberLinks),
    36     which_(NULL),
    37     sosType_(1)
    38 {
    39   id_=identifier;
     31CbcLink::CbcLink(CbcModel *model, int numberMembers,
     32  int numberLinks, int first, const double *weights, int identifier)
     33  : CbcObject(model)
     34  , numberMembers_(numberMembers)
     35  , numberLinks_(numberLinks)
     36  , which_(NULL)
     37  , sosType_(1)
     38{
     39  id_ = identifier;
    4040  if (numberMembers_) {
    4141    weights_ = new double[numberMembers_];
    42     which_ = new int[numberMembers_*numberLinks_];
     42    which_ = new int[numberMembers_ * numberLinks_];
    4343    if (weights) {
    44       memcpy(weights_,weights,numberMembers_*sizeof(double));
     44      memcpy(weights_, weights, numberMembers_ * sizeof(double));
    4545    } else {
    46       for (int i=0;i<numberMembers_;i++)
    47         weights_[i]=i;
     46      for (int i = 0; i < numberMembers_; i++)
     47        weights_[i] = i;
    4848    }
    4949    // weights must be increasing
    5050    int i;
    51     for (i=0;i<numberMembers_;i++)
    52       assert (i == 0 || weights_[i]>weights_[i-1]+1.0e-12);
    53     for (i=0;i<numberMembers_*numberLinks_;i++) {
    54       which_[i]=first+i;
     51    for (i = 0; i < numberMembers_; i++)
     52      assert(i == 0 || weights_[i] > weights_[i - 1] + 1.0e-12);
     53    for (i = 0; i < numberMembers_ * numberLinks_; i++) {
     54      which_[i] = first + i;
    5555    }
    5656  } else {
     
    6060
    6161// Useful constructor (which are indices)
    62 CbcLink::CbcLink (CbcModel * model, int numberMembers,
    63            int numberLinks, int sosType, const int * which , const double * weights, int identifier)
    64   : CbcObject(model),
    65     numberMembers_(numberMembers),
    66     numberLinks_(numberLinks),
    67     which_(NULL),
    68     sosType_(sosType)
    69 {
    70   id_=identifier;
     62CbcLink::CbcLink(CbcModel *model, int numberMembers,
     63  int numberLinks, int sosType, const int *which, const double *weights, int identifier)
     64  : CbcObject(model)
     65  , numberMembers_(numberMembers)
     66  , numberLinks_(numberLinks)
     67  , which_(NULL)
     68  , sosType_(sosType)
     69{
     70  id_ = identifier;
    7171  if (numberMembers_) {
    7272    weights_ = new double[numberMembers_];
    73     which_ = new int[numberMembers_*numberLinks_];
     73    which_ = new int[numberMembers_ * numberLinks_];
    7474    if (weights) {
    75       memcpy(weights_,weights,numberMembers_*sizeof(double));
     75      memcpy(weights_, weights, numberMembers_ * sizeof(double));
    7676    } else {
    77       for (int i=0;i<numberMembers_;i++)
    78         weights_[i]=i;
     77      for (int i = 0; i < numberMembers_; i++)
     78        weights_[i] = i;
    7979    }
    8080    // weights must be increasing
    8181    int i;
    82     for (i=0;i<numberMembers_;i++)
    83       assert (i == 0 || weights_[i]>weights_[i-1]+1.0e-12);
    84     for (i=0;i<numberMembers_*numberLinks_;i++) {
    85       which_[i]= which[i];
     82    for (i = 0; i < numberMembers_; i++)
     83      assert(i == 0 || weights_[i] > weights_[i - 1] + 1.0e-12);
     84    for (i = 0; i < numberMembers_ * numberLinks_; i++) {
     85      which_[i] = which[i];
    8686    }
    8787  } else {
     
    9090}
    9191
    92 // Copy constructor 
    93 CbcLink::CbcLink ( const CbcLink & rhs)
    94   :CbcObject(rhs)
     92// Copy constructor
     93CbcLink::CbcLink(const CbcLink &rhs)
     94  : CbcObject(rhs)
    9595{
    9696  numberMembers_ = rhs.numberMembers_;
     
    9898  sosType_ = rhs.sosType_;
    9999  if (numberMembers_) {
    100     weights_ = CoinCopyOfArray(rhs.weights_,numberMembers_);
    101     which_ = CoinCopyOfArray(rhs.which_,numberMembers_*numberLinks_);
     100    weights_ = CoinCopyOfArray(rhs.weights_, numberMembers_);
     101    which_ = CoinCopyOfArray(rhs.which_, numberMembers_ * numberLinks_);
    102102  } else {
    103103    weights_ = NULL;
     
    113113}
    114114
    115 // Assignment operator 
    116 CbcLink & 
    117 CbcLink::operator=( const CbcLink& rhs)
    118 {
    119   if (this!=&rhs) {
     115// Assignment operator
     116CbcLink &
     117CbcLink::operator=(const CbcLink &rhs)
     118{
     119  if (this != &rhs) {
    120120    CbcObject::operator=(rhs);
    121     delete [] weights_;
    122     delete [] which_;
     121    delete[] weights_;
     122    delete[] which_;
    123123    numberMembers_ = rhs.numberMembers_;
    124124    numberLinks_ = rhs.numberLinks_;
    125125    sosType_ = rhs.sosType_;
    126126    if (numberMembers_) {
    127       weights_ = CoinCopyOfArray(rhs.weights_,numberMembers_);
    128       which_ = CoinCopyOfArray(rhs.which_,numberMembers_*numberLinks_);
     127      weights_ = CoinCopyOfArray(rhs.weights_, numberMembers_);
     128      which_ = CoinCopyOfArray(rhs.which_, numberMembers_ * numberLinks_);
    129129    } else {
    130130      weights_ = NULL;
     
    135135}
    136136
    137 // Destructor 
    138 CbcLink::~CbcLink ()
    139 {
    140   delete [] weights_;
    141   delete [] which_;
     137// Destructor
     138CbcLink::~CbcLink()
     139{
     140  delete[] weights_;
     141  delete[] which_;
    142142}
    143143
    144144// Infeasibility - large is 0.5
    145 double 
    146 CbcLink::infeasibility(int & preferredWay) const
     145double
     146CbcLink::infeasibility(int &preferredWay) const
    147147{
    148148  int j;
    149   int firstNonZero=-1;
     149  int firstNonZero = -1;
    150150  int lastNonZero = -1;
    151   OsiSolverInterface * solver = model_->solver();
    152   const double * solution = model_->testSolution();
     151  OsiSolverInterface *solver = model_->solver();
     152  const double *solution = model_->testSolution();
    153153  //const double * lower = solver->getColLower();
    154   const double * upper = solver->getColUpper();
    155   double integerTolerance =
    156     model_->getDblParam(CbcModel::CbcIntegerTolerance);
     154  const double *upper = solver->getColUpper();
     155  double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);
    157156  double weight = 0.0;
    158   double sum =0.0;
     157  double sum = 0.0;
    159158
    160159  // check bounds etc
    161   double lastWeight=-1.0e100;
    162   int base=0;
    163   for (j=0;j<numberMembers_;j++) {
    164     for (int k=0;k<numberLinks_;k++) {
    165       int iColumn = which_[base+k];
     160  double lastWeight = -1.0e100;
     161  int base = 0;
     162  for (j = 0; j < numberMembers_; j++) {
     163    for (int k = 0; k < numberLinks_; k++) {
     164      int iColumn = which_[base + k];
    166165      //if (lower[iColumn])
    167166      //throw CoinError("Non zero lower bound in CBCLink","infeasibility","CbcLink");
    168       if (lastWeight>=weights_[j]-1.0e-7)
    169         throw CoinError("Weights too close together in CBCLink","infeasibility","CbcLink");
    170       double value = CoinMax(0.0,solution[iColumn]);
     167      if (lastWeight >= weights_[j] - 1.0e-7)
     168        throw CoinError("Weights too close together in CBCLink", "infeasibility", "CbcLink");
     169      double value = CoinMax(0.0, solution[iColumn]);
    171170      sum += value;
    172       if (value>integerTolerance&&upper[iColumn]) {
     171      if (value > integerTolerance && upper[iColumn]) {
    173172        // Possibly due to scaling a fixed variable might slip through
    174         if (value>upper[iColumn]+1.0e-8) {
     173        if (value > upper[iColumn] + 1.0e-8) {
    175174          // Could change to #ifdef CBC_DEBUG
    176175#ifndef NDEBUG
    177           if (model_->messageHandler()->logLevel()>1)
     176          if (model_->messageHandler()->logLevel() > 1)
    178177            printf("** Variable %d (%d) has value %g and upper bound of %g\n",
    179                    iColumn,j,value,upper[iColumn]);
     178              iColumn, j, value, upper[iColumn]);
    180179#endif
    181         } 
    182         value = CoinMin(value,upper[iColumn]);
    183         weight += weights_[j]*value;
    184         if (firstNonZero<0)
    185           firstNonZero=j;
    186         lastNonZero=j;
     180        }
     181        value = CoinMin(value, upper[iColumn]);
     182        weight += weights_[j] * value;
     183        if (firstNonZero < 0)
     184          firstNonZero = j;
     185        lastNonZero = j;
    187186      }
    188187    }
     
    190189  }
    191190  double valueInfeasibility;
    192   preferredWay=1;
    193   if (lastNonZero-firstNonZero>=sosType_) {
     191  preferredWay = 1;
     192  if (lastNonZero - firstNonZero >= sosType_) {
    194193    // find where to branch
    195     assert (sum>0.0);
     194    assert(sum > 0.0);
    196195    weight /= sum;
    197     valueInfeasibility = lastNonZero-firstNonZero+1;
    198     valueInfeasibility *= 0.5/((double) numberMembers_);
     196    valueInfeasibility = lastNonZero - firstNonZero + 1;
     197    valueInfeasibility *= 0.5 / ((double)numberMembers_);
    199198    //#define DISTANCE
    200199#ifdef DISTANCE
    201     assert (sosType_==1); // code up
     200    assert(sosType_ == 1); // code up
    202201    /* may still be satisfied.
    203202       For LOS type 2 we might wish to move coding around
     
    205204    */
    206205    int iWhere;
    207     bool possible=false;
    208     for (iWhere=firstNonZero;iWhere<=lastNonZero;iWhere++) {
    209       if (fabs(weight-weights_[iWhere])<1.0e-8) {
    210         possible=true;
    211         break;
     206    bool possible = false;
     207    for (iWhere = firstNonZero; iWhere <= lastNonZero; iWhere++) {
     208      if (fabs(weight - weights_[iWhere]) < 1.0e-8) {
     209        possible = true;
     210        break;
    212211      }
    213212    }
    214213    if (possible) {
    215214      // One could move some of this (+ arrays) into model_
    216       const CoinPackedMatrix * matrix = solver->getMatrixByCol();
    217       const double * element = matrix->getMutableElements();
    218       const int * row = matrix->getIndices();
    219       const CoinBigIndex * columnStart = matrix->getVectorStarts();
    220       const int * columnLength = matrix->getVectorLengths();
    221       const double * rowSolution = solver->getRowActivity();
    222       const double * rowLower = solver->getRowLower();
    223       const double * rowUpper = solver->getRowUpper();
     215      const CoinPackedMatrix *matrix = solver->getMatrixByCol();
     216      const double *element = matrix->getMutableElements();
     217      const int *row = matrix->getIndices();
     218      const CoinBigIndex *columnStart = matrix->getVectorStarts();
     219      const int *columnLength = matrix->getVectorLengths();
     220      const double *rowSolution = solver->getRowActivity();
     221      const double *rowLower = solver->getRowLower();
     222      const double *rowUpper = solver->getRowUpper();
    224223      int numberRows = matrix->getNumRows();
    225       double * array = new double [numberRows];
    226       CoinZeroN(array,numberRows);
    227       int * which = new int [numberRows];
    228       int n=0;
    229       int base=numberLinks_*firstNonZero;
    230       for (j=firstNonZero;j<=lastNonZero;j++) {
    231         for (int k=0;k<numberLinks_;k++) {
    232           int iColumn = which_[base+k];
    233           double value = CoinMax(0.0,solution[iColumn]);
    234           if (value>integerTolerance&&upper[iColumn]) {
    235             value = CoinMin(value,upper[iColumn]);
    236             for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
    237               int iRow = row[j];
    238               double a = array[iRow];
    239               if (a) {
    240                 a += value*element[j];
    241                 if (!a)
    242                   a = 1.0e-100;
    243               } else {
    244                 which[n++]=iRow;
    245                 a=value*element[j];
    246                 assert (a);
    247               }
    248               array[iRow]=a;
    249             }
    250           }
    251         }
    252         base += numberLinks_;
    253       }
    254       base=numberLinks_*iWhere;
    255       for (int k=0;k<numberLinks_;k++) {
    256         int iColumn = which_[base+k];
    257         const double value = 1.0;
    258         for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
    259           int iRow = row[j];
    260           double a = array[iRow];
    261           if (a) {
    262             a -= value*element[j];
    263             if (!a)
    264               a = 1.0e-100;
    265           } else {
    266             which[n++]=iRow;
    267             a=-value*element[j];
    268             assert (a);
    269           }
    270           array[iRow]=a;
    271         }
    272       }
    273       for (j=0;j<n;j++) {
    274         int iRow = which[j];
    275         // moving to point will increase row solution by this
    276         double distance = array[iRow];
    277         if (distance>1.0e-8) {
    278           if (distance+rowSolution[iRow]>rowUpper[iRow]+1.0e-8) {
    279             possible=false;
    280             break;
    281           }
    282         } else if (distance<-1.0e-8) {
    283           if (distance+rowSolution[iRow]<rowLower[iRow]-1.0e-8) {
    284             possible=false;
    285             break;
    286           }
    287         }
    288       }
    289       for (j=0;j<n;j++)
    290         array[which[j]]=0.0;
    291       delete [] array;
    292       delete [] which;
     224      double *array = new double[numberRows];
     225      CoinZeroN(array, numberRows);
     226      int *which = new int[numberRows];
     227      int n = 0;
     228      int base = numberLinks_ * firstNonZero;
     229      for (j = firstNonZero; j <= lastNonZero; j++) {
     230        for (int k = 0; k < numberLinks_; k++) {
     231          int iColumn = which_[base + k];
     232          double value = CoinMax(0.0, solution[iColumn]);
     233          if (value > integerTolerance && upper[iColumn]) {
     234            value = CoinMin(value, upper[iColumn]);
     235            for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     236              int iRow = row[j];
     237              double a = array[iRow];
     238              if (a) {
     239                a += value * element[j];
     240                if (!a)
     241                  a = 1.0e-100;
     242              } else {
     243                which[n++] = iRow;
     244                a = value * element[j];
     245                assert(a);
     246              }
     247              array[iRow] = a;
     248            }
     249          }
     250        }
     251        base += numberLinks_;
     252      }
     253      base = numberLinks_ * iWhere;
     254      for (int k = 0; k < numberLinks_; k++) {
     255        int iColumn = which_[base + k];
     256        const double value = 1.0;
     257        for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     258          int iRow = row[j];
     259          double a = array[iRow];
     260          if (a) {
     261            a -= value * element[j];
     262            if (!a)
     263              a = 1.0e-100;
     264          } else {
     265            which[n++] = iRow;
     266            a = -value * element[j];
     267            assert(a);
     268          }
     269          array[iRow] = a;
     270        }
     271      }
     272      for (j = 0; j < n; j++) {
     273        int iRow = which[j];
     274        // moving to point will increase row solution by this
     275        double distance = array[iRow];
     276        if (distance > 1.0e-8) {
     277          if (distance + rowSolution[iRow] > rowUpper[iRow] + 1.0e-8) {
     278            possible = false;
     279            break;
     280          }
     281        } else if (distance < -1.0e-8) {
     282          if (distance + rowSolution[iRow] < rowLower[iRow] - 1.0e-8) {
     283            possible = false;
     284            break;
     285          }
     286        }
     287      }
     288      for (j = 0; j < n; j++)
     289        array[which[j]] = 0.0;
     290      delete[] array;
     291      delete[] which;
    293292      if (possible) {
    294         valueInfeasibility=0.0;
    295         printf("possible %d %d %d\n",firstNonZero,lastNonZero,iWhere);
     293        valueInfeasibility = 0.0;
     294        printf("possible %d %d %d\n", firstNonZero, lastNonZero, iWhere);
    296295      }
    297296    }
     
    304303
    305304// This looks at solution and sets bounds to contain solution
    306 void
    307 CbcLink::feasibleRegion()
     305void CbcLink::feasibleRegion()
    308306{
    309307  int j;
    310   int firstNonZero=-1;
     308  int firstNonZero = -1;
    311309  int lastNonZero = -1;
    312   OsiSolverInterface * solver = model_->solver();
    313   const double * solution = model_->testSolution();
    314   const double * upper = solver->getColUpper();
    315   double integerTolerance =
    316     model_->getDblParam(CbcModel::CbcIntegerTolerance);
     310  OsiSolverInterface *solver = model_->solver();
     311  const double *solution = model_->testSolution();
     312  const double *upper = solver->getColUpper();
     313  double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);
    317314  double weight = 0.0;
    318   double sum =0.0;
    319 
    320   int base=0;
    321   for (j=0;j<numberMembers_;j++) {
    322     for (int k=0;k<numberLinks_;k++) {
    323       int iColumn = which_[base+k];
    324       double value = CoinMax(0.0,solution[iColumn]);
     315  double sum = 0.0;
     316
     317  int base = 0;
     318  for (j = 0; j < numberMembers_; j++) {
     319    for (int k = 0; k < numberLinks_; k++) {
     320      int iColumn = which_[base + k];
     321      double value = CoinMax(0.0, solution[iColumn]);
    325322      sum += value;
    326       if (value>integerTolerance&&upper[iColumn]) {
    327         weight += weights_[j]*value;
    328         if (firstNonZero<0)
    329           firstNonZero=j;
    330         lastNonZero=j;
     323      if (value > integerTolerance && upper[iColumn]) {
     324        weight += weights_[j] * value;
     325        if (firstNonZero < 0)
     326          firstNonZero = j;
     327        lastNonZero = j;
    331328      }
    332329    }
     
    334331  }
    335332#ifdef DISTANCE
    336   if (lastNonZero-firstNonZero>sosType_-1) {
     333  if (lastNonZero - firstNonZero > sosType_ - 1) {
    337334    /* may still be satisfied.
    338335       For LOS type 2 we might wish to move coding around
     
    340337    */
    341338    int iWhere;
    342     bool possible=false;
    343     for (iWhere=firstNonZero;iWhere<=lastNonZero;iWhere++) {
    344       if (fabs(weight-weights_[iWhere])<1.0e-8) {
    345         possible=true;
    346         break;
     339    bool possible = false;
     340    for (iWhere = firstNonZero; iWhere <= lastNonZero; iWhere++) {
     341      if (fabs(weight - weights_[iWhere]) < 1.0e-8) {
     342        possible = true;
     343        break;
    347344      }
    348345    }
    349346    if (possible) {
    350347      // One could move some of this (+ arrays) into model_
    351       const CoinPackedMatrix * matrix = solver->getMatrixByCol();
    352       const double * element = matrix->getMutableElements();
    353       const int * row = matrix->getIndices();
    354       const CoinBigIndex * columnStart = matrix->getVectorStarts();
    355       const int * columnLength = matrix->getVectorLengths();
    356       const double * rowSolution = solver->getRowActivity();
    357       const double * rowLower = solver->getRowLower();
    358       const double * rowUpper = solver->getRowUpper();
     348      const CoinPackedMatrix *matrix = solver->getMatrixByCol();
     349      const double *element = matrix->getMutableElements();
     350      const int *row = matrix->getIndices();
     351      const CoinBigIndex *columnStart = matrix->getVectorStarts();
     352      const int *columnLength = matrix->getVectorLengths();
     353      const double *rowSolution = solver->getRowActivity();
     354      const double *rowLower = solver->getRowLower();
     355      const double *rowUpper = solver->getRowUpper();
    359356      int numberRows = matrix->getNumRows();
    360       double * array = new double [numberRows];
    361       CoinZeroN(array,numberRows);
    362       int * which = new int [numberRows];
    363       int n=0;
    364       int base=numberLinks_*firstNonZero;
    365       for (j=firstNonZero;j<=lastNonZero;j++) {
    366         for (int k=0;k<numberLinks_;k++) {
    367           int iColumn = which_[base+k];
    368           double value = CoinMax(0.0,solution[iColumn]);
    369           if (value>integerTolerance&&upper[iColumn]) {
    370             value = CoinMin(value,upper[iColumn]);
    371             for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
    372               int iRow = row[j];
    373               double a = array[iRow];
    374               if (a) {
    375                 a += value*element[j];
    376                 if (!a)
    377                   a = 1.0e-100;
    378               } else {
    379                 which[n++]=iRow;
    380                 a=value*element[j];
    381                 assert (a);
    382               }
    383               array[iRow]=a;
    384             }
    385           }
    386         }
    387         base += numberLinks_;
    388       }
    389       base=numberLinks_*iWhere;
    390       for (int k=0;k<numberLinks_;k++) {
    391         int iColumn = which_[base+k];
    392         const double value = 1.0;
    393         for (int j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
    394           int iRow = row[j];
    395           double a = array[iRow];
    396           if (a) {
    397             a -= value*element[j];
    398             if (!a)
    399               a = 1.0e-100;
    400           } else {
    401             which[n++]=iRow;
    402             a=-value*element[j];
    403             assert (a);
    404           }
    405           array[iRow]=a;
    406         }
    407       }
    408       for (j=0;j<n;j++) {
    409         int iRow = which[j];
    410         // moving to point will increase row solution by this
    411         double distance = array[iRow];
    412         if (distance>1.0e-8) {
    413           if (distance+rowSolution[iRow]>rowUpper[iRow]+1.0e-8) {
    414             possible=false;
    415             break;
    416           }
    417         } else if (distance<-1.0e-8) {
    418           if (distance+rowSolution[iRow]<rowLower[iRow]-1.0e-8) {
    419             possible=false;
    420             break;
    421           }
    422         }
    423       }
    424       for (j=0;j<n;j++)
    425         array[which[j]]=0.0;
    426       delete [] array;
    427       delete [] which;
     357      double *array = new double[numberRows];
     358      CoinZeroN(array, numberRows);
     359      int *which = new int[numberRows];
     360      int n = 0;
     361      int base = numberLinks_ * firstNonZero;
     362      for (j = firstNonZero; j <= lastNonZero; j++) {
     363        for (int k = 0; k < numberLinks_; k++) {
     364          int iColumn = which_[base + k];
     365          double value = CoinMax(0.0, solution[iColumn]);
     366          if (value > integerTolerance && upper[iColumn]) {
     367            value = CoinMin(value, upper[iColumn]);
     368            for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     369              int iRow = row[j];
     370              double a = array[iRow];
     371              if (a) {
     372                a += value * element[j];
     373                if (!a)
     374                  a = 1.0e-100;
     375              } else {
     376                which[n++] = iRow;
     377                a = value * element[j];
     378                assert(a);
     379              }
     380              array[iRow] = a;
     381            }
     382          }
     383        }
     384        base += numberLinks_;
     385      }
     386      base = numberLinks_ * iWhere;
     387      for (int k = 0; k < numberLinks_; k++) {
     388        int iColumn = which_[base + k];
     389        const double value = 1.0;
     390        for (int j = columnStart[iColumn]; j < columnStart[iColumn] + columnLength[iColumn]; j++) {
     391          int iRow = row[j];
     392          double a = array[iRow];
     393          if (a) {
     394            a -= value * element[j];
     395            if (!a)
     396              a = 1.0e-100;
     397          } else {
     398            which[n++] = iRow;
     399            a = -value * element[j];
     400            assert(a);
     401          }
     402          array[iRow] = a;
     403        }
     404      }
     405      for (j = 0; j < n; j++) {
     406        int iRow = which[j];
     407        // moving to point will increase row solution by this
     408        double distance = array[iRow];
     409        if (distance > 1.0e-8) {
     410          if (distance + rowSolution[iRow] > rowUpper[iRow] + 1.0e-8) {
     411            possible = false;
     412            break;
     413          }
     414        } else if (distance < -1.0e-8) {
     415          if (distance + rowSolution[iRow] < rowLower[iRow] - 1.0e-8) {
     416            possible = false;
     417            break;
     418          }
     419        }
     420      }
     421      for (j = 0; j < n; j++)
     422        array[which[j]] = 0.0;
     423      delete[] array;
     424      delete[] which;
    428425      if (possible) {
    429         printf("possible feas region %d %d %d\n",firstNonZero,lastNonZero,iWhere);
    430         firstNonZero=iWhere;
    431         lastNonZero=iWhere;
     426        printf("possible feas region %d %d %d\n", firstNonZero, lastNonZero, iWhere);
     427        firstNonZero = iWhere;
     428        lastNonZero = iWhere;
    432429      }
    433430    }
    434431  }
    435432#else
    436   assert (lastNonZero-firstNonZero<sosType_) ;
     433  assert(lastNonZero - firstNonZero < sosType_);
    437434#endif
    438   base=0;
    439   for (j=0;j<firstNonZero;j++) {
    440     for (int k=0;k<numberLinks_;k++) {
    441       int iColumn = which_[base+k];
    442       solver->setColUpper(iColumn,0.0);
     435  base = 0;
     436  for (j = 0; j < firstNonZero; j++) {
     437    for (int k = 0; k < numberLinks_; k++) {
     438      int iColumn = which_[base + k];
     439      solver->setColUpper(iColumn, 0.0);
    443440    }
    444441    base += numberLinks_;
     
    446443  // skip
    447444  base += numberLinks_;
    448   for (j=lastNonZero+1;j<numberMembers_;j++) {
    449     for (int k=0;k<numberLinks_;k++) {
    450       int iColumn = which_[base+k];
    451       solver->setColUpper(iColumn,0.0);
     445  for (j = lastNonZero + 1; j < numberMembers_; j++) {
     446    for (int k = 0; k < numberLinks_; k++) {
     447      int iColumn = which_[base + k];
     448      solver->setColUpper(iColumn, 0.0);
    452449    }
    453450    base += numberLinks_;
     
    455452}
    456453
    457 
    458454// Creates a branching object
    459 CbcBranchingObject * 
    460 CbcLink::createCbcBranch(OsiSolverInterface * /*solver*/, const OsiBranchingInformation * /*info*/, int way) 
     455CbcBranchingObject *
     456CbcLink::createCbcBranch(OsiSolverInterface * /*solver*/, const OsiBranchingInformation * /*info*/, int way)
    461457{
    462458  int j;
    463   const double * solution = model_->testSolution();
    464   double integerTolerance =
    465       model_->getDblParam(CbcModel::CbcIntegerTolerance);
    466   OsiSolverInterface * solver = model_->solver();
    467   const double * upper = solver->getColUpper();
    468   int firstNonFixed=-1;
    469   int lastNonFixed=-1;
    470   int firstNonZero=-1;
     459  const double *solution = model_->testSolution();
     460  double integerTolerance = model_->getDblParam(CbcModel::CbcIntegerTolerance);
     461  OsiSolverInterface *solver = model_->solver();
     462  const double *upper = solver->getColUpper();
     463  int firstNonFixed = -1;
     464  int lastNonFixed = -1;
     465  int firstNonZero = -1;
    471466  int lastNonZero = -1;
    472467  double weight = 0.0;
    473   double sum =0.0;
    474   int base=0;
    475   for (j=0;j<numberMembers_;j++) {
    476     for (int k=0;k<numberLinks_;k++) {
    477       int iColumn = which_[base+k];
     468  double sum = 0.0;
     469  int base = 0;
     470  for (j = 0; j < numberMembers_; j++) {
     471    for (int k = 0; k < numberLinks_; k++) {
     472      int iColumn = which_[base + k];
    478473      if (upper[iColumn]) {
    479         double value = CoinMax(0.0,solution[iColumn]);
     474        double value = CoinMax(0.0, solution[iColumn]);
    480475        sum += value;
    481         if (firstNonFixed<0)
    482           firstNonFixed=j;
    483         lastNonFixed=j;
    484         if (value>integerTolerance) {
    485           weight += weights_[j]*value;
    486           if (firstNonZero<0)
    487             firstNonZero=j;
    488           lastNonZero=j;
     476        if (firstNonFixed < 0)
     477          firstNonFixed = j;
     478        lastNonFixed = j;
     479        if (value > integerTolerance) {
     480          weight += weights_[j] * value;
     481          if (firstNonZero < 0)
     482            firstNonZero = j;
     483          lastNonZero = j;
    489484        }
    490485      }
     
    492487    base += numberLinks_;
    493488  }
    494   assert (lastNonZero-firstNonZero>=sosType_) ;
     489  assert(lastNonZero - firstNonZero >= sosType_);
    495490  // find where to branch
    496   assert (sum>0.0);
     491  assert(sum > 0.0);
    497492  weight /= sum;
    498493  int iWhere;
    499   double separator=0.0;
    500   for (iWhere=firstNonZero;iWhere<lastNonZero;iWhere++)
    501     if (weight<weights_[iWhere+1])
     494  double separator = 0.0;
     495  for (iWhere = firstNonZero; iWhere < lastNonZero; iWhere++)
     496    if (weight < weights_[iWhere + 1])
    502497      break;
    503   if (sosType_==1) {
     498  if (sosType_ == 1) {
    504499    // SOS 1
    505     separator = 0.5 *(weights_[iWhere]+weights_[iWhere+1]);
     500    separator = 0.5 * (weights_[iWhere] + weights_[iWhere + 1]);
    506501  } else {
    507502    // SOS 2
    508     if (iWhere==firstNonFixed)
    509       iWhere++;;
    510     if (iWhere==lastNonFixed-1)
    511       iWhere = lastNonFixed-2;
    512     separator = weights_[iWhere+1];
     503    if (iWhere == firstNonFixed)
     504      iWhere++;
     505    ;
     506    if (iWhere == lastNonFixed - 1)
     507      iWhere = lastNonFixed - 2;
     508    separator = weights_[iWhere + 1];
    513509  }
    514510  // create object
    515   CbcBranchingObject * branch;
    516   branch = new CbcLinkBranchingObject(model_,this,way,separator);
     511  CbcBranchingObject *branch;
     512  branch = new CbcLinkBranchingObject(model_, this, way, separator);
    517513  return branch;
    518514}
    519515// Useful constructor
    520 CbcLinkBranchingObject::CbcLinkBranchingObject (CbcModel * model,
    521                                               const CbcLink * set,
    522                                               int way ,
    523                                               double separator)
    524   :CbcBranchingObject(model,set->id(),way,0.5)
     516CbcLinkBranchingObject::CbcLinkBranchingObject(CbcModel *model,
     517  const CbcLink *set,
     518  int way,
     519  double separator)
     520  : CbcBranchingObject(model, set->id(), way, 0.5)
    525521{
    526522  set_ = set;
     
    528524}
    529525
    530 // Copy constructor
    531 CbcLinkBranchingObject::CbcLinkBranchingObject ( const CbcLinkBranchingObject & rhs) :CbcBranchingObject(rhs)
    532 {
    533   set_=rhs.set_;
     526// Copy constructor
     527CbcLinkBranchingObject::CbcLinkBranchingObject(const CbcLinkBranchingObject &rhs)
     528  : CbcBranchingObject(rhs)
     529{
     530  set_ = rhs.set_;
    534531  separator_ = rhs.separator_;
    535532}
    536533
    537 // Assignment operator 
    538 CbcLinkBranchingObject & 
    539 CbcLinkBranchingObject::operator=( const CbcLinkBranchingObject& rhs)
     534// Assignment operator
     535CbcLinkBranchingObject &
     536CbcLinkBranchingObject::operator=(const CbcLinkBranchingObject &rhs)
    540537{
    541538  if (this != &rhs) {
    542539    CbcBranchingObject::operator=(rhs);
    543     set_=rhs.set_;
     540    set_ = rhs.set_;
    544541    separator_ = rhs.separator_;
    545542  }
    546543  return *this;
    547544}
    548 CbcBranchingObject * 
     545CbcBranchingObject *
    549546CbcLinkBranchingObject::clone() const
    550 { 
     547{
    551548  return (new CbcLinkBranchingObject(*this));
    552549}
    553550
    554 
    555 // Destructor
    556 CbcLinkBranchingObject::~CbcLinkBranchingObject ()
     551// Destructor
     552CbcLinkBranchingObject::~CbcLinkBranchingObject()
    557553{
    558554}
     
    563559  int numberMembers = set_->numberMembers();
    564560  int numberLinks = set_->numberLinks();
    565   const double * weights = set_->weights();
    566   const int * which = set_->which();
    567   OsiSolverInterface * solver = model_->solver();
     561  const double *weights = set_->weights();
     562  const int *which = set_->which();
     563  OsiSolverInterface *solver = model_->solver();
    568564  //const double * lower = solver->getColLower();
    569565  //const double * upper = solver->getColUpper();
    570566  // *** for way - up means fix all those in down section
    571   if (way_<0) {
     567  if (way_ < 0) {
    572568    int i;
    573     for ( i=0;i<numberMembers;i++) {
     569    for (i = 0; i < numberMembers; i++) {
    574570      if (weights[i] > separator_)
    575         break;
    576     }
    577     assert (i<numberMembers);
    578     int base=i*numberLinks;;
    579     for (;i<numberMembers;i++) {
    580       for (int k=0;k<numberLinks;k++) {
    581         int iColumn = which[base+k];
    582         solver->setColUpper(iColumn,0.0);
     571        break;
     572    }
     573    assert(i < numberMembers);
     574    int base = i * numberLinks;
     575    ;
     576    for (; i < numberMembers; i++) {
     577      for (int k = 0; k < numberLinks; k++) {
     578        int iColumn = which[base + k];
     579        solver->setColUpper(iColumn, 0.0);
    583580      }
    584581      base += numberLinks;
    585582    }
    586     way_=1;      // Swap direction
     583    way_ = 1; // Swap direction
    587584  } else {
    588585    int i;
    589     int base=0;
    590     for ( i=0;i<numberMembers;i++) {
     586    int base = 0;
     587    for (i = 0; i < numberMembers; i++) {
    591588      if (weights[i] >= separator_) {
    592         break;
     589        break;
    593590      } else {
    594         for (int k=0;k<numberLinks;k++) {
    595           int iColumn = which[base+k];
    596           solver->setColUpper(iColumn,0.0);
     591        for (int k = 0; k < numberLinks; k++) {
     592          int iColumn = which[base + k];
     593          solver->setColUpper(iColumn, 0.0);
    597594        }
    598595        base += numberLinks;
    599596      }
    600597    }
    601     assert (i<numberMembers);
    602     way_=-1;      // Swap direction
     598    assert(i < numberMembers);
     599    way_ = -1; // Swap direction
    603600  }
    604601  return 0.0;
    605602}
    606 // Print what would happen 
    607 void
    608 CbcLinkBranchingObject::print()
     603// Print what would happen
     604void CbcLinkBranchingObject::print()
    609605{
    610606  int numberMembers = set_->numberMembers();
    611607  int numberLinks = set_->numberLinks();
    612   const double * weights = set_->weights();
    613   const int * which = set_->which();
    614   OsiSolverInterface * solver = model_->solver();
    615   const double * upper = solver->getColUpper();
    616   int first=numberMembers;
    617   int last=-1;
    618   int numberFixed=0;
    619   int numberOther=0;
     608  const double *weights = set_->weights();
     609  const int *which = set_->which();
     610  OsiSolverInterface *solver = model_->solver();
     611  const double *upper = solver->getColUpper();
     612  int first = numberMembers;
     613  int last = -1;
     614  int numberFixed = 0;
     615  int numberOther = 0;
    620616  int i;
    621   int base=0;
    622   for ( i=0;i<numberMembers;i++) {
    623     for (int k=0;k<numberLinks;k++) {
    624       int iColumn = which[base+k];
     617  int base = 0;
     618  for (i = 0; i < numberMembers; i++) {
     619    for (int k = 0; k < numberLinks; k++) {
     620      int iColumn = which[base + k];
    625621      double bound = upper[iColumn];
    626622      if (bound) {
    627         first = CoinMin(first,i);
    628         last = CoinMax(last,i);
     623        first = CoinMin(first, i);
     624        last = CoinMax(last, i);
    629625      }
    630626    }
     
    632628  }
    633629  // *** for way - up means fix all those in down section
    634   base=0;
    635   if (way_<0) {
     630  base = 0;
     631  if (way_ < 0) {
    636632    printf("SOS Down");
    637     for ( i=0;i<numberMembers;i++) {
    638       if (weights[i] > separator_) 
    639         break;
    640       for (int k=0;k<numberLinks;k++) {
    641         int iColumn = which[base+k];
     633    for (i = 0; i < numberMembers; i++) {
     634      if (weights[i] > separator_)
     635        break;
     636      for (int k = 0; k < numberLinks; k++) {
     637        int iColumn = which[base + k];
    642638        double bound = upper[iColumn];
    643         if (bound)
     639        if (bound)
    644640          numberOther++;
    645641      }
    646642      base += numberLinks;
    647643    }
    648     assert (i<numberMembers);
    649     for (;i<numberMembers;i++) {
    650       for (int k=0;k<numberLinks;k++) {
    651         int iColumn = which[base+k];
     644    assert(i < numberMembers);
     645    for (; i < numberMembers; i++) {
     646      for (int k = 0; k < numberLinks; k++) {
     647        int iColumn = which[base + k];
    652648        double bound = upper[iColumn];
    653649        if (bound)
     
    658654  } else {
    659655    printf("SOS Up");
    660     for ( i=0;i<numberMembers;i++) {
     656    for (i = 0; i < numberMembers; i++) {
    661657      if (weights[i] >= separator_)
    662         break;
    663       for (int k=0;k<numberLinks;k++) {
    664         int iColumn = which[base+k];
     658        break;
     659      for (int k = 0; k < numberLinks; k++) {
     660        int iColumn = which[base + k];
    665661        double bound = upper[iColumn];
    666         if (bound)
     662        if (bound)
    667663          numberFixed++;
    668664      }
    669665      base += numberLinks;
    670666    }
    671     assert (i<numberMembers);
    672     for (;i<numberMembers;i++) {
    673       for (int k=0;k<numberLinks;k++) {
    674         int iColumn = which[base+k];
     667    assert(i < numberMembers);
     668    for (; i < numberMembers; i++) {
     669      for (int k = 0; k < numberLinks; k++) {
     670        int iColumn = which[base + k];
    675671        double bound = upper[iColumn];
    676672        if (bound)
     
    680676    }
    681677  }
    682   assert ((numberFixed%numberLinks)==0);
    683   assert ((numberOther%numberLinks)==0);
     678  assert((numberFixed % numberLinks) == 0);
     679  assert((numberOther % numberLinks) == 0);
    684680  printf(" - at %g, free range %d (%g) => %d (%g), %d would be fixed, %d other way\n",
    685          separator_,first,weights[first],last,weights[last],numberFixed/numberLinks,
    686          numberOther/numberLinks);
     681    separator_, first, weights[first], last, weights[last], numberFixed / numberLinks,
     682    numberOther / numberLinks);
    687683}
    688684/** Compare the \c this with \c brObj. \c this and \c brObj must be os the
     
    695691*/
    696692CbcRangeCompare
    697 CbcLinkBranchingObject::compareBranchingObject
    698 (const CbcBranchingObject* brObj, const bool replaceIfOverlap)
     693CbcLinkBranchingObject::compareBranchingObject(const CbcBranchingObject *brObj, const bool replaceIfOverlap)
    699694{
    700695  throw("must implement");
Note: See TracChangeset for help on using the changeset viewer.