source: trunk/Clp/src/ClpNonLinearCost.cpp @ 1264

Last change on this file since 1264 was 1264, checked in by forrest, 11 years ago

BSP changes from 1247 to 1259

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