Ignore:
Timestamp:
Jan 6, 2019 2:43:06 PM (4 months ago)
Author:
unxusr
Message:

formatting

File:
1 edited

Legend:

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

    r2235 r2385  
    2121// Default Constructor
    2222//-------------------------------------------------------------------
    23 ClpNonLinearCost::ClpNonLinearCost () :
    24      changeCost_(0.0),
    25      feasibleCost_(0.0),
    26      infeasibilityWeight_(-1.0),
    27      largestInfeasibility_(0.0),
    28      sumInfeasibilities_(0.0),
    29      averageTheta_(0.0),
    30      numberRows_(0),
    31      numberColumns_(0),
    32      start_(NULL),
    33      whichRange_(NULL),
    34      offset_(NULL),
    35      lower_(NULL),
    36      cost_(NULL),
    37      model_(NULL),
    38      infeasible_(NULL),
    39      numberInfeasibilities_(-1),
    40      status_(NULL),
    41      bound_(NULL),
    42      cost2_(NULL),
    43      method_(1),
    44      convex_(true),
    45      bothWays_(false)
    46 {
    47 
     23ClpNonLinearCost::ClpNonLinearCost()
     24  : changeCost_(0.0)
     25  , feasibleCost_(0.0)
     26  , infeasibilityWeight_(-1.0)
     27  , largestInfeasibility_(0.0)
     28  , sumInfeasibilities_(0.0)
     29  , averageTheta_(0.0)
     30  , numberRows_(0)
     31  , numberColumns_(0)
     32  , start_(NULL)
     33  , whichRange_(NULL)
     34  , offset_(NULL)
     35  , lower_(NULL)
     36  , cost_(NULL)
     37  , model_(NULL)
     38  , infeasible_(NULL)
     39  , numberInfeasibilities_(-1)
     40  , status_(NULL)
     41  , bound_(NULL)
     42  , cost2_(NULL)
     43  , method_(1)
     44  , convex_(true)
     45  , bothWays_(false)
     46{
    4847}
    4948//#define VALIDATE
    5049#ifdef VALIDATE
    51 static double * saveLowerV = NULL;
    52 static double * saveUpperV = NULL;
     50static double *saveLowerV = NULL;
     51static double *saveUpperV = NULL;
    5352#ifdef NDEBUG
    5453Validate sgould not be set if no debug
     
    5958   later may do dual analysis and even finding duplicate columns
    6059*/
    61 ClpNonLinearCost::ClpNonLinearCost ( ClpSimplex * model, int method)
    62 {
    63      method = 2;
    64      model_ = model;
    65      numberRows_ = model_->numberRows();
    66      //if (numberRows_==402) {
    67      //model_->setLogLevel(63);
    68      //model_->setMaximumIterations(30000);
    69      //}
    70      numberColumns_ = model_->numberColumns();
    71      // If gub then we need this extra
    72      int numberExtra = model_->numberExtraRows();
    73      if (numberExtra)
    74           method = 1;
    75      int numberTotal1 = numberRows_ + numberColumns_;
    76      int numberTotal = numberTotal1 + numberExtra;
    77      convex_ = true;
    78      bothWays_ = false;
    79      method_ = method;
    80      numberInfeasibilities_ = 0;
    81      changeCost_ = 0.0;
    82      feasibleCost_ = 0.0;
    83      infeasibilityWeight_ = -1.0;
    84      double * cost = model_->costRegion();
    85      // check if all 0
    86      int iSequence;
    87      bool allZero = true;
    88      for (iSequence = 0; iSequence < numberTotal1; iSequence++) {
    89           if (cost[iSequence]) {
    90                allZero = false;
    91                break;
    92           }
    93      }
    94      if (allZero&&model_->clpMatrix()->type()<15)
    95           model_->setInfeasibilityCost(1.0);
    96      double infeasibilityCost = model_->infeasibilityCost();
    97      sumInfeasibilities_ = 0.0;
    98      averageTheta_ = 0.0;
    99      largestInfeasibility_ = 0.0;
    100      // All arrays NULL to start
    101      status_ = NULL;
    102      bound_ = NULL;
    103      cost2_ = NULL;
    104      start_ = NULL;
    105      whichRange_ = NULL;
    106      offset_ = NULL;
    107      lower_ = NULL;
    108      cost_ = NULL;
    109      infeasible_ = NULL;
    110 
    111      double * upper = model_->upperRegion();
    112      double * lower = model_->lowerRegion();
    113 
    114      // See how we are storing things
    115      bool always4 = (model_->clpMatrix()->
    116                      generalExpanded(model_, 10, iSequence) != 0);
    117      if (always4)
    118           method_ = 1;
    119      if (CLP_METHOD1) {
    120           start_ = new int [numberTotal+1];
    121           whichRange_ = new int [numberTotal];
    122           offset_ = new int [numberTotal];
    123           memset(offset_, 0, numberTotal * sizeof(int));
    124 
    125 
    126           // First see how much space we need
    127           int put = 0;
    128 
    129           // For quadratic we need -inf,0,0,+inf
    130           for (iSequence = 0; iSequence < numberTotal1; iSequence++) {
    131                if (!always4) {
    132                     if (lower[iSequence] > -COIN_DBL_MAX)
    133                          put++;
    134                     if (upper[iSequence] < COIN_DBL_MAX)
    135                          put++;
    136                     put += 2;
    137                } else {
    138                     put += 4;
    139                }
    140           }
    141 
    142           // and for extra
    143           put += 4 * numberExtra;
     60ClpNonLinearCost::ClpNonLinearCost(ClpSimplex * model, int method)
     61{
     62  method = 2;
     63  model_ = model;
     64  numberRows_ = model_->numberRows();
     65  //if (numberRows_==402) {
     66  //model_->setLogLevel(63);
     67  //model_->setMaximumIterations(30000);
     68  //}
     69  numberColumns_ = model_->numberColumns();
     70  // If gub then we need this extra
     71  int numberExtra = model_->numberExtraRows();
     72  if (numberExtra)
     73    method = 1;
     74  int numberTotal1 = numberRows_ + numberColumns_;
     75  int numberTotal = numberTotal1 + numberExtra;
     76  convex_ = true;
     77  bothWays_ = false;
     78  method_ = method;
     79  numberInfeasibilities_ = 0;
     80  changeCost_ = 0.0;
     81  feasibleCost_ = 0.0;
     82  infeasibilityWeight_ = -1.0;
     83  double *cost = model_->costRegion();
     84  // check if all 0
     85  int iSequence;
     86  bool allZero = true;
     87  for (iSequence = 0; iSequence < numberTotal1; iSequence++) {
     88    if (cost[iSequence]) {
     89      allZero = false;
     90      break;
     91    }
     92  }
     93  if (allZero && model_->clpMatrix()->type() < 15)
     94    model_->setInfeasibilityCost(1.0);
     95  double infeasibilityCost = model_->infeasibilityCost();
     96  sumInfeasibilities_ = 0.0;
     97  averageTheta_ = 0.0;
     98  largestInfeasibility_ = 0.0;
     99  // All arrays NULL to start
     100  status_ = NULL;
     101  bound_ = NULL;
     102  cost2_ = NULL;
     103  start_ = NULL;
     104  whichRange_ = NULL;
     105  offset_ = NULL;
     106  lower_ = NULL;
     107  cost_ = NULL;
     108  infeasible_ = NULL;
     109
     110  double *upper = model_->upperRegion();
     111  double *lower = model_->lowerRegion();
     112
     113  // See how we are storing things
     114  bool always4 = (model_->clpMatrix()->generalExpanded(model_, 10, iSequence) != 0);
     115  if (always4)
     116    method_ = 1;
     117  if (CLP_METHOD1) {
     118    start_ = new int[numberTotal + 1];
     119    whichRange_ = new int[numberTotal];
     120    offset_ = new int[numberTotal];
     121    memset(offset_, 0, numberTotal * sizeof(int));
     122
     123    // First see how much space we need
     124    int put = 0;
     125
     126    // For quadratic we need -inf,0,0,+inf
     127    for (iSequence = 0; iSequence < numberTotal1; iSequence++) {
     128      if (!always4) {
     129        if (lower[iSequence] > -COIN_DBL_MAX)
     130          put++;
     131        if (upper[iSequence] < COIN_DBL_MAX)
     132          put++;
     133        put += 2;
     134      } else {
     135        put += 4;
     136      }
     137    }
     138
     139    // and for extra
     140    put += 4 * numberExtra;
    144141#ifndef NDEBUG
    145           int kPut = put;
    146 #endif
    147           lower_ = new double [put];
    148           cost_ = new double [put];
    149           infeasible_ = new unsigned int[(put+31)>>5];
    150           memset(infeasible_, 0, ((put + 31) >> 5)*sizeof(unsigned int));
    151 
    152           put = 0;
    153 
    154           start_[0] = 0;
    155 
    156           for (iSequence = 0; iSequence < numberTotal1; iSequence++) {
    157                if (!always4) {
    158                     if (lower[iSequence] > -COIN_DBL_MAX) {
    159                          lower_[put] = -COIN_DBL_MAX;
    160                          setInfeasible(put, true);
    161                          cost_[put++] = cost[iSequence] - infeasibilityCost;
    162                     }
    163                     whichRange_[iSequence] = put;
    164                     lower_[put] = lower[iSequence];
    165                     cost_[put++] = cost[iSequence];
    166                     lower_[put] = upper[iSequence];
    167                     cost_[put++] = cost[iSequence] + infeasibilityCost;
    168                     if (upper[iSequence] < COIN_DBL_MAX) {
    169                          lower_[put] = COIN_DBL_MAX;
    170                          setInfeasible(put - 1, true);
    171                          cost_[put++] = 1.0e50;
    172                     }
    173                } else {
    174                     lower_[put] = -COIN_DBL_MAX;
    175                     setInfeasible(put, true);
    176                     cost_[put++] = cost[iSequence] - infeasibilityCost;
    177                     whichRange_[iSequence] = put;
    178                     lower_[put] = lower[iSequence];
    179                     cost_[put++] = cost[iSequence];
    180                     lower_[put] = upper[iSequence];
    181                     cost_[put++] = cost[iSequence] + infeasibilityCost;
    182                     lower_[put] = COIN_DBL_MAX;
    183                     setInfeasible(put - 1, true);
    184                     cost_[put++] = 1.0e50;
    185                }
    186                start_[iSequence+1] = put;
    187           }
    188           for (; iSequence < numberTotal; iSequence++) {
    189 
    190                lower_[put] = -COIN_DBL_MAX;
    191                setInfeasible(put, true);
    192                put++;
    193                whichRange_[iSequence] = put;
    194                lower_[put] = 0.0;
    195                cost_[put++] = 0.0;
    196                lower_[put] = 0.0;
    197                cost_[put++] = 0.0;
    198                lower_[put] = COIN_DBL_MAX;
    199                setInfeasible(put - 1, true);
    200                cost_[put++] = 1.0e50;
    201                start_[iSequence+1] = put;
    202           }
    203           assert (put <= kPut);
    204      }
     142    int kPut = put;
     143#endif
     144    lower_ = new double[put];
     145    cost_ = new double[put];
     146    infeasible_ = new unsigned int[(put + 31) >> 5];
     147    memset(infeasible_, 0, ((put + 31) >> 5) * sizeof(unsigned int));
     148
     149    put = 0;
     150
     151    start_[0] = 0;
     152
     153    for (iSequence = 0; iSequence < numberTotal1; iSequence++) {
     154      if (!always4) {
     155        if (lower[iSequence] > -COIN_DBL_MAX) {
     156          lower_[put] = -COIN_DBL_MAX;
     157          setInfeasible(put, true);
     158          cost_[put++] = cost[iSequence] - infeasibilityCost;
     159        }
     160        whichRange_[iSequence] = put;
     161        lower_[put] = lower[iSequence];
     162        cost_[put++] = cost[iSequence];
     163        lower_[put] = upper[iSequence];
     164        cost_[put++] = cost[iSequence] + infeasibilityCost;
     165        if (upper[iSequence] < COIN_DBL_MAX) {
     166          lower_[put] = COIN_DBL_MAX;
     167          setInfeasible(put - 1, true);
     168          cost_[put++] = 1.0e50;
     169        }
     170      } else {
     171        lower_[put] = -COIN_DBL_MAX;
     172        setInfeasible(put, true);
     173        cost_[put++] = cost[iSequence] - infeasibilityCost;
     174        whichRange_[iSequence] = put;
     175        lower_[put] = lower[iSequence];
     176        cost_[put++] = cost[iSequence];
     177        lower_[put] = upper[iSequence];
     178        cost_[put++] = cost[iSequence] + infeasibilityCost;
     179        lower_[put] = COIN_DBL_MAX;
     180        setInfeasible(put - 1, true);
     181        cost_[put++] = 1.0e50;
     182      }
     183      start_[iSequence + 1] = put;
     184    }
     185    for (; iSequence < numberTotal; iSequence++) {
     186
     187      lower_[put] = -COIN_DBL_MAX;
     188      setInfeasible(put, true);
     189      put++;
     190      whichRange_[iSequence] = put;
     191      lower_[put] = 0.0;
     192      cost_[put++] = 0.0;
     193      lower_[put] = 0.0;
     194      cost_[put++] = 0.0;
     195      lower_[put] = COIN_DBL_MAX;
     196      setInfeasible(put - 1, true);
     197      cost_[put++] = 1.0e50;
     198      start_[iSequence + 1] = put;
     199    }
     200    assert(put <= kPut);
     201  }
    205202#ifdef FAST_CLPNON
    206      // See how we are storing things
    207      CoinAssert (model_->clpMatrix()->
    208                  generalExpanded(model_, 10, iSequence) == 0);
    209 #endif
    210      if (CLP_METHOD2) {
    211           assert (!numberExtra);
    212           bound_ = new double[numberTotal];
    213           cost2_ = new double[numberTotal];
    214           status_ = new unsigned char[numberTotal];
     203  // See how we are storing things
     204  CoinAssert(model_->clpMatrix()->generalExpanded(model_, 10, iSequence) == 0);
     205#endif
     206  if (CLP_METHOD2) {
     207    assert(!numberExtra);
     208    bound_ = new double[numberTotal];
     209    cost2_ = new double[numberTotal];
     210    status_ = new unsigned char[numberTotal];
    215211#ifdef VALIDATE
    216           delete [] saveLowerV;
    217           saveLowerV = CoinCopyOfArray(model_->lowerRegion(), numberTotal);
    218           delete [] saveUpperV;
    219           saveUpperV = CoinCopyOfArray(model_->upperRegion(), numberTotal);
    220 #endif
    221           for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    222                bound_[iSequence] = 0.0;
    223                cost2_[iSequence] = cost[iSequence];
    224                setInitialStatus(status_[iSequence]);
    225           }
    226      }
     212    delete[] saveLowerV;
     213    saveLowerV = CoinCopyOfArray(model_->lowerRegion(), numberTotal);
     214    delete[] saveUpperV;
     215    saveUpperV = CoinCopyOfArray(model_->upperRegion(), numberTotal);
     216#endif
     217    for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     218      bound_[iSequence] = 0.0;
     219      cost2_[iSequence] = cost[iSequence];
     220      setInitialStatus(status_[iSequence]);
     221    }
     222  }
    227223}
    228224#if 0
     
    281277#endif
    282278// Refresh one- assuming regions OK
    283 void
    284 ClpNonLinearCost::refresh(int iSequence)
    285 {
    286      double infeasibilityCost = model_->infeasibilityCost();
    287      double primalTolerance = model_->currentPrimalTolerance();
    288      double * cost = model_->costRegion();
    289      double * upper = model_->upperRegion();
    290      double * lower = model_->lowerRegion();
    291      double * solution = model_->solutionRegion();
    292      cost2_[iSequence] = cost[iSequence];
    293      double value = solution[iSequence];
    294      double lowerValue = lower[iSequence];
    295      double upperValue = upper[iSequence];
    296      if (value - upperValue <= primalTolerance) {
    297        if (value - lowerValue >= -primalTolerance) {
    298          // feasible
    299          status_[iSequence] = static_cast<unsigned char>(CLP_FEASIBLE | (CLP_SAME << 4));
    300          bound_[iSequence] = 0.0;
    301        } else {
    302          // below
    303          cost[iSequence] -= infeasibilityCost;
    304          status_[iSequence] = static_cast<unsigned char>(CLP_BELOW_LOWER | (CLP_SAME << 4));
    305          bound_[iSequence] = upperValue;
    306          upper[iSequence] = lowerValue;
    307          lower[iSequence] = -COIN_DBL_MAX;
    308        }
    309      } else {
    310        // above
    311        cost[iSequence] += infeasibilityCost;
    312        status_[iSequence] = static_cast<unsigned char>(CLP_ABOVE_UPPER | (CLP_SAME << 4));
    313        bound_[iSequence] = lowerValue;
    314        lower[iSequence] = upperValue;
    315        upper[iSequence] = COIN_DBL_MAX;
    316      }
    317      
     279void ClpNonLinearCost::refresh(int iSequence)
     280{
     281  double infeasibilityCost = model_->infeasibilityCost();
     282  double primalTolerance = model_->currentPrimalTolerance();
     283  double *cost = model_->costRegion();
     284  double *upper = model_->upperRegion();
     285  double *lower = model_->lowerRegion();
     286  double *solution = model_->solutionRegion();
     287  cost2_[iSequence] = cost[iSequence];
     288  double value = solution[iSequence];
     289  double lowerValue = lower[iSequence];
     290  double upperValue = upper[iSequence];
     291  if (value - upperValue <= primalTolerance) {
     292    if (value - lowerValue >= -primalTolerance) {
     293      // feasible
     294      status_[iSequence] = static_cast< unsigned char >(CLP_FEASIBLE | (CLP_SAME << 4));
     295      bound_[iSequence] = 0.0;
     296    } else {
     297      // below
     298      cost[iSequence] -= infeasibilityCost;
     299      status_[iSequence] = static_cast< unsigned char >(CLP_BELOW_LOWER | (CLP_SAME << 4));
     300      bound_[iSequence] = upperValue;
     301      upper[iSequence] = lowerValue;
     302      lower[iSequence] = -COIN_DBL_MAX;
     303    }
     304  } else {
     305    // above
     306    cost[iSequence] += infeasibilityCost;
     307    status_[iSequence] = static_cast< unsigned char >(CLP_ABOVE_UPPER | (CLP_SAME << 4));
     308    bound_[iSequence] = lowerValue;
     309    lower[iSequence] = upperValue;
     310    upper[iSequence] = COIN_DBL_MAX;
     311  }
    318312}
    319313// Refreshes costs always makes row costs zero
    320 void
    321 ClpNonLinearCost::refreshCosts(const double * columnCosts)
    322 {
    323      double * cost = model_->costRegion();
    324      // zero row costs
    325      memset(cost + numberColumns_, 0, numberRows_ * sizeof(double));
    326      // copy column costs
    327      CoinMemcpyN(columnCosts, numberColumns_, cost);
    328      if ((method_ & 1) != 0) {
    329           for (int iSequence = 0; iSequence < numberRows_ + numberColumns_; iSequence++) {
    330                int start = start_[iSequence];
    331                int end = start_[iSequence+1] - 1;
    332                double thisFeasibleCost = cost[iSequence];
    333                if (infeasible(start)) {
    334                     cost_[start] = thisFeasibleCost - infeasibilityWeight_;
    335                     cost_[start+1] = thisFeasibleCost;
    336                } else {
    337                     cost_[start] = thisFeasibleCost;
    338                }
    339                if (infeasible(end - 1)) {
    340                     cost_[end-1] = thisFeasibleCost + infeasibilityWeight_;
    341                }
    342           }
    343      }
    344      if (CLP_METHOD2) {
    345           for (int iSequence = 0; iSequence < numberRows_ + numberColumns_; iSequence++) {
    346                cost2_[iSequence] = cost[iSequence];
    347           }
    348      }
    349 }
    350 ClpNonLinearCost::ClpNonLinearCost(ClpSimplex * model, const int * starts,
    351                                    const double * lowerNon, const double * costNon)
     314void ClpNonLinearCost::refreshCosts(const double *columnCosts)
     315{
     316  double *cost = model_->costRegion();
     317  // zero row costs
     318  memset(cost + numberColumns_, 0, numberRows_ * sizeof(double));
     319  // copy column costs
     320  CoinMemcpyN(columnCosts, numberColumns_, cost);
     321  if ((method_ & 1) != 0) {
     322    for (int iSequence = 0; iSequence < numberRows_ + numberColumns_; iSequence++) {
     323      int start = start_[iSequence];
     324      int end = start_[iSequence + 1] - 1;
     325      double thisFeasibleCost = cost[iSequence];
     326      if (infeasible(start)) {
     327        cost_[start] = thisFeasibleCost - infeasibilityWeight_;
     328        cost_[start + 1] = thisFeasibleCost;
     329      } else {
     330        cost_[start] = thisFeasibleCost;
     331      }
     332      if (infeasible(end - 1)) {
     333        cost_[end - 1] = thisFeasibleCost + infeasibilityWeight_;
     334      }
     335    }
     336  }
     337  if (CLP_METHOD2) {
     338    for (int iSequence = 0; iSequence < numberRows_ + numberColumns_; iSequence++) {
     339      cost2_[iSequence] = cost[iSequence];
     340    }
     341  }
     342}
     343ClpNonLinearCost::ClpNonLinearCost(ClpSimplex *model, const int *starts,
     344  const double *lowerNon, const double *costNon)
    352345{
    353346#ifndef FAST_CLPNON
    354      // what about scaling? - only try without it initially
    355      assert(!model->scalingFlag());
    356      model_ = model;
    357      numberRows_ = model_->numberRows();
    358      numberColumns_ = model_->numberColumns();
    359      int numberTotal = numberRows_ + numberColumns_;
    360      convex_ = true;
    361      bothWays_ = true;
    362      start_ = new int [numberTotal+1];
    363      whichRange_ = new int [numberTotal];
    364      offset_ = new int [numberTotal];
    365      memset(offset_, 0, numberTotal * sizeof(int));
    366 
    367      double whichWay = model_->optimizationDirection();
    368      COIN_DETAIL_PRINT(printf("Direction %g\n", whichWay));
    369 
    370      numberInfeasibilities_ = 0;
    371      changeCost_ = 0.0;
    372      feasibleCost_ = 0.0;
    373      double infeasibilityCost = model_->infeasibilityCost();
    374      infeasibilityWeight_ = infeasibilityCost;;
    375      largestInfeasibility_ = 0.0;
    376      sumInfeasibilities_ = 0.0;
    377 
    378      int iSequence;
    379      assert (!model_->rowObjective());
    380      double * cost = model_->objective();
    381 
    382      // First see how much space we need
    383      // Also set up feasible bounds
    384      int put = starts[numberColumns_];
    385 
    386      double * columnUpper = model_->columnUpper();
    387      double * columnLower = model_->columnLower();
    388      for (iSequence = 0; iSequence < numberColumns_; iSequence++) {
    389           if (columnLower[iSequence] > -1.0e20)
    390                put++;
    391           if (columnUpper[iSequence] < 1.0e20)
    392                put++;
    393      }
    394 
    395      double * rowUpper = model_->rowUpper();
    396      double * rowLower = model_->rowLower();
    397      for (iSequence = 0; iSequence < numberRows_; iSequence++) {
    398           if (rowLower[iSequence] > -1.0e20)
    399                put++;
    400           if (rowUpper[iSequence] < 1.0e20)
    401                put++;
    402           put += 2;
    403      }
    404      lower_ = new double [put];
    405      cost_ = new double [put];
    406      infeasible_ = new unsigned int[(put+31)>>5];
    407      memset(infeasible_, 0, ((put + 31) >> 5)*sizeof(unsigned int));
    408 
    409      // now fill in
    410      put = 0;
    411 
    412      start_[0] = 0;
    413      for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    414           lower_[put] = -COIN_DBL_MAX;
    415           whichRange_[iSequence] = put + 1;
    416           double thisCost;
    417           double lowerValue;
    418           double upperValue;
    419           if (iSequence >= numberColumns_) {
    420                // rows
    421                lowerValue = rowLower[iSequence-numberColumns_];
    422                upperValue = rowUpper[iSequence-numberColumns_];
    423                if (lowerValue > -1.0e30) {
    424                     setInfeasible(put, true);
    425                     cost_[put++] = -infeasibilityCost;
    426                     lower_[put] = lowerValue;
    427                }
    428                cost_[put++] = 0.0;
    429                thisCost = 0.0;
    430           } else {
    431                // columns - move costs and see if convex
    432                lowerValue = columnLower[iSequence];
    433                upperValue = columnUpper[iSequence];
    434                if (lowerValue > -1.0e30) {
    435                     setInfeasible(put, true);
    436                     cost_[put++] = whichWay * cost[iSequence] - infeasibilityCost;
    437                     lower_[put] = lowerValue;
    438                }
    439                int iIndex = starts[iSequence];
    440                int end = starts[iSequence+1];
    441                assert (fabs(columnLower[iSequence] - lowerNon[iIndex]) < 1.0e-8);
    442                thisCost = -COIN_DBL_MAX;
    443                for (; iIndex < end; iIndex++) {
    444                     if (lowerNon[iIndex] < columnUpper[iSequence] - 1.0e-8) {
    445                          lower_[put] = lowerNon[iIndex];
    446                          cost_[put++] = whichWay * costNon[iIndex];
    447                          // check convexity
    448                          if (whichWay * costNon[iIndex] < thisCost - 1.0e-12)
    449                               convex_ = false;
    450                          thisCost = whichWay * costNon[iIndex];
    451                     } else {
    452                          break;
    453                     }
    454                }
    455           }
    456           lower_[put] = upperValue;
    457           setInfeasible(put, true);
    458           cost_[put++] = thisCost + infeasibilityCost;
    459           if (upperValue < 1.0e20) {
    460                lower_[put] = COIN_DBL_MAX;
    461                cost_[put++] = 1.0e50;
    462           }
    463           int iFirst = start_[iSequence];
    464           if (lower_[iFirst] != -COIN_DBL_MAX) {
    465                setInfeasible(iFirst, true);
    466                whichRange_[iSequence] = iFirst + 1;
    467           } else {
    468                whichRange_[iSequence] = iFirst;
    469           }
    470           start_[iSequence+1] = put;
    471      }
    472      // can't handle non-convex at present
    473      assert(convex_);
    474      status_ = NULL;
    475      bound_ = NULL;
    476      cost2_ = NULL;
    477      method_ = 1;
     347  // what about scaling? - only try without it initially
     348  assert(!model->scalingFlag());
     349  model_ = model;
     350  numberRows_ = model_->numberRows();
     351  numberColumns_ = model_->numberColumns();
     352  int numberTotal = numberRows_ + numberColumns_;
     353  convex_ = true;
     354  bothWays_ = true;
     355  start_ = new int[numberTotal + 1];
     356  whichRange_ = new int[numberTotal];
     357  offset_ = new int[numberTotal];
     358  memset(offset_, 0, numberTotal * sizeof(int));
     359
     360  double whichWay = model_->optimizationDirection();
     361  COIN_DETAIL_PRINT(printf("Direction %g\n", whichWay));
     362
     363  numberInfeasibilities_ = 0;
     364  changeCost_ = 0.0;
     365  feasibleCost_ = 0.0;
     366  double infeasibilityCost = model_->infeasibilityCost();
     367  infeasibilityWeight_ = infeasibilityCost;
     368  ;
     369  largestInfeasibility_ = 0.0;
     370  sumInfeasibilities_ = 0.0;
     371
     372  int iSequence;
     373  assert(!model_->rowObjective());
     374  double *cost = model_->objective();
     375
     376  // First see how much space we need
     377  // Also set up feasible bounds
     378  int put = starts[numberColumns_];
     379
     380  double *columnUpper = model_->columnUpper();
     381  double *columnLower = model_->columnLower();
     382  for (iSequence = 0; iSequence < numberColumns_; iSequence++) {
     383    if (columnLower[iSequence] > -1.0e20)
     384      put++;
     385    if (columnUpper[iSequence] < 1.0e20)
     386      put++;
     387  }
     388
     389  double *rowUpper = model_->rowUpper();
     390  double *rowLower = model_->rowLower();
     391  for (iSequence = 0; iSequence < numberRows_; iSequence++) {
     392    if (rowLower[iSequence] > -1.0e20)
     393      put++;
     394    if (rowUpper[iSequence] < 1.0e20)
     395      put++;
     396    put += 2;
     397  }
     398  lower_ = new double[put];
     399  cost_ = new double[put];
     400  infeasible_ = new unsigned int[(put + 31) >> 5];
     401  memset(infeasible_, 0, ((put + 31) >> 5) * sizeof(unsigned int));
     402
     403  // now fill in
     404  put = 0;
     405
     406  start_[0] = 0;
     407  for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     408    lower_[put] = -COIN_DBL_MAX;
     409    whichRange_[iSequence] = put + 1;
     410    double thisCost;
     411    double lowerValue;
     412    double upperValue;
     413    if (iSequence >= numberColumns_) {
     414      // rows
     415      lowerValue = rowLower[iSequence - numberColumns_];
     416      upperValue = rowUpper[iSequence - numberColumns_];
     417      if (lowerValue > -1.0e30) {
     418        setInfeasible(put, true);
     419        cost_[put++] = -infeasibilityCost;
     420        lower_[put] = lowerValue;
     421      }
     422      cost_[put++] = 0.0;
     423      thisCost = 0.0;
     424    } else {
     425      // columns - move costs and see if convex
     426      lowerValue = columnLower[iSequence];
     427      upperValue = columnUpper[iSequence];
     428      if (lowerValue > -1.0e30) {
     429        setInfeasible(put, true);
     430        cost_[put++] = whichWay * cost[iSequence] - infeasibilityCost;
     431        lower_[put] = lowerValue;
     432      }
     433      int iIndex = starts[iSequence];
     434      int end = starts[iSequence + 1];
     435      assert(fabs(columnLower[iSequence] - lowerNon[iIndex]) < 1.0e-8);
     436      thisCost = -COIN_DBL_MAX;
     437      for (; iIndex < end; iIndex++) {
     438        if (lowerNon[iIndex] < columnUpper[iSequence] - 1.0e-8) {
     439          lower_[put] = lowerNon[iIndex];
     440          cost_[put++] = whichWay * costNon[iIndex];
     441          // check convexity
     442          if (whichWay * costNon[iIndex] < thisCost - 1.0e-12)
     443            convex_ = false;
     444          thisCost = whichWay * costNon[iIndex];
     445        } else {
     446          break;
     447        }
     448      }
     449    }
     450    lower_[put] = upperValue;
     451    setInfeasible(put, true);
     452    cost_[put++] = thisCost + infeasibilityCost;
     453    if (upperValue < 1.0e20) {
     454      lower_[put] = COIN_DBL_MAX;
     455      cost_[put++] = 1.0e50;
     456    }
     457    int iFirst = start_[iSequence];
     458    if (lower_[iFirst] != -COIN_DBL_MAX) {
     459      setInfeasible(iFirst, true);
     460      whichRange_[iSequence] = iFirst + 1;
     461    } else {
     462      whichRange_[iSequence] = iFirst;
     463    }
     464    start_[iSequence + 1] = put;
     465  }
     466  // can't handle non-convex at present
     467  assert(convex_);
     468  status_ = NULL;
     469  bound_ = NULL;
     470  cost2_ = NULL;
     471  method_ = 1;
    478472#else
    479      printf("recompile ClpNonLinearCost without FAST_CLPNON\n");
    480      abort();
     473  printf("recompile ClpNonLinearCost without FAST_CLPNON\n");
     474  abort();
    481475#endif
    482476}
     
    484478// Copy constructor
    485479//-------------------------------------------------------------------
    486 ClpNonLinearCost::ClpNonLinearCost (const ClpNonLinearCost & rhs) :
    487      changeCost_(0.0),
    488      feasibleCost_(0.0),
    489      infeasibilityWeight_(-1.0),
    490      largestInfeasibility_(0.0),
    491      sumInfeasibilities_(0.0),
    492      averageTheta_(0.0),
    493      numberRows_(rhs.numberRows_),
    494      numberColumns_(rhs.numberColumns_),
    495      start_(NULL),
    496      whichRange_(NULL),
    497      offset_(NULL),
    498      lower_(NULL),
    499      cost_(NULL),
    500      model_(NULL),
    501      infeasible_(NULL),
    502      numberInfeasibilities_(-1),
    503      status_(NULL),
    504      bound_(NULL),
    505      cost2_(NULL),
    506      method_(rhs.method_),
    507      convex_(true),
    508     bothWays_(rhs.bothWays_)
    509 {
    510      if (numberRows_) {
    511           int numberTotal = numberRows_ + numberColumns_;
    512           model_ = rhs.model_;
    513           numberInfeasibilities_ = rhs.numberInfeasibilities_;
    514           changeCost_ = rhs.changeCost_;
    515           feasibleCost_ = rhs.feasibleCost_;
    516           infeasibilityWeight_ = rhs.infeasibilityWeight_;
    517           largestInfeasibility_ = rhs.largestInfeasibility_;
    518           sumInfeasibilities_ = rhs.sumInfeasibilities_;
    519           averageTheta_ = rhs.averageTheta_;
    520           convex_ = rhs.convex_;
    521           if (CLP_METHOD1) {
    522                start_ = new int [numberTotal+1];
    523                CoinMemcpyN(rhs.start_, (numberTotal + 1), start_);
    524                whichRange_ = new int [numberTotal];
    525                CoinMemcpyN(rhs.whichRange_, numberTotal, whichRange_);
    526                offset_ = new int [numberTotal];
    527                CoinMemcpyN(rhs.offset_, numberTotal, offset_);
    528                int numberEntries = start_[numberTotal];
    529                lower_ = new double [numberEntries];
    530                CoinMemcpyN(rhs.lower_, numberEntries, lower_);
    531                cost_ = new double [numberEntries];
    532                CoinMemcpyN(rhs.cost_, numberEntries, cost_);
    533                infeasible_ = new unsigned int[(numberEntries+31)>>5];
    534                CoinMemcpyN(rhs.infeasible_, ((numberEntries + 31) >> 5), infeasible_);
    535           }
    536           if (CLP_METHOD2) {
    537                bound_ = CoinCopyOfArray(rhs.bound_, numberTotal);
    538                cost2_ = CoinCopyOfArray(rhs.cost2_, numberTotal);
    539                status_ = CoinCopyOfArray(rhs.status_, numberTotal);
    540           }
    541      }
     480ClpNonLinearCost::ClpNonLinearCost(const ClpNonLinearCost &rhs)
     481  : changeCost_(0.0)
     482  , feasibleCost_(0.0)
     483  , infeasibilityWeight_(-1.0)
     484  , largestInfeasibility_(0.0)
     485  , sumInfeasibilities_(0.0)
     486  , averageTheta_(0.0)
     487  , numberRows_(rhs.numberRows_)
     488  , numberColumns_(rhs.numberColumns_)
     489  , start_(NULL)
     490  , whichRange_(NULL)
     491  , offset_(NULL)
     492  , lower_(NULL)
     493  , cost_(NULL)
     494  , model_(NULL)
     495  , infeasible_(NULL)
     496  , numberInfeasibilities_(-1)
     497  , status_(NULL)
     498  , bound_(NULL)
     499  , cost2_(NULL)
     500  , method_(rhs.method_)
     501  , convex_(true)
     502  , bothWays_(rhs.bothWays_)
     503{
     504  if (numberRows_) {
     505    int numberTotal = numberRows_ + numberColumns_;
     506    model_ = rhs.model_;
     507    numberInfeasibilities_ = rhs.numberInfeasibilities_;
     508    changeCost_ = rhs.changeCost_;
     509    feasibleCost_ = rhs.feasibleCost_;
     510    infeasibilityWeight_ = rhs.infeasibilityWeight_;
     511    largestInfeasibility_ = rhs.largestInfeasibility_;
     512    sumInfeasibilities_ = rhs.sumInfeasibilities_;
     513    averageTheta_ = rhs.averageTheta_;
     514    convex_ = rhs.convex_;
     515    if (CLP_METHOD1) {
     516      start_ = new int[numberTotal + 1];
     517      CoinMemcpyN(rhs.start_, (numberTotal + 1), start_);
     518      whichRange_ = new int[numberTotal];
     519      CoinMemcpyN(rhs.whichRange_, numberTotal, whichRange_);
     520      offset_ = new int[numberTotal];
     521      CoinMemcpyN(rhs.offset_, numberTotal, offset_);
     522      int numberEntries = start_[numberTotal];
     523      lower_ = new double[numberEntries];
     524      CoinMemcpyN(rhs.lower_, numberEntries, lower_);
     525      cost_ = new double[numberEntries];
     526      CoinMemcpyN(rhs.cost_, numberEntries, cost_);
     527      infeasible_ = new unsigned int[(numberEntries + 31) >> 5];
     528      CoinMemcpyN(rhs.infeasible_, ((numberEntries + 31) >> 5), infeasible_);
     529    }
     530    if (CLP_METHOD2) {
     531      bound_ = CoinCopyOfArray(rhs.bound_, numberTotal);
     532      cost2_ = CoinCopyOfArray(rhs.cost2_, numberTotal);
     533      status_ = CoinCopyOfArray(rhs.status_, numberTotal);
     534    }
     535  }
    542536}
    543537
     
    545539// Destructor
    546540//-------------------------------------------------------------------
    547 ClpNonLinearCost::~ClpNonLinearCost ()
    548 {
    549      delete [] start_;
    550      delete [] whichRange_;
    551      delete [] offset_;
    552      delete [] lower_;
    553      delete [] cost_;
    554      delete [] infeasible_;
    555      delete [] status_;
    556      delete [] bound_;
    557      delete [] cost2_;
     541ClpNonLinearCost::~ClpNonLinearCost()
     542{
     543  delete[] start_;
     544  delete[] whichRange_;
     545  delete[] offset_;
     546  delete[] lower_;
     547  delete[] cost_;
     548  delete[] infeasible_;
     549  delete[] status_;
     550  delete[] bound_;
     551  delete[] cost2_;
    558552}
    559553
     
    562556//-------------------------------------------------------------------
    563557ClpNonLinearCost &
    564 ClpNonLinearCost::operator=(const ClpNonLinearCost& rhs)
    565 {
    566      if (this != &rhs) {
    567           numberRows_ = rhs.numberRows_;
    568           numberColumns_ = rhs.numberColumns_;
    569           delete [] start_;
    570           delete [] whichRange_;
    571           delete [] offset_;
    572           delete [] lower_;
    573           delete [] cost_;
    574           delete [] infeasible_;
    575           delete [] status_;
    576           delete [] bound_;
    577           delete [] cost2_;
    578           start_ = NULL;
    579           whichRange_ = NULL;
    580           lower_ = NULL;
    581           cost_ = NULL;
    582           infeasible_ = NULL;
    583           status_ = NULL;
    584           bound_ = NULL;
    585           cost2_ = NULL;
    586           method_ = rhs.method_;
    587           if (numberRows_) {
    588                int numberTotal = numberRows_ + numberColumns_;
    589                if (CLP_METHOD1) {
    590                     start_ = new int [numberTotal+1];
    591                     CoinMemcpyN(rhs.start_, (numberTotal + 1), start_);
    592                     whichRange_ = new int [numberTotal];
    593                     CoinMemcpyN(rhs.whichRange_, numberTotal, whichRange_);
    594                     offset_ = new int [numberTotal];
    595                     CoinMemcpyN(rhs.offset_, numberTotal, offset_);
    596                     int numberEntries = start_[numberTotal];
    597                     lower_ = new double [numberEntries];
    598                     CoinMemcpyN(rhs.lower_, numberEntries, lower_);
    599                     cost_ = new double [numberEntries];
    600                     CoinMemcpyN(rhs.cost_, numberEntries, cost_);
    601                     infeasible_ = new unsigned int[(numberEntries+31)>>5];
    602                     CoinMemcpyN(rhs.infeasible_, ((numberEntries + 31) >> 5), infeasible_);
    603                }
    604                if (CLP_METHOD2) {
    605                     bound_ = CoinCopyOfArray(rhs.bound_, numberTotal);
    606                     cost2_ = CoinCopyOfArray(rhs.cost2_, numberTotal);
    607                     status_ = CoinCopyOfArray(rhs.status_, numberTotal);
    608                }
    609           }
    610           model_ = rhs.model_;
    611           numberInfeasibilities_ = rhs.numberInfeasibilities_;
    612           changeCost_ = rhs.changeCost_;
    613           feasibleCost_ = rhs.feasibleCost_;
    614           infeasibilityWeight_ = rhs.infeasibilityWeight_;
    615           largestInfeasibility_ = rhs.largestInfeasibility_;
    616           sumInfeasibilities_ = rhs.sumInfeasibilities_;
    617           averageTheta_ = rhs.averageTheta_;
    618           convex_ = rhs.convex_;
    619           bothWays_ = rhs.bothWays_;
    620      }
    621      return *this;
     558ClpNonLinearCost::operator=(const ClpNonLinearCost &rhs)
     559{
     560  if (this != &rhs) {
     561    numberRows_ = rhs.numberRows_;
     562    numberColumns_ = rhs.numberColumns_;
     563    delete[] start_;
     564    delete[] whichRange_;
     565    delete[] offset_;
     566    delete[] lower_;
     567    delete[] cost_;
     568    delete[] infeasible_;
     569    delete[] status_;
     570    delete[] bound_;
     571    delete[] cost2_;
     572    start_ = NULL;
     573    whichRange_ = NULL;
     574    lower_ = NULL;
     575    cost_ = NULL;
     576    infeasible_ = NULL;
     577    status_ = NULL;
     578    bound_ = NULL;
     579    cost2_ = NULL;
     580    method_ = rhs.method_;
     581    if (numberRows_) {
     582      int numberTotal = numberRows_ + numberColumns_;
     583      if (CLP_METHOD1) {
     584        start_ = new int[numberTotal + 1];
     585        CoinMemcpyN(rhs.start_, (numberTotal + 1), start_);
     586        whichRange_ = new int[numberTotal];
     587        CoinMemcpyN(rhs.whichRange_, numberTotal, whichRange_);
     588        offset_ = new int[numberTotal];
     589        CoinMemcpyN(rhs.offset_, numberTotal, offset_);
     590        int numberEntries = start_[numberTotal];
     591        lower_ = new double[numberEntries];
     592        CoinMemcpyN(rhs.lower_, numberEntries, lower_);
     593        cost_ = new double[numberEntries];
     594        CoinMemcpyN(rhs.cost_, numberEntries, cost_);
     595        infeasible_ = new unsigned int[(numberEntries + 31) >> 5];
     596        CoinMemcpyN(rhs.infeasible_, ((numberEntries + 31) >> 5), infeasible_);
     597      }
     598      if (CLP_METHOD2) {
     599        bound_ = CoinCopyOfArray(rhs.bound_, numberTotal);
     600        cost2_ = CoinCopyOfArray(rhs.cost2_, numberTotal);
     601        status_ = CoinCopyOfArray(rhs.status_, numberTotal);
     602      }
     603    }
     604    model_ = rhs.model_;
     605    numberInfeasibilities_ = rhs.numberInfeasibilities_;
     606    changeCost_ = rhs.changeCost_;
     607    feasibleCost_ = rhs.feasibleCost_;
     608    infeasibilityWeight_ = rhs.infeasibilityWeight_;
     609    largestInfeasibility_ = rhs.largestInfeasibility_;
     610    sumInfeasibilities_ = rhs.sumInfeasibilities_;
     611    averageTheta_ = rhs.averageTheta_;
     612    convex_ = rhs.convex_;
     613    bothWays_ = rhs.bothWays_;
     614  }
     615  return *this;
    622616}
    623617// Changes infeasible costs and computes number and cost of infeas
     
    625619// We will also need a 2 bit per variable array for some
    626620// purpose which will come to me later
    627 void
    628 ClpNonLinearCost::checkInfeasibilities(double oldTolerance)
    629 {
    630      numberInfeasibilities_ = 0;
    631      double infeasibilityCost = model_->infeasibilityCost();
    632      changeCost_ = 0.0;
    633      largestInfeasibility_ = 0.0;
    634      sumInfeasibilities_ = 0.0;
    635      double primalTolerance = model_->currentPrimalTolerance();
    636      int iSequence;
    637      double * solution = model_->solutionRegion();
    638      double * upper = model_->upperRegion();
    639      double * lower = model_->lowerRegion();
    640      double * cost = model_->costRegion();
    641      bool toNearest = oldTolerance <= 0.0;
    642      feasibleCost_ = 0.0;
    643      //bool checkCosts = (infeasibilityWeight_ != infeasibilityCost);
    644      infeasibilityWeight_ = infeasibilityCost;
    645      int numberTotal = numberColumns_ + numberRows_;
    646      //#define NONLIN_DEBUG
     621void ClpNonLinearCost::checkInfeasibilities(double oldTolerance)
     622{
     623  numberInfeasibilities_ = 0;
     624  double infeasibilityCost = model_->infeasibilityCost();
     625  changeCost_ = 0.0;
     626  largestInfeasibility_ = 0.0;
     627  sumInfeasibilities_ = 0.0;
     628  double primalTolerance = model_->currentPrimalTolerance();
     629  int iSequence;
     630  double *solution = model_->solutionRegion();
     631  double *upper = model_->upperRegion();
     632  double *lower = model_->lowerRegion();
     633  double *cost = model_->costRegion();
     634  bool toNearest = oldTolerance <= 0.0;
     635  feasibleCost_ = 0.0;
     636  //bool checkCosts = (infeasibilityWeight_ != infeasibilityCost);
     637  infeasibilityWeight_ = infeasibilityCost;
     638  int numberTotal = numberColumns_ + numberRows_;
     639  //#define NONLIN_DEBUG
    647640#ifdef NONLIN_DEBUG
    648      double * saveSolution = NULL;
    649      double * saveLower = NULL;
    650      double * saveUpper = NULL;
    651      unsigned char * saveStatus = NULL;
    652      if (method_ == 3) {
    653           // Save solution as we will be checking
    654           saveSolution = CoinCopyOfArray(solution, numberTotal);
    655           saveLower = CoinCopyOfArray(lower, numberTotal);
    656           saveUpper = CoinCopyOfArray(upper, numberTotal);
    657           saveStatus = CoinCopyOfArray(model_->statusArray(), numberTotal);
    658      }
     641  double *saveSolution = NULL;
     642  double *saveLower = NULL;
     643  double *saveUpper = NULL;
     644  unsigned char *saveStatus = NULL;
     645  if (method_ == 3) {
     646    // Save solution as we will be checking
     647    saveSolution = CoinCopyOfArray(solution, numberTotal);
     648    saveLower = CoinCopyOfArray(lower, numberTotal);
     649    saveUpper = CoinCopyOfArray(upper, numberTotal);
     650    saveStatus = CoinCopyOfArray(model_->statusArray(), numberTotal);
     651  }
    659652#else
    660      assert (method_ != 3);
    661 #endif
    662      if (CLP_METHOD1) {
    663           // nonbasic should be at a valid bound
    664           for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    665                double lowerValue;
    666                double upperValue;
    667                double value = solution[iSequence];
    668                int iRange;
    669                // get correct place
    670                int start = start_[iSequence];
    671                int end = start_[iSequence+1] - 1;
    672                // correct costs for this infeasibility weight
    673                // If free then true cost will be first
    674                double thisFeasibleCost = cost_[start];
    675                if (infeasible(start)) {
    676                     thisFeasibleCost = cost_[start+1];
    677                     cost_[start] = thisFeasibleCost - infeasibilityCost;
    678                }
    679                if (infeasible(end - 1)) {
    680                     thisFeasibleCost = cost_[end-2];
    681                     cost_[end-1] = thisFeasibleCost + infeasibilityCost;
    682                }
    683                for (iRange = start; iRange < end; iRange++) {
    684                     if (value < lower_[iRange+1] + primalTolerance) {
    685                          // put in better range if infeasible
    686                          if (value >= lower_[iRange+1] - primalTolerance && infeasible(iRange) && iRange == start)
    687                               iRange++;
    688                          whichRange_[iSequence] = iRange;
    689                          break;
    690                     }
    691                }
    692                assert(iRange < end);
    693                lowerValue = lower_[iRange];
    694                upperValue = lower_[iRange+1];
    695                ClpSimplex::Status status = model_->getStatus(iSequence);
    696                if (upperValue == lowerValue && status != ClpSimplex::isFixed) {
    697                     if (status != ClpSimplex::basic) {
    698                          model_->setStatus(iSequence, ClpSimplex::isFixed);
    699                          status = ClpSimplex::isFixed;
    700                     }
    701                }
    702                //#define PRINT_DETAIL7 2
    703 #if PRINT_DETAIL7>1
    704                printf("NL %d sol %g bounds %g %g\n",
    705                       iSequence,solution[iSequence],
    706                       lowerValue,upperValue);
    707 #endif
    708                switch(status) {
    709 
    710                case ClpSimplex::basic:
    711                case ClpSimplex::superBasic:
    712                     // iRange is in correct place
    713                     // slot in here
    714                     if (infeasible(iRange)) {
    715                          if (lower_[iRange] < -1.0e50) {
    716                               //cost_[iRange] = cost_[iRange+1]-infeasibilityCost;
    717                               // possibly below
    718                               lowerValue = lower_[iRange+1];
    719                               if (value - lowerValue < -primalTolerance) {
    720                                    value = lowerValue - value - primalTolerance;
     653  assert(method_ != 3);
     654#endif
     655  if (CLP_METHOD1) {
     656    // nonbasic should be at a valid bound
     657    for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     658      double lowerValue;
     659      double upperValue;
     660      double value = solution[iSequence];
     661      int iRange;
     662      // get correct place
     663      int start = start_[iSequence];
     664      int end = start_[iSequence + 1] - 1;
     665      // correct costs for this infeasibility weight
     666      // If free then true cost will be first
     667      double thisFeasibleCost = cost_[start];
     668      if (infeasible(start)) {
     669        thisFeasibleCost = cost_[start + 1];
     670        cost_[start] = thisFeasibleCost - infeasibilityCost;
     671      }
     672      if (infeasible(end - 1)) {
     673        thisFeasibleCost = cost_[end - 2];
     674        cost_[end - 1] = thisFeasibleCost + infeasibilityCost;
     675      }
     676      for (iRange = start; iRange < end; iRange++) {
     677        if (value < lower_[iRange + 1] + primalTolerance) {
     678          // put in better range if infeasible
     679          if (value >= lower_[iRange + 1] - primalTolerance && infeasible(iRange) && iRange == start)
     680            iRange++;
     681          whichRange_[iSequence] = iRange;
     682          break;
     683        }
     684      }
     685      assert(iRange < end);
     686      lowerValue = lower_[iRange];
     687      upperValue = lower_[iRange + 1];
     688      ClpSimplex::Status status = model_->getStatus(iSequence);
     689      if (upperValue == lowerValue && status != ClpSimplex::isFixed) {
     690        if (status != ClpSimplex::basic) {
     691          model_->setStatus(iSequence, ClpSimplex::isFixed);
     692          status = ClpSimplex::isFixed;
     693        }
     694      }
     695      //#define PRINT_DETAIL7 2
     696#if PRINT_DETAIL7 > 1
     697      printf("NL %d sol %g bounds %g %g\n",
     698        iSequence, solution[iSequence],
     699        lowerValue, upperValue);
     700#endif
     701      switch (status) {
     702
     703      case ClpSimplex::basic:
     704      case ClpSimplex::superBasic:
     705        // iRange is in correct place
     706        // slot in here
     707        if (infeasible(iRange)) {
     708          if (lower_[iRange] < -1.0e50) {
     709            //cost_[iRange] = cost_[iRange+1]-infeasibilityCost;
     710            // possibly below
     711            lowerValue = lower_[iRange + 1];
     712            if (value - lowerValue < -primalTolerance) {
     713              value = lowerValue - value - primalTolerance;
    721714#ifndef NDEBUG
    722                                    if(value > 1.0e15)
    723                                         printf("nonlincostb %d %g %g %g\n",
    724                                                iSequence, lowerValue, solution[iSequence], lower_[iRange+2]);
     715              if (value > 1.0e15)
     716                printf("nonlincostb %d %g %g %g\n",
     717                  iSequence, lowerValue, solution[iSequence], lower_[iRange + 2]);
    725718#endif
    726719#if PRINT_DETAIL7
    727                                    printf("**NL %d sol %g below %g\n",
    728                                           iSequence,solution[iSequence],
    729                                           lowerValue);
    730 #endif
    731                                    sumInfeasibilities_ += value;
    732                                    largestInfeasibility_ = CoinMax(largestInfeasibility_, value);
    733                                    changeCost_ -= lowerValue *
    734                                                   (cost_[iRange] - cost[iSequence]);
    735                                    numberInfeasibilities_++;
    736                               }
    737                          } else {
    738                               //cost_[iRange] = cost_[iRange-1]+infeasibilityCost;
    739                               // possibly above
    740                               upperValue = lower_[iRange];
    741                               if (value - upperValue > primalTolerance) {
    742                                    value = value - upperValue - primalTolerance;
     720              printf("**NL %d sol %g below %g\n",
     721                iSequence, solution[iSequence],
     722                lowerValue);
     723#endif
     724              sumInfeasibilities_ += value;
     725              largestInfeasibility_ = CoinMax(largestInfeasibility_, value);
     726              changeCost_ -= lowerValue * (cost_[iRange] - cost[iSequence]);
     727              numberInfeasibilities_++;
     728            }
     729          } else {
     730            //cost_[iRange] = cost_[iRange-1]+infeasibilityCost;
     731            // possibly above
     732            upperValue = lower_[iRange];
     733            if (value - upperValue > primalTolerance) {
     734              value = value - upperValue - primalTolerance;
    743735#ifndef NDEBUG
    744                                    if(value > 1.0e15)
    745                                         printf("nonlincostu %d %g %g %g\n",
    746                                                iSequence, lower_[iRange-1], solution[iSequence], upperValue);
     736              if (value > 1.0e15)
     737                printf("nonlincostu %d %g %g %g\n",
     738                  iSequence, lower_[iRange - 1], solution[iSequence], upperValue);
    747739#endif
    748740#if PRINT_DETAIL7
    749                                    printf("**NL %d sol %g above %g\n",
    750                                           iSequence,solution[iSequence],
    751                                           upperValue);
    752 #endif
    753                                    sumInfeasibilities_ += value;
    754                                    largestInfeasibility_ = CoinMax(largestInfeasibility_, value);
    755                                    changeCost_ -= upperValue *
    756                                                   (cost_[iRange] - cost[iSequence]);
    757                                    numberInfeasibilities_++;
    758                               }
    759                          }
    760                     }
    761                     //lower[iSequence] = lower_[iRange];
    762                     //upper[iSequence] = lower_[iRange+1];
    763                     //cost[iSequence] = cost_[iRange];
    764                     break;
    765                case ClpSimplex::isFree:
    766                     //if (toNearest)
    767                     //solution[iSequence] = 0.0;
    768                     break;
    769                case ClpSimplex::atUpperBound:
    770                     if (!toNearest) {
    771                          // With increasing tolerances - we may be at wrong place
    772                          if (fabs(value - upperValue) > oldTolerance * 1.0001) {
    773                               if (fabs(value - lowerValue) <= oldTolerance * 1.0001) {
    774                                    if  (fabs(value - lowerValue) > primalTolerance)
    775                                         solution[iSequence] = lowerValue;
    776                                    model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    777                               } else {
    778                                    model_->setStatus(iSequence, ClpSimplex::superBasic);
    779                               }
    780                          } else if  (fabs(value - upperValue) > primalTolerance) {
    781                               solution[iSequence] = upperValue;
    782                          }
    783                     } else {
    784                          // Set to nearest and make at upper bound
    785                          int kRange;
    786                          iRange = -1;
    787                          double nearest = COIN_DBL_MAX;
    788                          for (kRange = start; kRange < end; kRange++) {
    789                               if (fabs(lower_[kRange] - value) < nearest) {
    790                                    nearest = fabs(lower_[kRange] - value);
    791                                    iRange = kRange;
    792                               }
    793                          }
    794                          assert (iRange >= 0);
    795                          iRange--;
    796                          whichRange_[iSequence] = iRange;
    797                          solution[iSequence] = lower_[iRange+1];
    798                     }
    799                     break;
    800                case ClpSimplex::atLowerBound:
    801                     if (!toNearest) {
    802                          // With increasing tolerances - we may be at wrong place
    803                          // below stops compiler error with gcc 3.2!!!
    804                          if (iSequence == -119)
    805                               printf("ZZ %g %g %g %g\n", lowerValue, value, upperValue, oldTolerance);
    806                          if (fabs(value - lowerValue) > oldTolerance * 1.0001) {
    807                               if (fabs(value - upperValue) <= oldTolerance * 1.0001) {
    808                                    if  (fabs(value - upperValue) > primalTolerance)
    809                                         solution[iSequence] = upperValue;
    810                                    model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    811                               } else {
    812                                    model_->setStatus(iSequence, ClpSimplex::superBasic);
    813                               }
    814                          } else if  (fabs(value - lowerValue) > primalTolerance) {
    815                               solution[iSequence] = lowerValue;
    816                          }
    817                     } else {
    818                          // Set to nearest and make at lower bound
    819                          int kRange;
    820                          iRange = -1;
    821                          double nearest = COIN_DBL_MAX;
    822                          for (kRange = start; kRange < end; kRange++) {
    823                               if (fabs(lower_[kRange] - value) < nearest) {
    824                                    nearest = fabs(lower_[kRange] - value);
    825                                    iRange = kRange;
    826                               }
    827                          }
    828                          assert (iRange >= 0);
    829                          whichRange_[iSequence] = iRange;
    830                          solution[iSequence] = lower_[iRange];
    831                     }
    832                     break;
    833                case ClpSimplex::isFixed:
    834                     if (toNearest) {
    835                          // Set to true fixed
    836                          for (iRange = start; iRange < end; iRange++) {
    837                               if (lower_[iRange] == lower_[iRange+1])
    838                                    break;
    839                          }
    840                          if (iRange == end) {
    841                               // Odd - but make sensible
    842                               // Set to nearest and make at lower bound
    843                               int kRange;
    844                               iRange = -1;
    845                               double nearest = COIN_DBL_MAX;
    846                               for (kRange = start; kRange < end; kRange++) {
    847                                    if (fabs(lower_[kRange] - value) < nearest) {
    848                                         nearest = fabs(lower_[kRange] - value);
    849                                         iRange = kRange;
    850                                    }
    851                               }
    852                               assert (iRange >= 0);
    853                               whichRange_[iSequence] = iRange;
    854                               if (lower_[iRange] != lower_[iRange+1])
    855                                    model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    856                               else
    857                                    model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    858                          }
    859                          solution[iSequence] = lower_[iRange];
    860                     }
    861                     break;
    862                }
    863                lower[iSequence] = lower_[iRange];
    864                upper[iSequence] = lower_[iRange+1];
    865                cost[iSequence] = cost_[iRange];
    866                feasibleCost_ += thisFeasibleCost * solution[iSequence];
    867                //assert (iRange==whichRange_[iSequence]);
     741              printf("**NL %d sol %g above %g\n",
     742                iSequence, solution[iSequence],
     743                upperValue);
     744#endif
     745              sumInfeasibilities_ += value;
     746              largestInfeasibility_ = CoinMax(largestInfeasibility_, value);
     747              changeCost_ -= upperValue * (cost_[iRange] - cost[iSequence]);
     748              numberInfeasibilities_++;
     749            }
    868750          }
    869      }
     751        }
     752        //lower[iSequence] = lower_[iRange];
     753        //upper[iSequence] = lower_[iRange+1];
     754        //cost[iSequence] = cost_[iRange];
     755        break;
     756      case ClpSimplex::isFree:
     757        //if (toNearest)
     758        //solution[iSequence] = 0.0;
     759        break;
     760      case ClpSimplex::atUpperBound:
     761        if (!toNearest) {
     762          // With increasing tolerances - we may be at wrong place
     763          if (fabs(value - upperValue) > oldTolerance * 1.0001) {
     764            if (fabs(value - lowerValue) <= oldTolerance * 1.0001) {
     765              if (fabs(value - lowerValue) > primalTolerance)
     766                solution[iSequence] = lowerValue;
     767              model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     768            } else {
     769              model_->setStatus(iSequence, ClpSimplex::superBasic);
     770            }
     771          } else if (fabs(value - upperValue) > primalTolerance) {
     772            solution[iSequence] = upperValue;
     773          }
     774        } else {
     775          // Set to nearest and make at upper bound
     776          int kRange;
     777          iRange = -1;
     778          double nearest = COIN_DBL_MAX;
     779          for (kRange = start; kRange < end; kRange++) {
     780            if (fabs(lower_[kRange] - value) < nearest) {
     781              nearest = fabs(lower_[kRange] - value);
     782              iRange = kRange;
     783            }
     784          }
     785          assert(iRange >= 0);
     786          iRange--;
     787          whichRange_[iSequence] = iRange;
     788          solution[iSequence] = lower_[iRange + 1];
     789        }
     790        break;
     791      case ClpSimplex::atLowerBound:
     792        if (!toNearest) {
     793          // With increasing tolerances - we may be at wrong place
     794          // below stops compiler error with gcc 3.2!!!
     795          if (iSequence == -119)
     796            printf("ZZ %g %g %g %g\n", lowerValue, value, upperValue, oldTolerance);
     797          if (fabs(value - lowerValue) > oldTolerance * 1.0001) {
     798            if (fabs(value - upperValue) <= oldTolerance * 1.0001) {
     799              if (fabs(value - upperValue) > primalTolerance)
     800                solution[iSequence] = upperValue;
     801              model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     802            } else {
     803              model_->setStatus(iSequence, ClpSimplex::superBasic);
     804            }
     805          } else if (fabs(value - lowerValue) > primalTolerance) {
     806            solution[iSequence] = lowerValue;
     807          }
     808        } else {
     809          // Set to nearest and make at lower bound
     810          int kRange;
     811          iRange = -1;
     812          double nearest = COIN_DBL_MAX;
     813          for (kRange = start; kRange < end; kRange++) {
     814            if (fabs(lower_[kRange] - value) < nearest) {
     815              nearest = fabs(lower_[kRange] - value);
     816              iRange = kRange;
     817            }
     818          }
     819          assert(iRange >= 0);
     820          whichRange_[iSequence] = iRange;
     821          solution[iSequence] = lower_[iRange];
     822        }
     823        break;
     824      case ClpSimplex::isFixed:
     825        if (toNearest) {
     826          // Set to true fixed
     827          for (iRange = start; iRange < end; iRange++) {
     828            if (lower_[iRange] == lower_[iRange + 1])
     829              break;
     830          }
     831          if (iRange == end) {
     832            // Odd - but make sensible
     833            // Set to nearest and make at lower bound
     834            int kRange;
     835            iRange = -1;
     836            double nearest = COIN_DBL_MAX;
     837            for (kRange = start; kRange < end; kRange++) {
     838              if (fabs(lower_[kRange] - value) < nearest) {
     839                nearest = fabs(lower_[kRange] - value);
     840                iRange = kRange;
     841              }
     842            }
     843            assert(iRange >= 0);
     844            whichRange_[iSequence] = iRange;
     845            if (lower_[iRange] != lower_[iRange + 1])
     846              model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     847            else
     848              model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     849          }
     850          solution[iSequence] = lower_[iRange];
     851        }
     852        break;
     853      }
     854      lower[iSequence] = lower_[iRange];
     855      upper[iSequence] = lower_[iRange + 1];
     856      cost[iSequence] = cost_[iRange];
     857      feasibleCost_ += thisFeasibleCost * solution[iSequence];
     858      //assert (iRange==whichRange_[iSequence]);
     859    }
     860  }
    870861#ifdef NONLIN_DEBUG
    871      double saveCost = feasibleCost_;
    872      if (method_ == 3) {
    873           feasibleCost_ = 0.0;
    874           // Put back solution as we will be checking
    875           unsigned char * statusA = model_->statusArray();
    876           for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    877                double value = solution[iSequence];
    878                solution[iSequence] = saveSolution[iSequence];
    879                saveSolution[iSequence] = value;
    880                value = lower[iSequence];
    881                lower[iSequence] = saveLower[iSequence];
    882                saveLower[iSequence] = value;
    883                value = upper[iSequence];
    884                upper[iSequence] = saveUpper[iSequence];
    885                saveUpper[iSequence] = value;
    886                unsigned char value2 = statusA[iSequence];
    887                statusA[iSequence] = saveStatus[iSequence];
    888                saveStatus[iSequence] = value2;
     862  double saveCost = feasibleCost_;
     863  if (method_ == 3) {
     864    feasibleCost_ = 0.0;
     865    // Put back solution as we will be checking
     866    unsigned char *statusA = model_->statusArray();
     867    for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     868      double value = solution[iSequence];
     869      solution[iSequence] = saveSolution[iSequence];
     870      saveSolution[iSequence] = value;
     871      value = lower[iSequence];
     872      lower[iSequence] = saveLower[iSequence];
     873      saveLower[iSequence] = value;
     874      value = upper[iSequence];
     875      upper[iSequence] = saveUpper[iSequence];
     876      saveUpper[iSequence] = value;
     877      unsigned char value2 = statusA[iSequence];
     878      statusA[iSequence] = saveStatus[iSequence];
     879      saveStatus[iSequence] = value2;
     880    }
     881  }
     882#endif
     883  if (CLP_METHOD2) {
     884    //#define CLP_NON_JUST_BASIC
     885#ifndef CLP_NON_JUST_BASIC
     886    // nonbasic should be at a valid bound
     887    for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     888#else
     889    const int *pivotVariable = model_->pivotVariable();
     890    for (int i = 0; i < numberRows_; i++) {
     891      int iSequence = pivotVariable[i];
     892#endif
     893      double value = solution[iSequence];
     894      unsigned char iStatus = status_[iSequence];
     895      assert(currentStatus(iStatus) == CLP_SAME);
     896      double lowerValue = lower[iSequence];
     897      double upperValue = upper[iSequence];
     898      double costValue = cost2_[iSequence];
     899      double trueCost = costValue;
     900      int iWhere = originalStatus(iStatus);
     901      if (iWhere == CLP_BELOW_LOWER) {
     902        lowerValue = upperValue;
     903        upperValue = bound_[iSequence];
     904        costValue -= infeasibilityCost;
     905      } else if (iWhere == CLP_ABOVE_UPPER) {
     906        upperValue = lowerValue;
     907        lowerValue = bound_[iSequence];
     908        costValue += infeasibilityCost;
     909      }
     910      // get correct place
     911      int newWhere = CLP_FEASIBLE;
     912      ClpSimplex::Status status = model_->getStatus(iSequence);
     913      if (upperValue == lowerValue && status != ClpSimplex::isFixed) {
     914        if (status != ClpSimplex::basic) {
     915          model_->setStatus(iSequence, ClpSimplex::isFixed);
     916          status = ClpSimplex::isFixed;
     917        }
     918      }
     919      switch (status) {
     920
     921      case ClpSimplex::basic:
     922      case ClpSimplex::superBasic:
     923        if (value - upperValue <= primalTolerance) {
     924          if (value - lowerValue >= -primalTolerance) {
     925            // feasible
     926            //newWhere=CLP_FEASIBLE;
     927          } else {
     928            // below
     929            newWhere = CLP_BELOW_LOWER;
     930            assert(fabs(lowerValue) < 1.0e100);
     931            double infeasibility = lowerValue - value - primalTolerance;
     932            sumInfeasibilities_ += infeasibility;
     933            largestInfeasibility_ = CoinMax(largestInfeasibility_, infeasibility);
     934            costValue = trueCost - infeasibilityCost;
     935            changeCost_ -= lowerValue * (costValue - cost[iSequence]);
     936            numberInfeasibilities_++;
    889937          }
    890      }
    891 #endif
    892      if (CLP_METHOD2) {
    893        //#define CLP_NON_JUST_BASIC
    894 #ifndef CLP_NON_JUST_BASIC
    895           // nonbasic should be at a valid bound
    896           for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    897 #else
    898           const int * pivotVariable = model_->pivotVariable();
    899           for (int i=0;i<numberRows_;i++) {
    900             int iSequence = pivotVariable[i];
    901 #endif
    902                double value = solution[iSequence];
    903                unsigned char iStatus = status_[iSequence];
    904                assert (currentStatus(iStatus) == CLP_SAME);
    905                double lowerValue = lower[iSequence];
    906                double upperValue = upper[iSequence];
    907                double costValue = cost2_[iSequence];
    908                double trueCost = costValue;
    909                int iWhere = originalStatus(iStatus);
    910                if (iWhere == CLP_BELOW_LOWER) {
    911                     lowerValue = upperValue;
    912                     upperValue = bound_[iSequence];
    913                     costValue -= infeasibilityCost;
    914                } else if (iWhere == CLP_ABOVE_UPPER) {
    915                     upperValue = lowerValue;
    916                     lowerValue = bound_[iSequence];
    917                     costValue += infeasibilityCost;
    918                }
    919                // get correct place
    920                int newWhere = CLP_FEASIBLE;
    921                ClpSimplex::Status status = model_->getStatus(iSequence);
    922                if (upperValue == lowerValue && status != ClpSimplex::isFixed) {
    923                     if (status != ClpSimplex::basic) {
    924                          model_->setStatus(iSequence, ClpSimplex::isFixed);
    925                          status = ClpSimplex::isFixed;
    926                     }
    927                }
    928                switch(status) {
    929 
    930                case ClpSimplex::basic:
    931                case ClpSimplex::superBasic:
    932                     if (value - upperValue <= primalTolerance) {
    933                          if (value - lowerValue >= -primalTolerance) {
    934                               // feasible
    935                               //newWhere=CLP_FEASIBLE;
    936                          } else {
    937                               // below
    938                               newWhere = CLP_BELOW_LOWER;
    939                               assert (fabs(lowerValue) < 1.0e100);
    940                               double infeasibility = lowerValue - value - primalTolerance;
    941                               sumInfeasibilities_ += infeasibility;
    942                               largestInfeasibility_ = CoinMax(largestInfeasibility_, infeasibility);
    943                               costValue = trueCost - infeasibilityCost;
    944                               changeCost_ -= lowerValue * (costValue - cost[iSequence]);
    945                               numberInfeasibilities_++;
    946                          }
    947                     } else {
    948                          // above
    949                          newWhere = CLP_ABOVE_UPPER;
    950                          double infeasibility = value - upperValue - primalTolerance;
    951                          sumInfeasibilities_ += infeasibility;
    952                          largestInfeasibility_ = CoinMax(largestInfeasibility_, infeasibility);
    953                          costValue = trueCost + infeasibilityCost;
    954                          changeCost_ -= upperValue * (costValue - cost[iSequence]);
    955                          numberInfeasibilities_++;
    956                     }
    957                     break;
    958                case ClpSimplex::isFree:
    959                     break;
    960                case ClpSimplex::atUpperBound:
    961                     if (!toNearest) {
    962                          // With increasing tolerances - we may be at wrong place
    963                          if (fabs(value - upperValue) > oldTolerance * 1.0001) {
    964                               if (fabs(value - lowerValue) <= oldTolerance * 1.0001) {
    965                                    if  (fabs(value - lowerValue) > primalTolerance) {
    966                                         solution[iSequence] = lowerValue;
    967                                         value = lowerValue;
    968                                    }
    969                                    model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    970                               } else {
    971                                    if (value < upperValue) {
    972                                         if (value > lowerValue) {
    973                                              model_->setStatus(iSequence, ClpSimplex::superBasic);
    974                                         } else {
    975                                              // set to lower bound as infeasible
    976                                              solution[iSequence] = lowerValue;
    977                                              value = lowerValue;
    978                                              model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    979                                         }
    980                                    } else {
    981                                         // set to upper bound as infeasible
    982                                         solution[iSequence] = upperValue;
    983                                         value = upperValue;
    984                                    }
    985                               }
    986                          } else if  (fabs(value - upperValue) > primalTolerance) {
    987                               solution[iSequence] = upperValue;
    988                               value = upperValue;
    989                          }
    990                     } else {
    991                          // Set to nearest and make at bound
    992                          if (fabs(value - lowerValue) < fabs(value - upperValue)) {
    993                               solution[iSequence] = lowerValue;
    994                               value = lowerValue;
    995                               model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    996                          } else {
    997                               solution[iSequence] = upperValue;
    998                               value = upperValue;
    999                          }
    1000                     }
    1001                     break;
    1002                case ClpSimplex::atLowerBound:
    1003                     if (!toNearest) {
    1004                          // With increasing tolerances - we may be at wrong place
    1005                          if (fabs(value - lowerValue) > oldTolerance * 1.0001) {
    1006                               if (fabs(value - upperValue) <= oldTolerance * 1.0001) {
    1007                                    if  (fabs(value - upperValue) > primalTolerance) {
    1008                                         solution[iSequence] = upperValue;
    1009                                         value = upperValue;
    1010                                    }
    1011                                    model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    1012                               } else {
    1013                                    if (value < upperValue) {
    1014                                         if (value > lowerValue) {
    1015                                              model_->setStatus(iSequence, ClpSimplex::superBasic);
    1016                                         } else {
    1017                                              // set to lower bound as infeasible
    1018                                              solution[iSequence] = lowerValue;
    1019                                              value = lowerValue;
    1020                                         }
    1021                                    } else {
    1022                                         // set to upper bound as infeasible
    1023                                         solution[iSequence] = upperValue;
    1024                                         value = upperValue;
    1025                                         model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    1026                                    }
    1027                               }
    1028                          } else if  (fabs(value - lowerValue) > primalTolerance) {
    1029                               solution[iSequence] = lowerValue;
    1030                               value = lowerValue;
    1031                          }
    1032                     } else {
    1033                          // Set to nearest and make at bound
    1034                          if (fabs(value - lowerValue) < fabs(value - upperValue)) {
    1035                               solution[iSequence] = lowerValue;
    1036                               value = lowerValue;
    1037                          } else {
    1038                               solution[iSequence] = upperValue;
    1039                               value = upperValue;
    1040                               model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    1041                          }
    1042                     }
    1043                     break;
    1044                case ClpSimplex::isFixed:
    1045                     solution[iSequence] = lowerValue;
    1046                     value = lowerValue;
    1047                     break;
    1048                }
     938        } else {
     939          // above
     940          newWhere = CLP_ABOVE_UPPER;
     941          double infeasibility = value - upperValue - primalTolerance;
     942          sumInfeasibilities_ += infeasibility;
     943          largestInfeasibility_ = CoinMax(largestInfeasibility_, infeasibility);
     944          costValue = trueCost + infeasibilityCost;
     945          changeCost_ -= upperValue * (costValue - cost[iSequence]);
     946          numberInfeasibilities_++;
     947        }
     948        break;
     949      case ClpSimplex::isFree:
     950        break;
     951      case ClpSimplex::atUpperBound:
     952        if (!toNearest) {
     953          // With increasing tolerances - we may be at wrong place
     954          if (fabs(value - upperValue) > oldTolerance * 1.0001) {
     955            if (fabs(value - lowerValue) <= oldTolerance * 1.0001) {
     956              if (fabs(value - lowerValue) > primalTolerance) {
     957                solution[iSequence] = lowerValue;
     958                value = lowerValue;
     959              }
     960              model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     961            } else {
     962              if (value < upperValue) {
     963                if (value > lowerValue) {
     964                  model_->setStatus(iSequence, ClpSimplex::superBasic);
     965                } else {
     966                  // set to lower bound as infeasible
     967                  solution[iSequence] = lowerValue;
     968                  value = lowerValue;
     969                  model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     970                }
     971              } else {
     972                // set to upper bound as infeasible
     973                solution[iSequence] = upperValue;
     974                value = upperValue;
     975              }
     976            }
     977          } else if (fabs(value - upperValue) > primalTolerance) {
     978            solution[iSequence] = upperValue;
     979            value = upperValue;
     980          }
     981        } else {
     982          // Set to nearest and make at bound
     983          if (fabs(value - lowerValue) < fabs(value - upperValue)) {
     984            solution[iSequence] = lowerValue;
     985            value = lowerValue;
     986            model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     987          } else {
     988            solution[iSequence] = upperValue;
     989            value = upperValue;
     990          }
     991        }
     992        break;
     993      case ClpSimplex::atLowerBound:
     994        if (!toNearest) {
     995          // With increasing tolerances - we may be at wrong place
     996          if (fabs(value - lowerValue) > oldTolerance * 1.0001) {
     997            if (fabs(value - upperValue) <= oldTolerance * 1.0001) {
     998              if (fabs(value - upperValue) > primalTolerance) {
     999                solution[iSequence] = upperValue;
     1000                value = upperValue;
     1001              }
     1002              model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     1003            } else {
     1004              if (value < upperValue) {
     1005                if (value > lowerValue) {
     1006                  model_->setStatus(iSequence, ClpSimplex::superBasic);
     1007                } else {
     1008                  // set to lower bound as infeasible
     1009                  solution[iSequence] = lowerValue;
     1010                  value = lowerValue;
     1011                }
     1012              } else {
     1013                // set to upper bound as infeasible
     1014                solution[iSequence] = upperValue;
     1015                value = upperValue;
     1016                model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     1017              }
     1018            }
     1019          } else if (fabs(value - lowerValue) > primalTolerance) {
     1020            solution[iSequence] = lowerValue;
     1021            value = lowerValue;
     1022          }
     1023        } else {
     1024          // Set to nearest and make at bound
     1025          if (fabs(value - lowerValue) < fabs(value - upperValue)) {
     1026            solution[iSequence] = lowerValue;
     1027            value = lowerValue;
     1028          } else {
     1029            solution[iSequence] = upperValue;
     1030            value = upperValue;
     1031            model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     1032          }
     1033        }
     1034        break;
     1035      case ClpSimplex::isFixed:
     1036        solution[iSequence] = lowerValue;
     1037        value = lowerValue;
     1038        break;
     1039      }
    10491040#ifdef NONLIN_DEBUG
    1050                double lo = saveLower[iSequence];
    1051                double up = saveUpper[iSequence];
    1052                double cc = cost[iSequence];
    1053                unsigned char ss = saveStatus[iSequence];
    1054                unsigned char snow = model_->statusArray()[iSequence];
    1055 #endif
    1056                if (iWhere != newWhere) {
    1057                     setOriginalStatus(status_[iSequence], newWhere);
    1058                     if (newWhere == CLP_BELOW_LOWER) {
    1059                          bound_[iSequence] = upperValue;
    1060                          upperValue = lowerValue;
    1061                          lowerValue = -COIN_DBL_MAX;
    1062                          costValue = trueCost - infeasibilityCost;
    1063                     } else if (newWhere == CLP_ABOVE_UPPER) {
    1064                          bound_[iSequence] = lowerValue;
    1065                          lowerValue = upperValue;
    1066                          upperValue = COIN_DBL_MAX;
    1067                          costValue = trueCost + infeasibilityCost;
    1068                     } else {
    1069                          costValue = trueCost;
    1070                     }
    1071                     lower[iSequence] = lowerValue;
    1072                     upper[iSequence] = upperValue;
    1073                     assert (lowerValue<=upperValue);
    1074                }
    1075                // always do as other things may change
    1076                cost[iSequence] = costValue;
     1041      double lo = saveLower[iSequence];
     1042      double up = saveUpper[iSequence];
     1043      double cc = cost[iSequence];
     1044      unsigned char ss = saveStatus[iSequence];
     1045      unsigned char snow = model_->statusArray()[iSequence];
     1046#endif
     1047      if (iWhere != newWhere) {
     1048        setOriginalStatus(status_[iSequence], newWhere);
     1049        if (newWhere == CLP_BELOW_LOWER) {
     1050          bound_[iSequence] = upperValue;
     1051          upperValue = lowerValue;
     1052          lowerValue = -COIN_DBL_MAX;
     1053          costValue = trueCost - infeasibilityCost;
     1054        } else if (newWhere == CLP_ABOVE_UPPER) {
     1055          bound_[iSequence] = lowerValue;
     1056          lowerValue = upperValue;
     1057          upperValue = COIN_DBL_MAX;
     1058          costValue = trueCost + infeasibilityCost;
     1059        } else {
     1060          costValue = trueCost;
     1061        }
     1062        lower[iSequence] = lowerValue;
     1063        upper[iSequence] = upperValue;
     1064        assert(lowerValue <= upperValue);
     1065      }
     1066      // always do as other things may change
     1067      cost[iSequence] = costValue;
    10771068#ifdef NONLIN_DEBUG
    1078                if (method_ == 3) {
    1079                     assert (ss == snow);
    1080                     assert (cc == cost[iSequence]);
    1081                     assert (lo == lower[iSequence]);
    1082                     assert (up == upper[iSequence]);
    1083                     assert (value == saveSolution[iSequence]);
    1084                }
    1085 #endif
    1086                feasibleCost_ += trueCost * value;
    1087           }
    1088      }
     1069      if (method_ == 3) {
     1070        assert(ss == snow);
     1071        assert(cc == cost[iSequence]);
     1072        assert(lo == lower[iSequence]);
     1073        assert(up == upper[iSequence]);
     1074        assert(value == saveSolution[iSequence]);
     1075      }
     1076#endif
     1077      feasibleCost_ += trueCost * value;
     1078    }
     1079  }
    10891080#ifdef NONLIN_DEBUG
    1090      if (method_ == 3)
    1091           assert (fabs(saveCost - feasibleCost_) < 0.001 * (1.0 + CoinMax(fabs(saveCost), fabs(feasibleCost_))));
    1092      delete [] saveSolution;
    1093      delete [] saveLower;
    1094      delete [] saveUpper;
    1095      delete [] saveStatus;
    1096 #endif
    1097      //feasibleCost_ /= (model_->rhsScale()*model_->objScale());
     1081  if (method_ == 3)
     1082    assert(fabs(saveCost - feasibleCost_) < 0.001 * (1.0 + CoinMax(fabs(saveCost), fabs(feasibleCost_))));
     1083  delete[] saveSolution;
     1084  delete[] saveLower;
     1085  delete[] saveUpper;
     1086  delete[] saveStatus;
     1087#endif
     1088  //feasibleCost_ /= (model_->rhsScale()*model_->objScale());
    10981089}
    10991090// Puts feasible bounds into lower and upper
    1100 void
    1101 ClpNonLinearCost::feasibleBounds()
    1102 {
    1103      if (CLP_METHOD2) {
    1104           int iSequence;
    1105           double * upper = model_->upperRegion();
    1106           double * lower = model_->lowerRegion();
    1107           double * cost = model_->costRegion();
    1108           int numberTotal = numberColumns_ + numberRows_;
    1109           for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    1110                unsigned char iStatus = status_[iSequence];
    1111                assert (currentStatus(iStatus) == CLP_SAME);
    1112                double lowerValue = lower[iSequence];
    1113                double upperValue = upper[iSequence];
    1114                double costValue = cost2_[iSequence];
    1115                int iWhere = originalStatus(iStatus);
    1116                if (iWhere == CLP_BELOW_LOWER) {
    1117                     lowerValue = upperValue;
    1118                     upperValue = bound_[iSequence];
    1119                     assert (fabs(lowerValue) < 1.0e100);
    1120                } else if (iWhere == CLP_ABOVE_UPPER) {
    1121                     upperValue = lowerValue;
    1122                     lowerValue = bound_[iSequence];
    1123                }
    1124                setOriginalStatus(status_[iSequence], CLP_FEASIBLE);
    1125                lower[iSequence] = lowerValue;
    1126                upper[iSequence] = upperValue;
    1127                cost[iSequence] = costValue;
    1128           }
    1129      }
     1091void ClpNonLinearCost::feasibleBounds()
     1092{
     1093  if (CLP_METHOD2) {
     1094    int iSequence;
     1095    double *upper = model_->upperRegion();
     1096    double *lower = model_->lowerRegion();
     1097    double *cost = model_->costRegion();
     1098    int numberTotal = numberColumns_ + numberRows_;
     1099    for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     1100      unsigned char iStatus = status_[iSequence];
     1101      assert(currentStatus(iStatus) == CLP_SAME);
     1102      double lowerValue = lower[iSequence];
     1103      double upperValue = upper[iSequence];
     1104      double costValue = cost2_[iSequence];
     1105      int iWhere = originalStatus(iStatus);
     1106      if (iWhere == CLP_BELOW_LOWER) {
     1107        lowerValue = upperValue;
     1108        upperValue = bound_[iSequence];
     1109        assert(fabs(lowerValue) < 1.0e100);
     1110      } else if (iWhere == CLP_ABOVE_UPPER) {
     1111        upperValue = lowerValue;
     1112        lowerValue = bound_[iSequence];
     1113      }
     1114      setOriginalStatus(status_[iSequence], CLP_FEASIBLE);
     1115      lower[iSequence] = lowerValue;
     1116      upper[iSequence] = upperValue;
     1117      cost[iSequence] = costValue;
     1118    }
     1119  }
    11301120}
    11311121/* Goes through one bound for each variable.
     
    11331123   The indices are row indices and need converting to sequences
    11341124*/
    1135 void
    1136 ClpNonLinearCost::goThru(int numberInArray, double multiplier,
    1137                          const int * index, const double * array,
    1138                          double * rhs)
    1139 {
    1140      assert (model_ != NULL);
    1141      abort();
    1142      const int * pivotVariable = model_->pivotVariable();
    1143      if (CLP_METHOD1) {
    1144           for (int i = 0; i < numberInArray; i++) {
    1145                int iRow = index[i];
    1146                int iSequence = pivotVariable[iRow];
    1147                double alpha = multiplier * array[iRow];
    1148                // get where in bound sequence
    1149                int iRange = whichRange_[iSequence];
    1150                iRange += offset_[iSequence]; //add temporary bias
    1151                double value = model_->solution(iSequence);
    1152                if (alpha > 0.0) {
    1153                     // down one
    1154                     iRange--;
    1155                     assert(iRange >= start_[iSequence]);
    1156                     rhs[iRow] = value - lower_[iRange];
    1157                } else {
    1158                     // up one
    1159                     iRange++;
    1160                     assert(iRange < start_[iSequence+1] - 1);
    1161                     rhs[iRow] = lower_[iRange+1] - value;
    1162                }
    1163                offset_[iSequence] = iRange - whichRange_[iSequence];
    1164           }
    1165      }
     1125void ClpNonLinearCost::goThru(int numberInArray, double multiplier,
     1126  const int *index, const double *array,
     1127  double *rhs)
     1128{
     1129  assert(model_ != NULL);
     1130  abort();
     1131  const int *pivotVariable = model_->pivotVariable();
     1132  if (CLP_METHOD1) {
     1133    for (int i = 0; i < numberInArray; i++) {
     1134      int iRow = index[i];
     1135      int iSequence = pivotVariable[iRow];
     1136      double alpha = multiplier * array[iRow];
     1137      // get where in bound sequence
     1138      int iRange = whichRange_[iSequence];
     1139      iRange += offset_[iSequence]; //add temporary bias
     1140      double value = model_->solution(iSequence);
     1141      if (alpha > 0.0) {
     1142        // down one
     1143        iRange--;
     1144        assert(iRange >= start_[iSequence]);
     1145        rhs[iRow] = value - lower_[iRange];
     1146      } else {
     1147        // up one
     1148        iRange++;
     1149        assert(iRange < start_[iSequence + 1] - 1);
     1150        rhs[iRow] = lower_[iRange + 1] - value;
     1151      }
     1152      offset_[iSequence] = iRange - whichRange_[iSequence];
     1153    }
     1154  }
    11661155#ifdef NONLIN_DEBUG
    1167      double * saveRhs = NULL;
    1168      if (method_ == 3) {
    1169           int numberRows = model_->numberRows();
    1170           saveRhs = CoinCopyOfArray(rhs, numberRows);
    1171      }
    1172 #endif
    1173      if (CLP_METHOD2) {
    1174           const double * solution = model_->solutionRegion();
    1175           const double * upper = model_->upperRegion();
    1176           const double * lower = model_->lowerRegion();
    1177           for (int i = 0; i < numberInArray; i++) {
    1178                int iRow = index[i];
    1179                int iSequence = pivotVariable[iRow];
    1180                double alpha = multiplier * array[iRow];
    1181                double value = solution[iSequence];
    1182                unsigned char iStatus = status_[iSequence];
    1183                int iWhere = currentStatus(iStatus);
    1184                if (iWhere == CLP_SAME)
    1185                     iWhere = originalStatus(iStatus);
    1186                if (iWhere == CLP_FEASIBLE) {
    1187                     if (alpha > 0.0) {
    1188                          // going below
    1189                          iWhere = CLP_BELOW_LOWER;
    1190                          rhs[iRow] = value - lower[iSequence];
    1191                     } else {
    1192                          // going above
    1193                          iWhere = CLP_ABOVE_UPPER;
    1194                          rhs[iRow] = upper[iSequence] - value;
    1195                     }
    1196                } else if(iWhere == CLP_BELOW_LOWER) {
    1197                     assert (alpha < 0);
    1198                     // going feasible
    1199                     iWhere = CLP_FEASIBLE;
    1200                     rhs[iRow] = upper[iSequence] - value;
    1201                } else {
    1202                     assert (iWhere == CLP_ABOVE_UPPER);
    1203                     // going feasible
    1204                     iWhere = CLP_FEASIBLE;
    1205                     rhs[iRow] = value - lower[iSequence];
    1206                }
     1156  double *saveRhs = NULL;
     1157  if (method_ == 3) {
     1158    int numberRows = model_->numberRows();
     1159    saveRhs = CoinCopyOfArray(rhs, numberRows);
     1160  }
     1161#endif
     1162  if (CLP_METHOD2) {
     1163    const double *solution = model_->solutionRegion();
     1164    const double *upper = model_->upperRegion();
     1165    const double *lower = model_->lowerRegion();
     1166    for (int i = 0; i < numberInArray; i++) {
     1167      int iRow = index[i];
     1168      int iSequence = pivotVariable[iRow];
     1169      double alpha = multiplier * array[iRow];
     1170      double value = solution[iSequence];
     1171      unsigned char iStatus = status_[iSequence];
     1172      int iWhere = currentStatus(iStatus);
     1173      if (iWhere == CLP_SAME)
     1174        iWhere = originalStatus(iStatus);
     1175      if (iWhere == CLP_FEASIBLE) {
     1176        if (alpha > 0.0) {
     1177          // going below
     1178          iWhere = CLP_BELOW_LOWER;
     1179          rhs[iRow] = value - lower[iSequence];
     1180        } else {
     1181          // going above
     1182          iWhere = CLP_ABOVE_UPPER;
     1183          rhs[iRow] = upper[iSequence] - value;
     1184        }
     1185      } else if (iWhere == CLP_BELOW_LOWER) {
     1186        assert(alpha < 0);
     1187        // going feasible
     1188        iWhere = CLP_FEASIBLE;
     1189        rhs[iRow] = upper[iSequence] - value;
     1190      } else {
     1191        assert(iWhere == CLP_ABOVE_UPPER);
     1192        // going feasible
     1193        iWhere = CLP_FEASIBLE;
     1194        rhs[iRow] = value - lower[iSequence];
     1195      }
    12071196#ifdef NONLIN_DEBUG
    1208                if (method_ == 3)
    1209                     assert (rhs[iRow] == saveRhs[iRow]);
    1210 #endif
    1211                setCurrentStatus(status_[iSequence], iWhere);
    1212           }
    1213      }
     1197      if (method_ == 3)
     1198        assert(rhs[iRow] == saveRhs[iRow]);
     1199#endif
     1200      setCurrentStatus(status_[iSequence], iWhere);
     1201    }
     1202  }
    12141203#ifdef NONLIN_DEBUG
    1215      delete [] saveRhs;
     1204  delete[] saveRhs;
    12161205#endif
    12171206}
    12181207/* Takes off last iteration (i.e. offsets closer to 0)
    12191208 */
    1220 void
    1221 ClpNonLinearCost::goBack(int numberInArray, const int * index,
    1222                          double * rhs)
    1223 {
    1224      assert (model_ != NULL);
    1225      abort();
    1226      const int * pivotVariable = model_->pivotVariable();
    1227      if (CLP_METHOD1) {
    1228           for (int i = 0; i < numberInArray; i++) {
    1229                int iRow = index[i];
    1230                int iSequence = pivotVariable[iRow];
    1231                // get where in bound sequence
    1232                int iRange = whichRange_[iSequence];
    1233                // get closer to original
    1234                if (offset_[iSequence] > 0) {
    1235                     offset_[iSequence]--;
    1236                     assert (offset_[iSequence] >= 0);
    1237                     iRange += offset_[iSequence]; //add temporary bias
    1238                     double value = model_->solution(iSequence);
    1239                     // up one
    1240                     assert(iRange < start_[iSequence+1] - 1);
    1241                     rhs[iRow] = lower_[iRange+1] - value; // was earlier lower_[iRange]
    1242                } else {
    1243                     offset_[iSequence]++;
    1244                     assert (offset_[iSequence] <= 0);
    1245                     iRange += offset_[iSequence]; //add temporary bias
    1246                     double value = model_->solution(iSequence);
    1247                     // down one
    1248                     assert(iRange >= start_[iSequence]);
    1249                     rhs[iRow] = value - lower_[iRange]; // was earlier lower_[iRange+1]
    1250                }
    1251           }
    1252      }
     1209void ClpNonLinearCost::goBack(int numberInArray, const int *index,
     1210  double *rhs)
     1211{
     1212  assert(model_ != NULL);
     1213  abort();
     1214  const int *pivotVariable = model_->pivotVariable();
     1215  if (CLP_METHOD1) {
     1216    for (int i = 0; i < numberInArray; i++) {
     1217      int iRow = index[i];
     1218      int iSequence = pivotVariable[iRow];
     1219      // get where in bound sequence
     1220      int iRange = whichRange_[iSequence];
     1221      // get closer to original
     1222      if (offset_[iSequence] > 0) {
     1223        offset_[iSequence]--;
     1224        assert(offset_[iSequence] >= 0);
     1225        iRange += offset_[iSequence]; //add temporary bias
     1226        double value = model_->solution(iSequence);
     1227        // up one
     1228        assert(iRange < start_[iSequence + 1] - 1);
     1229        rhs[iRow] = lower_[iRange + 1] - value; // was earlier lower_[iRange]
     1230      } else {
     1231        offset_[iSequence]++;
     1232        assert(offset_[iSequence] <= 0);
     1233        iRange += offset_[iSequence]; //add temporary bias
     1234        double value = model_->solution(iSequence);
     1235        // down one
     1236        assert(iRange >= start_[iSequence]);
     1237        rhs[iRow] = value - lower_[iRange]; // was earlier lower_[iRange+1]
     1238      }
     1239    }
     1240  }
    12531241#ifdef NONLIN_DEBUG
    1254      double * saveRhs = NULL;
    1255      if (method_ == 3) {
    1256           int numberRows = model_->numberRows();
    1257           saveRhs = CoinCopyOfArray(rhs, numberRows);
    1258      }
    1259 #endif
    1260      if (CLP_METHOD2) {
    1261           const double * solution = model_->solutionRegion();
    1262           const double * upper = model_->upperRegion();
    1263           const double * lower = model_->lowerRegion();
    1264           for (int i = 0; i < numberInArray; i++) {
    1265                int iRow = index[i];
    1266                int iSequence = pivotVariable[iRow];
    1267                double value = solution[iSequence];
    1268                unsigned char iStatus = status_[iSequence];
    1269                int iWhere = currentStatus(iStatus);
    1270                int original = originalStatus(iStatus);
    1271                assert (iWhere != CLP_SAME);
    1272                if (iWhere == CLP_FEASIBLE) {
    1273                     if (original == CLP_BELOW_LOWER) {
    1274                          // going below
    1275                          iWhere = CLP_BELOW_LOWER;
    1276                          rhs[iRow] = lower[iSequence] - value;
    1277                     } else {
    1278                          // going above
    1279                          iWhere = CLP_ABOVE_UPPER;
    1280                          rhs[iRow] = value - upper[iSequence];
    1281                     }
    1282                } else if(iWhere == CLP_BELOW_LOWER) {
    1283                     // going feasible
    1284                     iWhere = CLP_FEASIBLE;
    1285                     rhs[iRow] = value - upper[iSequence];
    1286                } else {
    1287                     // going feasible
    1288                     iWhere = CLP_FEASIBLE;
    1289                     rhs[iRow] = lower[iSequence] - value;
    1290                }
     1242  double *saveRhs = NULL;
     1243  if (method_ == 3) {
     1244    int numberRows = model_->numberRows();
     1245    saveRhs = CoinCopyOfArray(rhs, numberRows);
     1246  }
     1247#endif
     1248  if (CLP_METHOD2) {
     1249    const double *solution = model_->solutionRegion();
     1250    const double *upper = model_->upperRegion();
     1251    const double *lower = model_->lowerRegion();
     1252    for (int i = 0; i < numberInArray; i++) {
     1253      int iRow = index[i];
     1254      int iSequence = pivotVariable[iRow];
     1255      double value = solution[iSequence];
     1256      unsigned char iStatus = status_[iSequence];
     1257      int iWhere = currentStatus(iStatus);
     1258      int original = originalStatus(iStatus);
     1259      assert(iWhere != CLP_SAME);
     1260      if (iWhere == CLP_FEASIBLE) {
     1261        if (original == CLP_BELOW_LOWER) {
     1262          // going below
     1263          iWhere = CLP_BELOW_LOWER;
     1264          rhs[iRow] = lower[iSequence] - value;
     1265        } else {
     1266          // going above
     1267          iWhere = CLP_ABOVE_UPPER;
     1268          rhs[iRow] = value - upper[iSequence];
     1269        }
     1270      } else if (iWhere == CLP_BELOW_LOWER) {
     1271        // going feasible
     1272        iWhere = CLP_FEASIBLE;
     1273        rhs[iRow] = value - upper[iSequence];
     1274      } else {
     1275        // going feasible
     1276        iWhere = CLP_FEASIBLE;
     1277        rhs[iRow] = lower[iSequence] - value;
     1278      }
    12911279#ifdef NONLIN_DEBUG
    1292                if (method_ == 3)
    1293                     assert (rhs[iRow] == saveRhs[iRow]);
    1294 #endif
    1295                setCurrentStatus(status_[iSequence], iWhere);
    1296           }
    1297      }
     1280      if (method_ == 3)
     1281        assert(rhs[iRow] == saveRhs[iRow]);
     1282#endif
     1283      setCurrentStatus(status_[iSequence], iWhere);
     1284    }
     1285  }
    12981286#ifdef NONLIN_DEBUG
    1299      delete [] saveRhs;
    1300 #endif
    1301 }
    1302 void
    1303 ClpNonLinearCost::goBackAll(const CoinIndexedVector * update)
    1304 {
    1305      assert (model_ != NULL);
    1306      const int * pivotVariable = model_->pivotVariable();
    1307      int number = update->getNumElements();
    1308      const int * index = update->getIndices();
    1309      if (CLP_METHOD1) {
    1310           for (int i = 0; i < number; i++) {
    1311                int iRow = index[i];
    1312                int iSequence = pivotVariable[iRow];
    1313                offset_[iSequence] = 0;
    1314           }
     1287  delete[] saveRhs;
     1288#endif
     1289}
     1290void ClpNonLinearCost::goBackAll(const CoinIndexedVector *update)
     1291{
     1292  assert(model_ != NULL);
     1293  const int *pivotVariable = model_->pivotVariable();
     1294  int number = update->getNumElements();
     1295  const int *index = update->getIndices();
     1296  if (CLP_METHOD1) {
     1297    for (int i = 0; i < number; i++) {
     1298      int iRow = index[i];
     1299      int iSequence = pivotVariable[iRow];
     1300      offset_[iSequence] = 0;
     1301    }
    13151302#ifdef CLP_DEBUG
    1316           for (i = 0; i < numberRows_ + numberColumns_; i++)
    1317                assert(!offset_[i]);
    1318 #endif
    1319      }
    1320      if (CLP_METHOD2) {
    1321           for (int i = 0; i < number; i++) {
    1322                int iRow = index[i];
    1323                int iSequence = pivotVariable[iRow];
    1324                setSameStatus(status_[iSequence]);
    1325           }
    1326      }
    1327 }
    1328 void
    1329 ClpNonLinearCost::checkInfeasibilities(int numberInArray, const int * index)
    1330 {
    1331      assert (model_ != NULL);
    1332      double primalTolerance = model_->currentPrimalTolerance();
    1333      const int * pivotVariable = model_->pivotVariable();
    1334      if (CLP_METHOD1) {
    1335           for (int i = 0; i < numberInArray; i++) {
    1336                int iRow = index[i];
    1337                int iSequence = pivotVariable[iRow];
    1338                // get where in bound sequence
    1339                int iRange;
    1340                int currentRange = whichRange_[iSequence];
    1341                double value = model_->solution(iSequence);
    1342                int start = start_[iSequence];
    1343                int end = start_[iSequence+1] - 1;
    1344                for (iRange = start; iRange < end; iRange++) {
    1345                     if (value < lower_[iRange+1] + primalTolerance) {
    1346                          // put in better range
    1347                          if (value >= lower_[iRange+1] - primalTolerance && infeasible(iRange) && iRange == start)
    1348                               iRange++;
    1349                          break;
    1350                     }
    1351                }
    1352                assert(iRange < end);
    1353                assert(model_->getStatus(iSequence) == ClpSimplex::basic);
    1354                double & lower = model_->lowerAddress(iSequence);
    1355                double & upper = model_->upperAddress(iSequence);
    1356                double & cost = model_->costAddress(iSequence);
    1357                whichRange_[iSequence] = iRange;
    1358                if (iRange != currentRange) {
    1359                     if (infeasible(iRange))
    1360                          numberInfeasibilities_++;
    1361                     if (infeasible(currentRange))
    1362                          numberInfeasibilities_--;
    1363                }
    1364                lower = lower_[iRange];
    1365                upper = lower_[iRange+1];
    1366                cost = cost_[iRange];
    1367           }
    1368      }
    1369      if (CLP_METHOD2) {
    1370           double * solution = model_->solutionRegion();
    1371           double * upper = model_->upperRegion();
    1372           double * lower = model_->lowerRegion();
    1373           double * cost = model_->costRegion();
    1374           for (int i = 0; i < numberInArray; i++) {
    1375                int iRow = index[i];
    1376                int iSequence = pivotVariable[iRow];
    1377                double value = solution[iSequence];
    1378                unsigned char iStatus = status_[iSequence];
    1379                assert (currentStatus(iStatus) == CLP_SAME);
    1380                double lowerValue = lower[iSequence];
    1381                double upperValue = upper[iSequence];
    1382                double costValue = cost2_[iSequence];
    1383                int iWhere = originalStatus(iStatus);
    1384                if (iWhere == CLP_BELOW_LOWER) {
    1385                     lowerValue = upperValue;
    1386                     upperValue = bound_[iSequence];
    1387                     numberInfeasibilities_--;
    1388                     assert (fabs(lowerValue) < 1.0e100);
    1389                } else if (iWhere == CLP_ABOVE_UPPER) {
    1390                     upperValue = lowerValue;
    1391                     lowerValue = bound_[iSequence];
    1392                     numberInfeasibilities_--;
    1393                }
    1394                // get correct place
    1395                int newWhere = CLP_FEASIBLE;
    1396                if (value - upperValue <= primalTolerance) {
    1397                     if (value - lowerValue >= -primalTolerance) {
    1398                          // feasible
    1399                          //newWhere=CLP_FEASIBLE;
    1400                     } else {
    1401                          // below
    1402                          newWhere = CLP_BELOW_LOWER;
    1403                          assert (fabs(lowerValue) < 1.0e100);
    1404                          costValue -= infeasibilityWeight_;
    1405                          numberInfeasibilities_++;
    1406                     }
    1407                } else {
    1408                     // above
    1409                     newWhere = CLP_ABOVE_UPPER;
    1410                     costValue += infeasibilityWeight_;
    1411                     numberInfeasibilities_++;
    1412                }
    1413                if (iWhere != newWhere) {
    1414                     setOriginalStatus(status_[iSequence], newWhere);
    1415                     if (newWhere == CLP_BELOW_LOWER) {
    1416                          bound_[iSequence] = upperValue;
    1417                          upperValue = lowerValue;
    1418                          lowerValue = -COIN_DBL_MAX;
    1419                     } else if (newWhere == CLP_ABOVE_UPPER) {
    1420                          bound_[iSequence] = lowerValue;
    1421                          lowerValue = upperValue;
    1422                          upperValue = COIN_DBL_MAX;
    1423                     }
    1424                     lower[iSequence] = lowerValue;
    1425                     upper[iSequence] = upperValue;
    1426                     cost[iSequence] = costValue;
    1427                }
    1428           }
    1429      }
     1303    for (i = 0; i < numberRows_ + numberColumns_; i++)
     1304      assert(!offset_[i]);
     1305#endif
     1306  }
     1307  if (CLP_METHOD2) {
     1308    for (int i = 0; i < number; i++) {
     1309      int iRow = index[i];
     1310      int iSequence = pivotVariable[iRow];
     1311      setSameStatus(status_[iSequence]);
     1312    }
     1313  }
     1314}
     1315void ClpNonLinearCost::checkInfeasibilities(int numberInArray, const int *index)
     1316{
     1317  assert(model_ != NULL);
     1318  double primalTolerance = model_->currentPrimalTolerance();
     1319  const int *pivotVariable = model_->pivotVariable();
     1320  if (CLP_METHOD1) {
     1321    for (int i = 0; i < numberInArray; i++) {
     1322      int iRow = index[i];
     1323      int iSequence = pivotVariable[iRow];
     1324      // get where in bound sequence
     1325      int iRange;
     1326      int currentRange = whichRange_[iSequence];
     1327      double value = model_->solution(iSequence);
     1328      int start = start_[iSequence];
     1329      int end = start_[iSequence + 1] - 1;
     1330      for (iRange = start; iRange < end; iRange++) {
     1331        if (value < lower_[iRange + 1] + primalTolerance) {
     1332          // put in better range
     1333          if (value >= lower_[iRange + 1] - primalTolerance && infeasible(iRange) && iRange == start)
     1334            iRange++;
     1335          break;
     1336        }
     1337      }
     1338      assert(iRange < end);
     1339      assert(model_->getStatus(iSequence) == ClpSimplex::basic);
     1340      double &lower = model_->lowerAddress(iSequence);
     1341      double &upper = model_->upperAddress(iSequence);
     1342      double &cost = model_->costAddress(iSequence);
     1343      whichRange_[iSequence] = iRange;
     1344      if (iRange != currentRange) {
     1345        if (infeasible(iRange))
     1346          numberInfeasibilities_++;
     1347        if (infeasible(currentRange))
     1348          numberInfeasibilities_--;
     1349      }
     1350      lower = lower_[iRange];
     1351      upper = lower_[iRange + 1];
     1352      cost = cost_[iRange];
     1353    }
     1354  }
     1355  if (CLP_METHOD2) {
     1356    double *solution = model_->solutionRegion();
     1357    double *upper = model_->upperRegion();
     1358    double *lower = model_->lowerRegion();
     1359    double *cost = model_->costRegion();
     1360    for (int i = 0; i < numberInArray; i++) {
     1361      int iRow = index[i];
     1362      int iSequence = pivotVariable[iRow];
     1363      double value = solution[iSequence];
     1364      unsigned char iStatus = status_[iSequence];
     1365      assert(currentStatus(iStatus) == CLP_SAME);
     1366      double lowerValue = lower[iSequence];
     1367      double upperValue = upper[iSequence];
     1368      double costValue = cost2_[iSequence];
     1369      int iWhere = originalStatus(iStatus);
     1370      if (iWhere == CLP_BELOW_LOWER) {
     1371        lowerValue = upperValue;
     1372        upperValue = bound_[iSequence];
     1373        numberInfeasibilities_--;
     1374        assert(fabs(lowerValue) < 1.0e100);
     1375      } else if (iWhere == CLP_ABOVE_UPPER) {
     1376        upperValue = lowerValue;
     1377        lowerValue = bound_[iSequence];
     1378        numberInfeasibilities_--;
     1379      }
     1380      // get correct place
     1381      int newWhere = CLP_FEASIBLE;
     1382      if (value - upperValue <= primalTolerance) {
     1383        if (value - lowerValue >= -primalTolerance) {
     1384          // feasible
     1385          //newWhere=CLP_FEASIBLE;
     1386        } else {
     1387          // below
     1388          newWhere = CLP_BELOW_LOWER;
     1389          assert(fabs(lowerValue) < 1.0e100);
     1390          costValue -= infeasibilityWeight_;
     1391          numberInfeasibilities_++;
     1392        }
     1393      } else {
     1394        // above
     1395        newWhere = CLP_ABOVE_UPPER;
     1396        costValue += infeasibilityWeight_;
     1397        numberInfeasibilities_++;
     1398      }
     1399      if (iWhere != newWhere) {
     1400        setOriginalStatus(status_[iSequence], newWhere);
     1401        if (newWhere == CLP_BELOW_LOWER) {
     1402          bound_[iSequence] = upperValue;
     1403          upperValue = lowerValue;
     1404          lowerValue = -COIN_DBL_MAX;
     1405        } else if (newWhere == CLP_ABOVE_UPPER) {
     1406          bound_[iSequence] = lowerValue;
     1407          lowerValue = upperValue;
     1408          upperValue = COIN_DBL_MAX;
     1409        }
     1410        lower[iSequence] = lowerValue;
     1411        upper[iSequence] = upperValue;
     1412        cost[iSequence] = costValue;
     1413      }
     1414    }
     1415  }
    14301416}
    14311417/* Puts back correct infeasible costs for each variable
     
    14351421   changed costs will be stored as normal CoinIndexedVector
    14361422*/
    1437 void
    1438 ClpNonLinearCost::checkChanged(int numberInArray, CoinIndexedVector * update)
    1439 {
    1440      assert (model_ != NULL);
    1441      double primalTolerance = model_->currentPrimalTolerance();
    1442      const int * pivotVariable = model_->pivotVariable();
    1443      int number = 0;
    1444      int * index = update->getIndices();
    1445      double * work = update->denseVector();
    1446      if (CLP_METHOD1) {
    1447           for (int i = 0; i < numberInArray; i++) {
    1448                int iRow = index[i];
    1449                int iSequence = pivotVariable[iRow];
    1450                // get where in bound sequence
    1451                int iRange;
    1452                double value = model_->solution(iSequence);
    1453                int start = start_[iSequence];
    1454                int end = start_[iSequence+1] - 1;
    1455                for (iRange = start; iRange < end; iRange++) {
    1456                     if (value < lower_[iRange+1] + primalTolerance) {
    1457                          // put in better range
    1458                          if (value >= lower_[iRange+1] - primalTolerance && infeasible(iRange) && iRange == start)
    1459                               iRange++;
    1460                          break;
    1461                     }
    1462                }
    1463                assert(iRange < end);
    1464                assert(model_->getStatus(iSequence) == ClpSimplex::basic);
    1465                int jRange = whichRange_[iSequence];
    1466                if (iRange != jRange) {
    1467                     // changed
    1468                     work[iRow] = cost_[jRange] - cost_[iRange];
    1469                     index[number++] = iRow;
    1470                     double & lower = model_->lowerAddress(iSequence);
    1471                     double & upper = model_->upperAddress(iSequence);
    1472                     double & cost = model_->costAddress(iSequence);
    1473                     whichRange_[iSequence] = iRange;
    1474                     if (infeasible(iRange))
    1475                          numberInfeasibilities_++;
    1476                     if (infeasible(jRange))
    1477                          numberInfeasibilities_--;
    1478                     lower = lower_[iRange];
    1479                     upper = lower_[iRange+1];
    1480                     cost = cost_[iRange];
    1481                }
    1482           }
    1483      }
    1484      if (CLP_METHOD2) {
    1485           double * solution = model_->solutionRegion();
    1486           double * upper = model_->upperRegion();
    1487           double * lower = model_->lowerRegion();
    1488           double * cost = model_->costRegion();
    1489           for (int i = 0; i < numberInArray; i++) {
    1490                int iRow = index[i];
    1491                int iSequence = pivotVariable[iRow];
    1492                double value = solution[iSequence];
    1493                unsigned char iStatus = status_[iSequence];
    1494                assert (currentStatus(iStatus) == CLP_SAME);
    1495                double lowerValue = lower[iSequence];
    1496                double upperValue = upper[iSequence];
    1497                double costValue = cost2_[iSequence];
    1498                int iWhere = originalStatus(iStatus);
    1499                if (iWhere == CLP_BELOW_LOWER) {
    1500                     lowerValue = upperValue;
    1501                     upperValue = bound_[iSequence];
    1502                     numberInfeasibilities_--;
    1503                     assert (fabs(lowerValue) < 1.0e100);
    1504                } else if (iWhere == CLP_ABOVE_UPPER) {
    1505                     upperValue = lowerValue;
    1506                     lowerValue = bound_[iSequence];
    1507                     numberInfeasibilities_--;
    1508                }
    1509                // get correct place
    1510                int newWhere = CLP_FEASIBLE;
    1511                if (value - upperValue <= primalTolerance) {
    1512                     if (value - lowerValue >= -primalTolerance) {
    1513                          // feasible
    1514                          //newWhere=CLP_FEASIBLE;
    1515                     } else {
    1516                          // below
    1517                          newWhere = CLP_BELOW_LOWER;
    1518                          costValue -= infeasibilityWeight_;
    1519                          numberInfeasibilities_++;
    1520                          assert (fabs(lowerValue) < 1.0e100);
    1521                     }
    1522                } else {
    1523                     // above
    1524                     newWhere = CLP_ABOVE_UPPER;
    1525                     costValue += infeasibilityWeight_;
    1526                     numberInfeasibilities_++;
    1527                }
    1528                if (iWhere != newWhere) {
    1529                     work[iRow] = cost[iSequence] - costValue;
    1530                     index[number++] = iRow;
    1531                     setOriginalStatus(status_[iSequence], newWhere);
    1532                     if (newWhere == CLP_BELOW_LOWER) {
    1533                          bound_[iSequence] = upperValue;
    1534                          upperValue = lowerValue;
    1535                          lowerValue = -COIN_DBL_MAX;
    1536                     } else if (newWhere == CLP_ABOVE_UPPER) {
    1537                          bound_[iSequence] = lowerValue;
    1538                          lowerValue = upperValue;
    1539                          upperValue = COIN_DBL_MAX;
    1540                     }
    1541                     lower[iSequence] = lowerValue;
    1542                     upper[iSequence] = upperValue;
    1543                     cost[iSequence] = costValue;
    1544                }
    1545           }
    1546      }
    1547      update->setNumElements(number);
     1423void ClpNonLinearCost::checkChanged(int numberInArray, CoinIndexedVector *update)
     1424{
     1425  assert(model_ != NULL);
     1426  double primalTolerance = model_->currentPrimalTolerance();
     1427  const int *pivotVariable = model_->pivotVariable();
     1428  int number = 0;
     1429  int *index = update->getIndices();
     1430  double *work = update->denseVector();
     1431  if (CLP_METHOD1) {
     1432    for (int i = 0; i < numberInArray; i++) {
     1433      int iRow = index[i];
     1434      int iSequence = pivotVariable[iRow];
     1435      // get where in bound sequence
     1436      int iRange;
     1437      double value = model_->solution(iSequence);
     1438      int start = start_[iSequence];
     1439      int end = start_[iSequence + 1] - 1;
     1440      for (iRange = start; iRange < end; iRange++) {
     1441        if (value < lower_[iRange + 1] + primalTolerance) {
     1442          // put in better range
     1443          if (value >= lower_[iRange + 1] - primalTolerance && infeasible(iRange) && iRange == start)
     1444            iRange++;
     1445          break;
     1446        }
     1447      }
     1448      assert(iRange < end);
     1449      assert(model_->getStatus(iSequence) == ClpSimplex::basic);
     1450      int jRange = whichRange_[iSequence];
     1451      if (iRange != jRange) {
     1452        // changed
     1453        work[iRow] = cost_[jRange] - cost_[iRange];
     1454        index[number++] = iRow;
     1455        double &lower = model_->lowerAddress(iSequence);
     1456        double &upper = model_->upperAddress(iSequence);
     1457        double &cost = model_->costAddress(iSequence);
     1458        whichRange_[iSequence] = iRange;
     1459        if (infeasible(iRange))
     1460          numberInfeasibilities_++;
     1461        if (infeasible(jRange))
     1462          numberInfeasibilities_--;
     1463        lower = lower_[iRange];
     1464        upper = lower_[iRange + 1];
     1465        cost = cost_[iRange];
     1466      }
     1467    }
     1468  }
     1469  if (CLP_METHOD2) {
     1470    double *solution = model_->solutionRegion();
     1471    double *upper = model_->upperRegion();
     1472    double *lower = model_->lowerRegion();
     1473    double *cost = model_->costRegion();
     1474    for (int i = 0; i < numberInArray; i++) {
     1475      int iRow = index[i];
     1476      int iSequence = pivotVariable[iRow];
     1477      double value = solution[iSequence];
     1478      unsigned char iStatus = status_[iSequence];
     1479      assert(currentStatus(iStatus) == CLP_SAME);
     1480      double lowerValue = lower[iSequence];
     1481      double upperValue = upper[iSequence];
     1482      double costValue = cost2_[iSequence];
     1483      int iWhere = originalStatus(iStatus);
     1484      if (iWhere == CLP_BELOW_LOWER) {
     1485        lowerValue = upperValue;
     1486        upperValue = bound_[iSequence];
     1487        numberInfeasibilities_--;
     1488        assert(fabs(lowerValue) < 1.0e100);
     1489      } else if (iWhere == CLP_ABOVE_UPPER) {
     1490        upperValue = lowerValue;
     1491        lowerValue = bound_[iSequence];
     1492        numberInfeasibilities_--;
     1493      }
     1494      // get correct place
     1495      int newWhere = CLP_FEASIBLE;
     1496      if (value - upperValue <= primalTolerance) {
     1497        if (value - lowerValue >= -primalTolerance) {
     1498          // feasible
     1499          //newWhere=CLP_FEASIBLE;
     1500        } else {
     1501          // below
     1502          newWhere = CLP_BELOW_LOWER;
     1503          costValue -= infeasibilityWeight_;
     1504          numberInfeasibilities_++;
     1505          assert(fabs(lowerValue) < 1.0e100);
     1506        }
     1507      } else {
     1508        // above
     1509        newWhere = CLP_ABOVE_UPPER;
     1510        costValue += infeasibilityWeight_;
     1511        numberInfeasibilities_++;
     1512      }
     1513      if (iWhere != newWhere) {
     1514        work[iRow] = cost[iSequence] - costValue;
     1515        index[number++] = iRow;
     1516        setOriginalStatus(status_[iSequence], newWhere);
     1517        if (newWhere == CLP_BELOW_LOWER) {
     1518          bound_[iSequence] = upperValue;
     1519          upperValue = lowerValue;
     1520          lowerValue = -COIN_DBL_MAX;
     1521        } else if (newWhere == CLP_ABOVE_UPPER) {
     1522          bound_[iSequence] = lowerValue;
     1523          lowerValue = upperValue;
     1524          upperValue = COIN_DBL_MAX;
     1525        }
     1526        lower[iSequence] = lowerValue;
     1527        upper[iSequence] = upperValue;
     1528        cost[iSequence] = costValue;
     1529      }
     1530    }
     1531  }
     1532  update->setNumElements(number);
    15481533}
    15491534/* Sets bounds and cost for one variable - returns change in cost*/
     
    15511536ClpNonLinearCost::setOne(int iSequence, double value)
    15521537{
    1553      assert (model_ != NULL);
    1554      double primalTolerance = model_->currentPrimalTolerance();
    1555      // difference in cost
    1556      double difference = 0.0;
    1557      if (CLP_METHOD1) {
    1558           // get where in bound sequence
    1559           int iRange;
    1560           int currentRange = whichRange_[iSequence];
    1561           int start = start_[iSequence];
    1562           int end = start_[iSequence+1] - 1;
    1563           if (!bothWays_) {
    1564                // If fixed try and get feasible
    1565                if (lower_[start+1] == lower_[start+2] && fabs(value - lower_[start+1]) < 1.001 * primalTolerance) {
    1566                     iRange = start + 1;
    1567                } else {
    1568                     for (iRange = start; iRange < end; iRange++) {
    1569                          if (value <= lower_[iRange+1] + primalTolerance) {
    1570                               // put in better range
    1571                               if (value >= lower_[iRange+1] - primalTolerance && infeasible(iRange) && iRange == start)
    1572                                    iRange++;
    1573                               break;
    1574                          }
    1575                     }
    1576                }
    1577           } else {
    1578                // leave in current if possible
    1579                iRange = whichRange_[iSequence];
    1580                if (value < lower_[iRange] - primalTolerance || value > lower_[iRange+1] + primalTolerance) {
    1581                     for (iRange = start; iRange < end; iRange++) {
    1582                          if (value < lower_[iRange+1] + primalTolerance) {
    1583                               // put in better range
    1584                               if (value >= lower_[iRange+1] - primalTolerance && infeasible(iRange) && iRange == start)
    1585                                    iRange++;
    1586                               break;
    1587                          }
    1588                     }
    1589                }
     1538  assert(model_ != NULL);
     1539  double primalTolerance = model_->currentPrimalTolerance();
     1540  // difference in cost
     1541  double difference = 0.0;
     1542  if (CLP_METHOD1) {
     1543    // get where in bound sequence
     1544    int iRange;
     1545    int currentRange = whichRange_[iSequence];
     1546    int start = start_[iSequence];
     1547    int end = start_[iSequence + 1] - 1;
     1548    if (!bothWays_) {
     1549      // If fixed try and get feasible
     1550      if (lower_[start + 1] == lower_[start + 2] && fabs(value - lower_[start + 1]) < 1.001 * primalTolerance) {
     1551        iRange = start + 1;
     1552      } else {
     1553        for (iRange = start; iRange < end; iRange++) {
     1554          if (value <= lower_[iRange + 1] + primalTolerance) {
     1555            // put in better range
     1556            if (value >= lower_[iRange + 1] - primalTolerance && infeasible(iRange) && iRange == start)
     1557              iRange++;
     1558            break;
    15901559          }
    1591           assert(iRange < end);
    1592           whichRange_[iSequence] = iRange;
    1593           if (iRange != currentRange) {
    1594                if (infeasible(iRange))
    1595                     numberInfeasibilities_++;
    1596                if (infeasible(currentRange))
    1597                     numberInfeasibilities_--;
     1560        }
     1561      }
     1562    } else {
     1563      // leave in current if possible
     1564      iRange = whichRange_[iSequence];
     1565      if (value < lower_[iRange] - primalTolerance || value > lower_[iRange + 1] + primalTolerance) {
     1566        for (iRange = start; iRange < end; iRange++) {
     1567          if (value < lower_[iRange + 1] + primalTolerance) {
     1568            // put in better range
     1569            if (value >= lower_[iRange + 1] - primalTolerance && infeasible(iRange) && iRange == start)
     1570              iRange++;
     1571            break;
    15981572          }
    1599           double & lower = model_->lowerAddress(iSequence);
    1600           double & upper = model_->upperAddress(iSequence);
    1601           double & cost = model_->costAddress(iSequence);
    1602           lower = lower_[iRange];
    1603           upper = lower_[iRange+1];
    1604           ClpSimplex::Status status = model_->getStatus(iSequence);
    1605           if (upper == lower) {
    1606                if (status != ClpSimplex::basic) {
    1607                     model_->setStatus(iSequence, ClpSimplex::isFixed);
    1608                     status = ClpSimplex::basic; // so will skip
    1609                }
    1610           }
    1611           switch(status) {
    1612 
    1613           case ClpSimplex::basic:
    1614           case ClpSimplex::superBasic:
    1615           case ClpSimplex::isFree:
    1616                break;
    1617           case ClpSimplex::atUpperBound:
    1618           case ClpSimplex::atLowerBound:
    1619           case ClpSimplex::isFixed:
    1620                // set correctly
    1621                if (fabs(value - lower) <= primalTolerance * 1.001) {
    1622                     model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    1623                } else if (fabs(value - upper) <= primalTolerance * 1.001) {
    1624                     model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    1625                } else {
    1626                     // set superBasic
    1627                     model_->setStatus(iSequence, ClpSimplex::superBasic);
    1628                }
    1629                break;
    1630           }
    1631           difference = cost - cost_[iRange];
    1632           cost = cost_[iRange];
    1633      }
    1634      if (CLP_METHOD2) {
    1635           double * upper = model_->upperRegion();
    1636           double * lower = model_->lowerRegion();
    1637           double * cost = model_->costRegion();
    1638           unsigned char iStatus = status_[iSequence];
    1639           assert (currentStatus(iStatus) == CLP_SAME);
    1640           double lowerValue = lower[iSequence];
    1641           double upperValue = upper[iSequence];
    1642           double costValue = cost2_[iSequence];
    1643           int iWhere = originalStatus(iStatus);
     1573        }
     1574      }
     1575    }
     1576    assert(iRange < end);
     1577    whichRange_[iSequence] = iRange;
     1578    if (iRange != currentRange) {
     1579      if (infeasible(iRange))
     1580        numberInfeasibilities_++;
     1581      if (infeasible(currentRange))
     1582        numberInfeasibilities_--;
     1583    }
     1584    double &lower = model_->lowerAddress(iSequence);
     1585    double &upper = model_->upperAddress(iSequence);
     1586    double &cost = model_->costAddress(iSequence);
     1587    lower = lower_[iRange];
     1588    upper = lower_[iRange + 1];
     1589    ClpSimplex::Status status = model_->getStatus(iSequence);
     1590    if (upper == lower) {
     1591      if (status != ClpSimplex::basic) {
     1592        model_->setStatus(iSequence, ClpSimplex::isFixed);
     1593        status = ClpSimplex::basic; // so will skip
     1594      }
     1595    }
     1596    switch (status) {
     1597
     1598    case ClpSimplex::basic:
     1599    case ClpSimplex::superBasic:
     1600    case ClpSimplex::isFree:
     1601      break;
     1602    case ClpSimplex::atUpperBound:
     1603    case ClpSimplex::atLowerBound:
     1604    case ClpSimplex::isFixed:
     1605      // set correctly
     1606      if (fabs(value - lower) <= primalTolerance * 1.001) {
     1607        model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     1608      } else if (fabs(value - upper) <= primalTolerance * 1.001) {
     1609        model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     1610      } else {
     1611        // set superBasic
     1612        model_->setStatus(iSequence, ClpSimplex::superBasic);
     1613      }
     1614      break;
     1615    }
     1616    difference = cost - cost_[iRange];
     1617    cost = cost_[iRange];
     1618  }
     1619  if (CLP_METHOD2) {
     1620    double *upper = model_->upperRegion();
     1621    double *lower = model_->lowerRegion();
     1622    double *cost = model_->costRegion();
     1623    unsigned char iStatus = status_[iSequence];
     1624    assert(currentStatus(iStatus) == CLP_SAME);
     1625    double lowerValue = lower[iSequence];
     1626    double upperValue = upper[iSequence];
     1627    double costValue = cost2_[iSequence];
     1628    int iWhere = originalStatus(iStatus);
    16441629#undef CLP_USER_DRIVEN
    16451630#ifdef CLP_USER_DRIVEN
    1646           clpUserStruct info;
    1647           info.type=3;
    1648           info.sequence=iSequence;
    1649           info.value=value;
    1650           info.change=0;
    1651           info.lower=lowerValue;
    1652           info.upper=upperValue;
    1653           info.cost=costValue;
    1654 #endif             
    1655           if (iWhere == CLP_BELOW_LOWER) {
    1656                lowerValue = upperValue;
    1657                upperValue = bound_[iSequence];
    1658                numberInfeasibilities_--;
    1659                assert (fabs(lowerValue) < 1.0e100);
    1660           } else if (iWhere == CLP_ABOVE_UPPER) {
    1661                upperValue = lowerValue;
    1662                lowerValue = bound_[iSequence];
    1663                numberInfeasibilities_--;
    1664           }
    1665           // get correct place
    1666           int newWhere = CLP_FEASIBLE;
    1667           if (value - upperValue <= primalTolerance) {
    1668                if (value - lowerValue >= -primalTolerance) {
    1669                     // feasible
    1670                     //newWhere=CLP_FEASIBLE;
    1671                } else {
    1672                     // below
    1673                     newWhere = CLP_BELOW_LOWER;
    1674                     costValue -= infeasibilityWeight_;
    1675                     numberInfeasibilities_++;
    1676                     assert (fabs(lowerValue) < 1.0e100);
    1677                }
    1678           } else {
    1679                // above
    1680                newWhere = CLP_ABOVE_UPPER;
    1681                costValue += infeasibilityWeight_;
    1682                numberInfeasibilities_++;
    1683           }
    1684           if (iWhere != newWhere) {
    1685                difference = cost[iSequence] - costValue;
    1686                setOriginalStatus(status_[iSequence], newWhere);
    1687                if (newWhere == CLP_BELOW_LOWER) {
    1688                     bound_[iSequence] = upperValue;
    1689                     upperValue = lowerValue;
    1690                     lowerValue = -COIN_DBL_MAX;
    1691                } else if (newWhere == CLP_ABOVE_UPPER) {
    1692                     bound_[iSequence] = lowerValue;
    1693                     lowerValue = upperValue;
    1694                     upperValue = COIN_DBL_MAX;
    1695                }
    1696                lower[iSequence] = lowerValue;
    1697                upper[iSequence] = upperValue;
    1698                cost[iSequence] = costValue;
    1699           }
    1700           ClpSimplex::Status status = model_->getStatus(iSequence);
    1701           if (upperValue == lowerValue) {
    1702                if (status != ClpSimplex::basic) {
    1703                     model_->setStatus(iSequence, ClpSimplex::isFixed);
    1704                     status = ClpSimplex::basic; // so will skip
    1705                }
    1706           }
    1707           switch(status) {
    1708 
    1709           case ClpSimplex::basic:
    1710           case ClpSimplex::superBasic:
    1711           case ClpSimplex::isFree:
    1712                break;
    1713           case ClpSimplex::atUpperBound:
    1714           case ClpSimplex::atLowerBound:
    1715           case ClpSimplex::isFixed:
    1716                // set correctly
    1717                if (fabs(value - lowerValue) <= primalTolerance * 1.001) {
    1718                     model_->setStatus(iSequence, ClpSimplex::atLowerBound);
    1719                } else if (fabs(value - upperValue) <= primalTolerance * 1.001) {
    1720                     model_->setStatus(iSequence, ClpSimplex::atUpperBound);
    1721                } else {
    1722                     // set superBasic
    1723                     model_->setStatus(iSequence, ClpSimplex::superBasic);
    1724                }
    1725                break;
    1726           }
     1631    clpUserStruct info;
     1632    info.type = 3;
     1633    info.sequence = iSequence;
     1634    info.value = value;
     1635    info.change = 0;
     1636    info.lower = lowerValue;
     1637    info.upper = upperValue;
     1638    info.cost = costValue;
     1639#endif
     1640    if (iWhere == CLP_BELOW_LOWER) {
     1641      lowerValue = upperValue;
     1642      upperValue = bound_[iSequence];
     1643      numberInfeasibilities_--;
     1644      assert(fabs(lowerValue) < 1.0e100);
     1645    } else if (iWhere == CLP_ABOVE_UPPER) {
     1646      upperValue = lowerValue;
     1647      lowerValue = bound_[iSequence];
     1648      numberInfeasibilities_--;
     1649    }
     1650    // get correct place
     1651    int newWhere = CLP_FEASIBLE;
     1652    if (value - upperValue <= primalTolerance) {
     1653      if (value - lowerValue >= -primalTolerance) {
     1654        // feasible
     1655        //newWhere=CLP_FEASIBLE;
     1656      } else {
     1657        // below
     1658        newWhere = CLP_BELOW_LOWER;
     1659        costValue -= infeasibilityWeight_;
     1660        numberInfeasibilities_++;
     1661        assert(fabs(lowerValue) < 1.0e100);
     1662      }
     1663    } else {
     1664      // above
     1665      newWhere = CLP_ABOVE_UPPER;
     1666      costValue += infeasibilityWeight_;
     1667      numberInfeasibilities_++;
     1668    }
     1669    if (iWhere != newWhere) {
     1670      difference = cost[iSequence] - costValue;
     1671      setOriginalStatus(status_[iSequence], newWhere);
     1672      if (newWhere == CLP_BELOW_LOWER) {
     1673        bound_[iSequence] = upperValue;
     1674        upperValue = lowerValue;
     1675        lowerValue = -COIN_DBL_MAX;
     1676      } else if (newWhere == CLP_ABOVE_UPPER) {
     1677        bound_[iSequence] = lowerValue;
     1678        lowerValue = upperValue;
     1679        upperValue = COIN_DBL_MAX;
     1680      }
     1681      lower[iSequence] = lowerValue;
     1682      upper[iSequence] = upperValue;
     1683      cost[iSequence] = costValue;
     1684    }
     1685    ClpSimplex::Status status = model_->getStatus(iSequence);
     1686    if (upperValue == lowerValue) {
     1687      if (status != ClpSimplex::basic) {
     1688        model_->setStatus(iSequence, ClpSimplex::isFixed);
     1689        status = ClpSimplex::basic; // so will skip
     1690      }
     1691    }
     1692    switch (status) {
     1693
     1694    case ClpSimplex::basic:
     1695    case ClpSimplex::superBasic:
     1696    case ClpSimplex::isFree:
     1697      break;
     1698    case ClpSimplex::atUpperBound:
     1699    case ClpSimplex::atLowerBound:
     1700    case ClpSimplex::isFixed:
     1701      // set correctly
     1702      if (fabs(value - lowerValue) <= primalTolerance * 1.001) {
     1703        model_->setStatus(iSequence, ClpSimplex::atLowerBound);
     1704      } else if (fabs(value - upperValue) <= primalTolerance * 1.001) {
     1705        model_->setStatus(iSequence, ClpSimplex::atUpperBound);
     1706      } else {
     1707        // set superBasic
     1708        model_->setStatus(iSequence, ClpSimplex::superBasic);
     1709      }
     1710      break;
     1711    }
    17271712#ifdef CLP_USER_DRIVEN
    1728           model_->eventHandler()->eventWithInfo(ClpEventHandler::pivotRow,
    1729                                                            &info);
    1730           if (info.change) {
    1731           }
    1732 #endif             
    1733      }
    1734      changeCost_ += value * difference;
    1735      return difference;
     1713    model_->eventHandler()->eventWithInfo(ClpEventHandler::pivotRow,
     1714      &info);
     1715    if (info.change) {
     1716    }
     1717#endif
     1718  }
     1719  changeCost_ += value * difference;
     1720  return difference;
    17361721}
    17371722/* Sets bounds and infeasible cost and true cost for one variable
    17381723   This is for gub and column generation etc */
    1739 void
    1740 ClpNonLinearCost::setOne(int sequence, double solutionValue, double lowerValue, double upperValue,
    1741                          double costValue)
    1742 {
    1743      if (CLP_METHOD1) {
    1744           int iRange = -1;
    1745           int start = start_[sequence];
    1746           double infeasibilityCost = model_->infeasibilityCost();
    1747           cost_[start] = costValue - infeasibilityCost;
    1748           lower_[start+1] = lowerValue;
    1749           cost_[start+1] = costValue;
    1750           lower_[start+2] = upperValue;
    1751           cost_[start+2] = costValue + infeasibilityCost;
    1752           double primalTolerance = model_->currentPrimalTolerance();
    1753           if (solutionValue - lowerValue >= -primalTolerance) {
    1754                if (solutionValue - upperValue <= primalTolerance) {
    1755                     iRange = start + 1;
    1756                } else {
    1757                     iRange = start + 2;
    1758                }
    1759           } else {
    1760                iRange = start;
    1761           }
    1762           model_->costRegion()[sequence] = cost_[iRange];
    1763           whichRange_[sequence] = iRange;
    1764      }
    1765      if (CLP_METHOD2) {
    1766        bound_[sequence]=0.0;
    1767        cost2_[sequence]=costValue;
    1768        setInitialStatus(status_[sequence]);
    1769      }
    1770 
     1724void ClpNonLinearCost::setOne(int sequence, double solutionValue, double lowerValue, double upperValue,
     1725  double costValue)
     1726{
     1727  if (CLP_METHOD1) {
     1728    int iRange = -1;
     1729    int start = start_[sequence];
     1730    double infeasibilityCost = model_->infeasibilityCost();
     1731    cost_[start] = costValue - infeasibilityCost;
     1732    lower_[start + 1] = lowerValue;
     1733    cost_[start + 1] = costValue;
     1734    lower_[start + 2] = upperValue;
     1735    cost_[start + 2] = costValue + infeasibilityCost;
     1736    double primalTolerance = model_->currentPrimalTolerance();
     1737    if (solutionValue - lowerValue >= -primalTolerance) {
     1738      if (solutionValue - upperValue <= primalTolerance) {
     1739        iRange = start + 1;
     1740      } else {
     1741        iRange = start + 2;
     1742      }
     1743    } else {
     1744      iRange = start;
     1745    }
     1746    model_->costRegion()[sequence] = cost_[iRange];
     1747    whichRange_[sequence] = iRange;
     1748  }
     1749  if (CLP_METHOD2) {
     1750    bound_[sequence] = 0.0;
     1751    cost2_[sequence] = costValue;
     1752    setInitialStatus(status_[sequence]);
     1753  }
    17711754}
    17721755/* Sets bounds and cost for outgoing variable
    17731756   may change value
    17741757   Returns direction */
    1775 int
    1776 ClpNonLinearCost::setOneOutgoing(int iSequence, double & value)
    1777 {
    1778      assert (model_ != NULL);
    1779      double primalTolerance = model_->currentPrimalTolerance();
    1780      // difference in cost
    1781      double difference = 0.0;
    1782      int direction = 0;
     1758int ClpNonLinearCost::setOneOutgoing(int iSequence, double &value)
     1759{
     1760  assert(model_ != NULL);
     1761  double primalTolerance = model_->currentPrimalTolerance();
     1762  // difference in cost
     1763  double difference = 0.0;
     1764  int direction = 0;
    17831765#ifdef CLP_USER_DRIVEN
    1784      double saveLower = model_->lowerRegion()[iSequence];
    1785      double saveUpper = model_->upperRegion()[iSequence];
    1786      double saveCost = model_->costRegion()[iSequence];
    1787 #endif
    1788      if (CLP_METHOD1) {
    1789           // get where in bound sequence
    1790           int iRange;
    1791           int currentRange = whichRange_[iSequence];
    1792           int start = start_[iSequence];
    1793           int end = start_[iSequence+1] - 1;
    1794           // Set perceived direction out
    1795           if (value <= lower_[currentRange] + 1.001 * primalTolerance) {
    1796                direction = 1;
    1797           } else if (value >= lower_[currentRange+1] - 1.001 * primalTolerance) {
    1798                direction = -1;
    1799           } else {
    1800                // odd
    1801                direction = 0;
     1766  double saveLower = model_->lowerRegion()[iSequence];
     1767  double saveUpper = model_->upperRegion()[iSequence];
     1768  double saveCost = model_->costRegion()[iSequence];
     1769#endif
     1770  if (CLP_METHOD1) {
     1771    // get where in bound sequence
     1772    int iRange;
     1773    int currentRange = whichRange_[iSequence];
     1774    int start = start_[iSequence];
     1775    int end = start_[iSequence + 1] - 1;
     1776    // Set perceived direction out
     1777    if (value <= lower_[currentRange] + 1.001 * primalTolerance) {
     1778      direction = 1;
     1779    } else if (value >= lower_[currentRange + 1] - 1.001 * primalTolerance) {
     1780      direction = -1;
     1781    } else {
     1782      // odd
     1783      direction = 0;
     1784    }
     1785    // If fixed try and get feasible
     1786    if (lower_[start + 1] == lower_[start + 2] && fabs(value - lower_[start + 1]) < 1.001 * primalTolerance) {
     1787      iRange = start + 1;
     1788    } else {
     1789      // See if exact
     1790      for (iRange = start; iRange < end; iRange++) {
     1791        if (value == lower_[iRange + 1]) {
     1792          // put in better range
     1793          if (infeasible(iRange) && iRange == start)
     1794            iRange++;
     1795          break;
     1796        }
     1797      }
     1798      if (iRange == end) {
     1799        // not exact
     1800        for (iRange = start; iRange < end; iRange++) {
     1801          if (value <= lower_[iRange + 1] + primalTolerance) {
     1802            // put in better range
     1803            if (value >= lower_[iRange + 1] - primalTolerance && infeasible(iRange) && iRange == start)
     1804              iRange++;
     1805            break;
    18021806          }
    1803           // If fixed try and get feasible
    1804           if (lower_[start+1] == lower_[start+2] && fabs(value - lower_[start+1]) < 1.001 * primalTolerance) {
    1805                iRange = start + 1;
    1806           } else {
    1807                // See if exact
    1808                for (iRange = start; iRange < end; iRange++) {
    1809                     if (value == lower_[iRange+1]) {
    1810                          // put in better range
    1811                          if (infeasible(iRange) && iRange == start)
    1812                               iRange++;
    1813                          break;
    1814                     }
    1815                }
    1816                if (iRange == end) {
    1817                     // not exact
    1818                     for (iRange = start; iRange < end; iRange++) {
    1819                          if (value <= lower_[iRange+1] + primalTolerance) {
    1820                               // put in better range
    1821                               if (value >= lower_[iRange+1] - primalTolerance && infeasible(iRange) && iRange == start)
    1822                                    iRange++;
    1823                               break;
    1824                          }
    1825                     }
    1826                }
    1827           }
    1828           assert(iRange < end);
    1829           whichRange_[iSequence] = iRange;
    1830           if (iRange != currentRange) {
    1831                if (infeasible(iRange))
    1832                     numberInfeasibilities_++;
    1833                if (infeasible(currentRange))
    1834                     numberInfeasibilities_--;
    1835           }
    1836           double & lower = model_->lowerAddress(iSequence);
    1837           double & upper = model_->upperAddress(iSequence);
    1838           double & cost = model_->costAddress(iSequence);
    1839           lower = lower_[iRange];
    1840           upper = lower_[iRange+1];
    1841           if (upper == lower) {
    1842                value = upper;
    1843           } else {
    1844                // set correctly
    1845                if (fabs(value - lower) <= primalTolerance * 1.001) {
    1846                     value = CoinMin(value, lower + primalTolerance);
    1847                } else if (fabs(value - upper) <= primalTolerance * 1.001) {
    1848                     value = CoinMax(value, upper - primalTolerance);
    1849                } else {
    1850                     //printf("*** variable wandered off bound %g %g %g!\n",
    1851                     //     lower,value,upper);
    1852                     if (value - lower <= upper - value)
    1853                          value = lower + primalTolerance;
    1854                     else
    1855                          value = upper - primalTolerance;
    1856                }
    1857           }
    1858           difference = cost - cost_[iRange];
    1859           cost = cost_[iRange];
    1860      }
    1861      if (CLP_METHOD2) {
    1862           double * upper = model_->upperRegion();
    1863           double * lower = model_->lowerRegion();
    1864           double * cost = model_->costRegion();
    1865           unsigned char iStatus = status_[iSequence];
    1866           assert (currentStatus(iStatus) == CLP_SAME);
    1867           double lowerValue = lower[iSequence];
    1868           double upperValue = upper[iSequence];
    1869           double costValue = cost2_[iSequence];
    1870           // Set perceived direction out
    1871           if (value <= lowerValue + 1.001 * primalTolerance) {
    1872                direction = 1;
    1873           } else if (value >= upperValue - 1.001 * primalTolerance) {
    1874                direction = -1;
    1875           } else {
    1876                // odd
    1877                direction = 0;
    1878           }
    1879           int iWhere = originalStatus(iStatus);
    1880           if (iWhere == CLP_BELOW_LOWER) {
    1881                lowerValue = upperValue;
    1882                upperValue = bound_[iSequence];
    1883                numberInfeasibilities_--;
    1884                assert (fabs(lowerValue) < 1.0e100);
    1885           } else if (iWhere == CLP_ABOVE_UPPER) {
    1886                upperValue = lowerValue;
    1887                lowerValue = bound_[iSequence];
    1888                numberInfeasibilities_--;
    1889           }
    1890           // get correct place
    1891           // If fixed give benefit of doubt
    1892           if (lowerValue == upperValue)
    1893                value = lowerValue;
    1894           int newWhere = CLP_FEASIBLE;
    1895           if (value - upperValue <= primalTolerance) {
    1896                if (value - lowerValue >= -primalTolerance) {
    1897                     // feasible
    1898                     //newWhere=CLP_FEASIBLE;
    1899                } else {
    1900                     // below
    1901                     newWhere = CLP_BELOW_LOWER;
    1902                     costValue -= infeasibilityWeight_;
    1903                     numberInfeasibilities_++;
    1904                     assert (fabs(lowerValue) < 1.0e100);
    1905                }
    1906           } else {
    1907                // above
    1908                newWhere = CLP_ABOVE_UPPER;
    1909                costValue += infeasibilityWeight_;
    1910                numberInfeasibilities_++;
    1911           }
    1912           if (iWhere != newWhere) {
    1913                difference = cost[iSequence] - costValue;
    1914                setOriginalStatus(status_[iSequence], newWhere);
    1915                if (newWhere == CLP_BELOW_LOWER) {
    1916                     bound_[iSequence] = upperValue;
    1917                     upper[iSequence] = lowerValue;
    1918                     lower[iSequence] = -COIN_DBL_MAX;
    1919                } else if (newWhere == CLP_ABOVE_UPPER) {
    1920                     bound_[iSequence] = lowerValue;
    1921                     lower[iSequence] = upperValue;
    1922                     upper[iSequence] = COIN_DBL_MAX;
    1923                } else {
    1924                     lower[iSequence] = lowerValue;
    1925                     upper[iSequence] = upperValue;
    1926                }
    1927                cost[iSequence] = costValue;
    1928           }
    1929           // set correctly
    1930           if (fabs(value - lowerValue) <= primalTolerance * 1.001) {
    1931                value = CoinMin(value, lowerValue + primalTolerance);
    1932           } else if (fabs(value - upperValue) <= primalTolerance * 1.001) {
    1933                value = CoinMax(value, upperValue - primalTolerance);
    1934           } else {
    1935                //printf("*** variable wandered off bound %g %g %g!\n",
    1936                //     lowerValue,value,upperValue);
    1937                if (value - lowerValue <= upperValue - value)
    1938                     value = lowerValue + primalTolerance;
    1939                else
    1940                     value = upperValue - primalTolerance;
    1941           }
    1942      }
    1943      changeCost_ += value * difference;
     1807        }
     1808      }
     1809    }
     1810    assert(iRange < end);
     1811    whichRange_[iSequence] = iRange;
     1812    if (iRange != currentRange) {
     1813      if (infeasible(iRange))
     1814        numberInfeasibilities_++;
     1815      if (infeasible(currentRange))
     1816        numberInfeasibilities_--;
     1817    }
     1818    double &lower = model_->lowerAddress(iSequence);
     1819    double &upper = model_->upperAddress(iSequence);
     1820    double &cost = model_->costAddress(iSequence);
     1821    lower = lower_[iRange];
     1822    upper = lower_[iRange + 1];
     1823    if (upper == lower) {
     1824      value = upper;
     1825    } else {
     1826      // set correctly
     1827      if (fabs(value - lower) <= primalTolerance * 1.001) {
     1828        value = CoinMin(value, lower + primalTolerance);
     1829      } else if (fabs(value - upper) <= primalTolerance * 1.001) {
     1830        value = CoinMax(value, upper - primalTolerance);
     1831      } else {
     1832        //printf("*** variable wandered off bound %g %g %g!\n",
     1833        //     lower,value,upper);
     1834        if (value - lower <= upper - value)
     1835          value = lower + primalTolerance;
     1836        else
     1837          value = upper - primalTolerance;
     1838      }
     1839    }
     1840    difference = cost - cost_[iRange];
     1841    cost = cost_[iRange];
     1842  }
     1843  if (CLP_METHOD2) {
     1844    double *upper = model_->upperRegion();
     1845    double *lower = model_->lowerRegion();
     1846    double *cost = model_->costRegion();
     1847    unsigned char iStatus = status_[iSequence];
     1848    assert(currentStatus(iStatus) == CLP_SAME);
     1849    double lowerValue = lower[iSequence];
     1850    double upperValue = upper[iSequence];
     1851    double costValue = cost2_[iSequence];
     1852    // Set perceived direction out
     1853    if (value <= lowerValue + 1.001 * primalTolerance) {
     1854      direction = 1;
     1855    } else if (value >= upperValue - 1.001 * primalTolerance) {
     1856      direction = -1;
     1857    } else {
     1858      // odd
     1859      direction = 0;
     1860    }
     1861    int iWhere = originalStatus(iStatus);
     1862    if (iWhere == CLP_BELOW_LOWER) {
     1863      lowerValue = upperValue;
     1864      upperValue = bound_[iSequence];
     1865      numberInfeasibilities_--;
     1866      assert(fabs(lowerValue) < 1.0e100);
     1867    } else if (iWhere == CLP_ABOVE_UPPER) {
     1868      upperValue = lowerValue;
     1869      lowerValue = bound_[iSequence];
     1870      numberInfeasibilities_--;
     1871    }
     1872    // get correct place
     1873    // If fixed give benefit of doubt
     1874    if (lowerValue == upperValue)
     1875      value = lowerValue;
     1876    int newWhere = CLP_FEASIBLE;
     1877    if (value - upperValue <= primalTolerance) {
     1878      if (value - lowerValue >= -primalTolerance) {
     1879        // feasible
     1880        //newWhere=CLP_FEASIBLE;
     1881      } else {
     1882        // below
     1883        newWhere = CLP_BELOW_LOWER;
     1884        costValue -= infeasibilityWeight_;
     1885        numberInfeasibilities_++;
     1886        assert(fabs(lowerValue) < 1.0e100);
     1887      }
     1888    } else {
     1889      // above
     1890      newWhere = CLP_ABOVE_UPPER;
     1891      costValue += infeasibilityWeight_;
     1892      numberInfeasibilities_++;
     1893    }
     1894    if (iWhere != newWhere) {
     1895      difference = cost[iSequence] - costValue;
     1896      setOriginalStatus(status_[iSequence], newWhere);
     1897      if (newWhere == CLP_BELOW_LOWER) {
     1898        bound_[iSequence] = upperValue;
     1899        upper[iSequence] = lowerValue;
     1900        lower[iSequence] = -COIN_DBL_MAX;
     1901      } else if (newWhere == CLP_ABOVE_UPPER) {
     1902        bound_[iSequence] = lowerValue;
     1903        lower[iSequence] = upperValue;
     1904        upper[iSequence] = COIN_DBL_MAX;
     1905      } else {
     1906        lower[iSequence] = lowerValue;
     1907        upper[iSequence] = upperValue;
     1908      }
     1909      cost[iSequence] = costValue;
     1910    }
     1911    // set correctly
     1912    if (fabs(value - lowerValue) <= primalTolerance * 1.001) {
     1913      value = CoinMin(value, lowerValue + primalTolerance);
     1914    } else if (fabs(value - upperValue) <= primalTolerance * 1.001) {
     1915      value = CoinMax(value, upperValue - primalTolerance);
     1916    } else {
     1917      //printf("*** variable wandered off bound %g %g %g!\n",
     1918      //     lowerValue,value,upperValue);
     1919      if (value - lowerValue <= upperValue - value)
     1920        value = lowerValue + primalTolerance;
     1921      else
     1922        value = upperValue - primalTolerance;
     1923    }
     1924  }
     1925  changeCost_ += value * difference;
    19441926#ifdef CLP_USER_DRIVEN
    1945      //if (iSequence>=model_->numberColumns)
    1946      if (difference||saveCost!=model_->costRegion()[iSequence]) {
    1947        printf("Out %d (%d) %g <= %g <= %g cost %g x=>x %g < %g cost %g\n",
    1948               iSequence,iSequence-model_->numberColumns(),
    1949               saveLower,value,saveUpper,saveCost,
    1950               model_->lowerRegion()[iSequence],
    1951               model_->upperRegion()[iSequence],
    1952               model_->costRegion()[iSequence]);
    1953      }
    1954 #endif
    1955      return direction;
     1927  //if (iSequence>=model_->numberColumns)
     1928  if (difference || saveCost != model_->costRegion()[iSequence]) {
     1929    printf("Out %d (%d) %g <= %g <= %g cost %g x=>x %g < %g cost %g\n",
     1930      iSequence, iSequence - model_->numberColumns(),
     1931      saveLower, value, saveUpper, saveCost,
     1932      model_->lowerRegion()[iSequence],
     1933      model_->upperRegion()[iSequence],
     1934      model_->costRegion()[iSequence]);
     1935  }
     1936#endif
     1937  return direction;
    19561938}
    19571939// Returns nearest bound
     
    19591941ClpNonLinearCost::nearest(int iSequence, double solutionValue)
    19601942{
    1961      assert (model_ != NULL);
    1962      double nearest = 0.0;
    1963      if (CLP_METHOD1) {
    1964           // get where in bound sequence
    1965           int iRange;
    1966           int start = start_[iSequence];
    1967           int end = start_[iSequence+1];
    1968           int jRange = -1;
    1969           nearest = COIN_DBL_MAX;
    1970           for (iRange = start; iRange < end; iRange++) {
    1971                if (fabs(solutionValue - lower_[iRange]) < nearest) {
    1972                     jRange = iRange;
    1973                     nearest = fabs(solutionValue - lower_[iRange]);
    1974                }
    1975           }
    1976           assert(jRange < end);
    1977           nearest = lower_[jRange];
    1978      }
    1979      if (CLP_METHOD2) {
    1980           const double * upper = model_->upperRegion();
    1981           const double * lower = model_->lowerRegion();
    1982           double lowerValue = lower[iSequence];
    1983           double upperValue = upper[iSequence];
    1984           int iWhere = originalStatus(status_[iSequence]);
    1985           if (iWhere == CLP_BELOW_LOWER) {
    1986                lowerValue = upperValue;
    1987                upperValue = bound_[iSequence];
    1988                assert (fabs(lowerValue) < 1.0e100);
    1989           } else if (iWhere == CLP_ABOVE_UPPER) {
    1990                upperValue = lowerValue;
    1991                lowerValue = bound_[iSequence];
    1992           }
    1993           if (fabs(solutionValue - lowerValue) < fabs(solutionValue - upperValue))
    1994                nearest = lowerValue;
    1995           else
    1996                nearest = upperValue;
    1997      }
    1998      return nearest;
     1943  assert(model_ != NULL);
     1944  double nearest = 0.0;
     1945  if (CLP_METHOD1) {
     1946    // get where in bound sequence
     1947    int iRange;
     1948    int start = start_[iSequence];
     1949    int end = start_[iSequence + 1];
     1950    int jRange = -1;
     1951    nearest = COIN_DBL_MAX;
     1952    for (iRange = start; iRange < end; iRange++) {
     1953      if (fabs(solutionValue - lower_[iRange]) < nearest) {
     1954        jRange = iRange;
     1955        nearest = fabs(solutionValue - lower_[iRange]);
     1956      }
     1957    }
     1958    assert(jRange < end);
     1959    nearest = lower_[jRange];
     1960  }
     1961  if (CLP_METHOD2) {
     1962    const double *upper = model_->upperRegion();
     1963    const double *lower = model_->lowerRegion();
     1964    double lowerValue = lower[iSequence];
     1965    double upperValue = upper[iSequence];
     1966    int iWhere = originalStatus(status_[iSequence]);
     1967    if (iWhere == CLP_BELOW_LOWER) {
     1968      lowerValue = upperValue;
     1969      upperValue = bound_[iSequence];
     1970      assert(fabs(lowerValue) < 1.0e100);
     1971    } else if (iWhere == CLP_ABOVE_UPPER) {
     1972      upperValue = lowerValue;
     1973      lowerValue = bound_[iSequence];
     1974    }
     1975    if (fabs(solutionValue - lowerValue) < fabs(solutionValue - upperValue))
     1976      nearest = lowerValue;
     1977    else
     1978      nearest = upperValue;
     1979  }
     1980  return nearest;
    19991981}
    20001982/// Feasible cost with offset and direction (i.e. for reporting)
     
    20021984ClpNonLinearCost::feasibleReportCost() const
    20031985{
    2004      double value;
    2005      model_->getDblParam(ClpObjOffset, value);
    2006      return (feasibleCost_ + model_->objectiveAsObject()->nonlinearOffset()) * model_->optimizationDirection() /
    2007             (model_->objectiveScale() * model_->rhsScale()) - value;
     1986  double value;
     1987  model_->getDblParam(ClpObjOffset, value);
     1988  return (feasibleCost_ + model_->objectiveAsObject()->nonlinearOffset()) * model_->optimizationDirection() / (model_->objectiveScale() * model_->rhsScale()) - value;
    20081989}
    20091990// Get rid of real costs (just for moment)
    2010 void
    2011 ClpNonLinearCost::zapCosts()
    2012 {
    2013      int iSequence;
    2014      double infeasibilityCost = model_->infeasibilityCost();
    2015      // zero out all costs
    2016      int numberTotal = numberColumns_ + numberRows_;
    2017      if (CLP_METHOD1) {
    2018           int n = start_[numberTotal];
    2019           memset(cost_, 0, n * sizeof(double));
    2020           for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    2021                int start = start_[iSequence];
    2022                int end = start_[iSequence+1] - 1;
    2023                // correct costs for this infeasibility weight
    2024                if (infeasible(start)) {
    2025                     cost_[start] = -infeasibilityCost;
    2026                }
    2027                if (infeasible(end - 1)) {
    2028                     cost_[end-1] = infeasibilityCost;
    2029                }
    2030           }
    2031      }
    2032      if (CLP_METHOD2) {
    2033      }
     1991void ClpNonLinearCost::zapCosts()
     1992{
     1993  int iSequence;
     1994  double infeasibilityCost = model_->infeasibilityCost();
     1995  // zero out all costs
     1996  int numberTotal = numberColumns_ + numberRows_;
     1997  if (CLP_METHOD1) {
     1998    int n = start_[numberTotal];
     1999    memset(cost_, 0, n * sizeof(double));
     2000    for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     2001      int start = start_[iSequence];
     2002      int end = start_[iSequence + 1] - 1;
     2003      // correct costs for this infeasibility weight
     2004      if (infeasible(start)) {
     2005        cost_[start] = -infeasibilityCost;
     2006      }
     2007      if (infeasible(end - 1)) {
     2008        cost_[end - 1] = infeasibilityCost;
     2009      }
     2010    }
     2011  }
     2012  if (CLP_METHOD2) {
     2013  }
    20342014}
    20352015#ifdef VALIDATE
    20362016// For debug
    2037 void
    2038 ClpNonLinearCost::validate()
    2039 {
    2040      double primalTolerance = model_->currentPrimalTolerance();
    2041      int iSequence;
    2042      const double * solution = model_->solutionRegion();
    2043      const double * upper = model_->upperRegion();
    2044      const double * lower = model_->lowerRegion();
    2045      const double * cost = model_->costRegion();
    2046      double infeasibilityCost = model_->infeasibilityCost();
    2047      int numberTotal = numberRows_ + numberColumns_;
    2048      int numberInfeasibilities = 0;
    2049      double sumInfeasibilities = 0.0;
    2050 
    2051      for (iSequence = 0; iSequence < numberTotal; iSequence++) {
    2052           double value = solution[iSequence];
    2053           int iStatus = status_[iSequence];
    2054           assert (currentStatus(iStatus) == CLP_SAME);
    2055           double lowerValue = lower[iSequence];
    2056           double upperValue = upper[iSequence];
    2057           double costValue = cost2_[iSequence];
    2058           int iWhere = originalStatus(iStatus);
    2059           if (iWhere == CLP_BELOW_LOWER) {
    2060                lowerValue = upperValue;
    2061                upperValue = bound_[iSequence];
    2062                assert (fabs(lowerValue) < 1.0e100);
    2063                costValue -= infeasibilityCost;
    2064                assert (value <= lowerValue - primalTolerance);
    2065                numberInfeasibilities++;
    2066                sumInfeasibilities += lowerValue - value - primalTolerance;
    2067                assert (model_->getStatus(iSequence) == ClpSimplex::basic);
    2068           } else if (iWhere == CLP_ABOVE_UPPER) {
    2069                upperValue = lowerValue;
    2070                lowerValue = bound_[iSequence];
    2071                costValue += infeasibilityCost;
    2072                assert (value >= upperValue + primalTolerance);
    2073                numberInfeasibilities++;
    2074                sumInfeasibilities += value - upperValue - primalTolerance;
    2075                assert (model_->getStatus(iSequence) == ClpSimplex::basic);
    2076           } else {
    2077                assert (value >= lowerValue - primalTolerance && value <= upperValue + primalTolerance);
    2078           }
    2079           assert (lowerValue == saveLowerV[iSequence]);
    2080           assert (upperValue == saveUpperV[iSequence]);
    2081           assert (costValue == cost[iSequence]);
    2082      }
    2083      if (numberInfeasibilities)
    2084           printf("JJ %d infeasibilities summing to %g\n",
    2085                  numberInfeasibilities, sumInfeasibilities);
    2086 }
    2087 #endif
     2017void ClpNonLinearCost::validate()
     2018{
     2019  double primalTolerance = model_->currentPrimalTolerance();
     2020  int iSequence;
     2021  const double *solution = model_->solutionRegion();
     2022  const double *upper = model_->upperRegion();
     2023  const double *lower = model_->lowerRegion();
     2024  const double *cost = model_->costRegion();
     2025  double infeasibilityCost = model_->infeasibilityCost();
     2026  int numberTotal = numberRows_ + numberColumns_;
     2027  int numberInfeasibilities = 0;
     2028  double sumInfeasibilities = 0.0;
     2029
     2030  for (iSequence = 0; iSequence < numberTotal; iSequence++) {
     2031    double value = solution[iSequence];
     2032    int iStatus = status_[iSequence];
     2033    assert(currentStatus(iStatus) == CLP_SAME);
     2034    double lowerValue = lower[iSequence];
     2035    double upperValue = upper[iSequence];
     2036    double costValue = cost2_[iSequence];
     2037    int iWhere = originalStatus(iStatus);
     2038    if (iWhere == CLP_BELOW_LOWER) {
     2039      lowerValue = upperValue;
     2040      upperValue = bound_[iSequence];
     2041      assert(fabs(lowerValue) < 1.0e100);
     2042      costValue -= infeasibilityCost;
     2043      assert(value <= lowerValue - primalTolerance);
     2044      numberInfeasibilities++;
     2045      sumInfeasibilities += lowerValue - value - primalTolerance;
     2046      assert(model_->getStatus(iSequence) == ClpSimplex::basic);
     2047    } else if (iWhere == CLP_ABOVE_UPPER) {
     2048      upperValue = lowerValue;
     2049      lowerValue = bound_[iSequence];
     2050      costValue += infeasibilityCost;
     2051      assert(value >= upperValue + primalTolerance);
     2052      numberInfeasibilities++;
     2053      sumInfeasibilities += value - upperValue - primalTolerance;
     2054      assert(model_->getStatus(iSequence) == ClpSimplex::basic);
     2055    } else {
     2056      assert(value >= lowerValue - primalTolerance && value <= upperValue + primalTolerance);
     2057    }
     2058    assert(lowerValue == saveLowerV[iSequence]);
     2059    assert(upperValue == saveUpperV[iSequence]);
     2060    assert(costValue == cost[iSequence]);
     2061  }
     2062  if (numberInfeasibilities)
     2063    printf("JJ %d infeasibilities summing to %g\n",
     2064      numberInfeasibilities, sumInfeasibilities);
     2065}
     2066#endif
     2067
     2068/* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
     2069*/
Note: See TracChangeset for help on using the changeset viewer.