source: trunk/Clp/src/OsiClp/OsiClpSolverInterface.cpp @ 2030

Last change on this file since 2030 was 2030, checked in by forrest, 6 years ago

fix some ampl stuff, add ClpSolver? and a few fixes

File size: 302.6 KB
Line 
1// $Id$
2// Copyright (C) 2002, International Business Machines
3// Corporation and others.  All Rights Reserved.
4// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6#include <cassert>
7#ifdef CBC_STATISTICS
8extern int osi_crunch;
9extern int osi_primal;
10extern int osi_dual;
11extern int osi_hot;
12#endif
13#include "CoinTime.hpp"
14#include "CoinHelperFunctions.hpp"
15#include "CoinIndexedVector.hpp"
16#include "CoinModel.hpp"
17#include "CoinMpsIO.hpp"
18#include "CoinSort.hpp"
19#include "ClpDualRowSteepest.hpp"
20#include "ClpPrimalColumnSteepest.hpp"
21#include "ClpPackedMatrix.hpp"
22#include "ClpDualRowDantzig.hpp"
23#include "ClpPrimalColumnDantzig.hpp"
24#include "ClpFactorization.hpp"
25#include "ClpObjective.hpp"
26#include "ClpSimplex.hpp"
27#include "ClpSimplexOther.hpp"
28#include "ClpSimplexPrimal.hpp"
29#include "ClpSimplexDual.hpp"
30#include "ClpNonLinearCost.hpp"
31#include "OsiClpSolverInterface.hpp"
32#include "OsiBranchingObject.hpp"
33#include "OsiCuts.hpp"
34#include "OsiRowCut.hpp"
35#include "OsiColCut.hpp"
36#include "ClpPresolve.hpp"
37#include "CoinLpIO.hpp"
38static double totalTime=0.0;
39//#define SAVE_MODEL 1
40#ifdef SAVE_MODEL
41static int resolveTry=0;
42static int loResolveTry=0;
43static int hiResolveTry=9999999;
44#endif
45//#############################################################################
46// Solve methods
47//#############################################################################
48void OsiClpSolverInterface::initialSolve()
49{
50#define KEEP_SMALL
51#ifdef KEEP_SMALL
52  if (smallModel_) {
53    delete [] spareArrays_;
54    spareArrays_ = NULL;
55    delete smallModel_;
56    smallModel_=NULL;
57  }
58#endif
59  if ((specialOptions_&2097152)!=0||(specialOptions_&4194304)!=0) {
60    bool takeHint;
61    OsiHintStrength strength;
62    int algorithm = 0;
63    getHintParam(OsiDoDualInInitial,takeHint,strength);
64    if (strength!=OsiHintIgnore)
65      algorithm = takeHint ? -1 : 1;
66    if (algorithm>0||(specialOptions_&4194304)!=0) {
67      // Gub
68      resolveGub((9*modelPtr_->numberRows())/10);
69      return;
70    }
71  }
72  bool deleteSolver;
73  ClpSimplex * solver;
74  double time1 = CoinCpuTime();
75  int userFactorizationFrequency = modelPtr_->factorization()->maximumPivots();
76  int totalIterations=0;
77  bool abortSearch=false;
78  ClpObjective * savedObjective=NULL;
79  double savedDualLimit=modelPtr_->dblParam_[ClpDualObjectiveLimit];
80  if (fakeObjective_) {
81    // Clear (no objective, 0-1 and in B&B)
82    modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~128));
83    // See if all with costs fixed
84    int numberColumns = modelPtr_->numberColumns_;
85    const double * obj = modelPtr_->objective();
86    const double * lower = modelPtr_->columnLower();
87    const double * upper = modelPtr_->columnUpper();
88    int i;
89    for (i=0;i<numberColumns;i++) {
90      double objValue = obj[i];
91      if (objValue) {
92        if (lower[i]!=upper[i])
93          break;
94      }
95    }
96    if (i==numberColumns) {
97      // Check (Clp fast dual)
98      if ((specialOptions_&524288)==0) {
99        // Set fake
100        savedObjective=modelPtr_->objective_;
101        modelPtr_->objective_=fakeObjective_;
102        modelPtr_->dblParam_[ClpDualObjectiveLimit]=COIN_DBL_MAX;
103      } else {
104        // Set (no objective, 0-1 and in B&B)
105        modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()|128);
106      }
107    }
108  }
109  // Check (in branch and bound)
110  if ((specialOptions_&1024)==0) {
111    solver = new ClpSimplex(true);
112    deleteSolver=true;
113    solver->borrowModel(*modelPtr_);
114    // See if user set factorization frequency
115    // borrowModel does not move
116    solver->factorization()->maximumPivots(userFactorizationFrequency);
117  } else {
118    solver = modelPtr_;
119    deleteSolver=false;
120  }
121  // Treat as if user simplex not enabled
122  int saveSolveType=solver->solveType();
123  bool doingPrimal = solver->algorithm()>0;
124  if (saveSolveType==2) {
125    disableSimplexInterface();
126    solver->setSolveType(1);
127  }
128  int saveOptions = solver->specialOptions();
129  solver->setSpecialOptions(saveOptions|64|32768); // go as far as possible
130  // get original log levels
131  int saveMessageLevel=modelPtr_->logLevel();
132  int messageLevel=messageHandler()->logLevel();
133  int saveMessageLevel2 = messageLevel;
134  // Set message handler
135  if (!defaultHandler_)
136    solver->passInMessageHandler(handler_);
137  // But keep log level
138  solver->messageHandler()->setLogLevel(saveMessageLevel);
139  // set reasonable defaults
140  bool takeHint;
141  OsiHintStrength strength;
142  // Switch off printing if asked to
143  bool gotHint = (getHintParam(OsiDoReducePrint,takeHint,strength));
144  assert (gotHint);
145  if (strength!=OsiHintIgnore&&takeHint) {
146    if (messageLevel>0)
147      messageLevel--;
148  }
149  if (messageLevel<saveMessageLevel)
150    solver->messageHandler()->setLogLevel(messageLevel);
151  // Allow for specialOptions_==1+8 forcing saving factorization
152  int startFinishOptions=0;
153  if ((specialOptions_&9)==(1+8)) {
154    startFinishOptions =1+2+4; // allow re-use of factorization
155  }
156  bool defaultHints=true;
157  {
158    int hint;
159    for (hint=OsiDoPresolveInInitial;hint<OsiLastHintParam;hint++) {
160      if (hint!=OsiDoReducePrint&&
161          hint!=OsiDoInBranchAndCut) {
162        bool yesNo;
163        OsiHintStrength strength;
164        getHintParam(static_cast<OsiHintParam> (hint),yesNo,strength);
165        if (yesNo) {
166          defaultHints=false;
167          break;
168        }
169        if (strength != OsiHintIgnore) {
170          defaultHints=false;
171          break;
172        }
173      }
174    }
175  }
176  ClpPresolve * pinfo = NULL;
177  /*
178    If basis then do primal (as user could do dual with resolve)
179    If not then see if dual feasible (and allow for gubs etc?)
180  */
181  bool doPrimal = (basis_.numberBasicStructurals()>0);
182  setBasis(basis_,solver);
183  bool inCbcOrOther = (modelPtr_->specialOptions()&0x03000000)!=0;
184  if ((!defaultHints||doPrimal)&&!solveOptions_.getSpecialOption(6)) {
185    // scaling
186    // save initial state
187    const double * rowScale1 = solver->rowScale();
188    if (modelPtr_->solveType()==1) {
189      gotHint = (getHintParam(OsiDoScale,takeHint,strength));
190      assert (gotHint);
191      if (strength==OsiHintIgnore||takeHint) {
192        if (!solver->scalingFlag())
193          solver->scaling(3);
194      } else {
195        solver->scaling(0);
196      }
197    } else {
198      solver->scaling(0);
199    }
200    //solver->setDualBound(1.0e6);
201    //solver->setDualTolerance(1.0e-7);
202   
203    //ClpDualRowSteepest steep;
204    //solver->setDualRowPivotAlgorithm(steep);
205    //solver->setPrimalTolerance(1.0e-8);
206    //ClpPrimalColumnSteepest steepP;
207    //solver->setPrimalColumnPivotAlgorithm(steepP);
208   
209    // sort out hints;
210    // algorithm 0 whatever, -1 force dual, +1 force primal
211    int algorithm = 0;
212    gotHint = (getHintParam(OsiDoDualInInitial,takeHint,strength));
213    assert (gotHint);
214    if (strength!=OsiHintIgnore)
215      algorithm = takeHint ? -1 : 1;
216    // crash 0 do lightweight if all slack, 1 do, -1 don't
217    int doCrash=0;
218    gotHint = (getHintParam(OsiDoCrash,takeHint,strength));
219    assert (gotHint);
220    if (strength!=OsiHintIgnore)
221      doCrash = takeHint ? 1 : -1;
222    // doPrimal set true if any structurals in basis so switch off crash
223    if (doPrimal)
224      doCrash = -1;
225   
226    // presolve
227    gotHint = (getHintParam(OsiDoPresolveInInitial,takeHint,strength));
228    assert (gotHint);
229    if (strength!=OsiHintIgnore&&takeHint) {
230      pinfo = new ClpPresolve();
231      ClpSimplex * model2 = pinfo->presolvedModel(*solver,1.0e-8);
232      if (!model2) {
233        // problem found to be infeasible - whats best?
234        model2 = solver;
235        delete pinfo;
236        pinfo = NULL;
237      } else {
238        model2->setSpecialOptions(solver->specialOptions());
239      }
240     
241      // change from 200 (unless changed)
242      if (modelPtr_->factorization()->maximumPivots()==200)
243        model2->factorization()->maximumPivots(100+model2->numberRows()/50);
244      else
245        model2->factorization()->maximumPivots(userFactorizationFrequency);
246      int savePerturbation = model2->perturbation();
247      if (savePerturbation==100)
248        model2->setPerturbation(50);
249      if (!doPrimal) {
250        // faster if bounds tightened
251        //int numberInfeasibilities = model2->tightenPrimalBounds();
252        model2->tightenPrimalBounds();
253        // look further
254        bool crashResult=false;
255        if (doCrash>0)
256          crashResult =  (solver->crash(1000.0,1)>0);
257        else if (doCrash==0&&algorithm>0)
258          crashResult =  (solver->crash(1000.0,1)>0);
259        doPrimal=crashResult;
260      }
261      if (algorithm<0)
262        doPrimal=false;
263      else if (algorithm>0)
264        doPrimal=true;
265      if (!doPrimal) {
266        //if (numberInfeasibilities)
267        //std::cout<<"** Analysis indicates model infeasible"
268        //       <<std::endl;
269        // up dual bound for safety
270        //model2->setDualBound(1.0e11);
271        disasterHandler_->setOsiModel(this);
272        if (inCbcOrOther) {
273          disasterHandler_->setSimplex(model2);
274          disasterHandler_->setWhereFrom(4);
275          model2->setDisasterHandler(disasterHandler_);
276        }
277        model2->dual(0);
278        totalIterations += model2->numberIterations();
279        if (inCbcOrOther) {
280          if(disasterHandler_->inTrouble()) {
281#ifdef COIN_DEVELOP
282            printf("dual trouble a\n");
283#endif
284            if (disasterHandler_->typeOfDisaster()) {
285              // We want to abort
286              abortSearch=true;
287              goto disaster;
288            }
289            // try just going back in
290            disasterHandler_->setPhase(1);
291            model2->dual();
292            totalIterations += model2->numberIterations();
293            if (disasterHandler_->inTrouble()) {
294#ifdef COIN_DEVELOP
295              printf("dual trouble b\n");
296#endif
297              if (disasterHandler_->typeOfDisaster()) {
298                // We want to abort
299                abortSearch=true;
300                goto disaster;
301              }
302              // try primal with original basis
303              disasterHandler_->setPhase(2);
304              setBasis(basis_,model2);
305              model2->primal();
306              totalIterations += model2->numberIterations();
307            }
308            if(disasterHandler_->inTrouble()) {
309#ifdef COIN_DEVELOP
310              printf("disaster - treat as infeasible\n");
311#endif
312              if (disasterHandler_->typeOfDisaster()) {
313                // We want to abort
314                abortSearch=true;
315                goto disaster;
316              }
317              model2->setProblemStatus(1);
318            }
319          }
320          // reset
321          model2->setDisasterHandler(NULL);
322        }
323        // check if clp thought it was in a loop
324        if (model2->status()==3&&!model2->hitMaximumIterations()) {
325          // switch algorithm
326          disasterHandler_->setOsiModel(this);
327          if (inCbcOrOther) {
328            disasterHandler_->setSimplex(model2);
329            disasterHandler_->setWhereFrom(6);
330            model2->setDisasterHandler(disasterHandler_);
331          }
332          model2->primal();
333          totalIterations += model2->numberIterations();
334          if (inCbcOrOther) {
335            if(disasterHandler_->inTrouble()) {
336#ifdef COIN_DEVELOP
337              printf("primal trouble a\n");
338#endif
339              if (disasterHandler_->typeOfDisaster()) {
340                // We want to abort
341                abortSearch=true;
342                goto disaster;
343              }
344              // try just going back in (but with dual)
345              disasterHandler_->setPhase(1);
346              model2->dual();
347              totalIterations += model2->numberIterations();
348              if (disasterHandler_->inTrouble()) {
349#ifdef COIN_DEVELOP
350                printf("primal trouble b\n");
351#endif
352                if (disasterHandler_->typeOfDisaster()) {
353                  // We want to abort
354                  abortSearch=true;
355                  goto disaster;
356                }
357                // try primal with original basis
358                disasterHandler_->setPhase(2);
359                setBasis(basis_,model2);
360                model2->dual();
361                totalIterations += model2->numberIterations();
362              }
363              if(disasterHandler_->inTrouble()) {
364#ifdef COIN_DEVELOP
365                printf("disaster - treat as infeasible\n");
366#endif
367                if (disasterHandler_->typeOfDisaster()) {
368                  // We want to abort
369                  abortSearch=true;
370                  goto disaster;
371                }
372                model2->setProblemStatus(1);
373              }
374            }
375            // reset
376            model2->setDisasterHandler(NULL);
377          }
378        }
379      } else {
380        // up infeasibility cost for safety
381        //model2->setInfeasibilityCost(1.0e10);
382        disasterHandler_->setOsiModel(this);
383        if (inCbcOrOther) {
384          disasterHandler_->setSimplex(model2);
385          disasterHandler_->setWhereFrom(6);
386          model2->setDisasterHandler(disasterHandler_);
387        }
388        model2->primal(1);
389        totalIterations += model2->numberIterations();
390        if (inCbcOrOther) {
391          if(disasterHandler_->inTrouble()) {
392#ifdef COIN_DEVELOP
393            printf("primal trouble a\n");
394#endif
395            if (disasterHandler_->typeOfDisaster()) {
396              // We want to abort
397              abortSearch=true;
398              goto disaster;
399            }
400            // try just going back in (but with dual)
401            disasterHandler_->setPhase(1);
402            model2->dual();
403            totalIterations += model2->numberIterations();
404            if (disasterHandler_->inTrouble()) {
405#ifdef COIN_DEVELOP
406              printf("primal trouble b\n");
407#endif
408              if (disasterHandler_->typeOfDisaster()) {
409                // We want to abort
410                abortSearch=true;
411                goto disaster;
412              }
413              // try primal with original basis
414              disasterHandler_->setPhase(2);
415              setBasis(basis_,model2);
416              model2->dual();
417              totalIterations += model2->numberIterations();
418            }
419            if(disasterHandler_->inTrouble()) {
420#ifdef COIN_DEVELOP
421              printf("disaster - treat as infeasible\n");
422#endif
423              if (disasterHandler_->typeOfDisaster()) {
424                // We want to abort
425                abortSearch=true;
426                goto disaster;
427              }
428              model2->setProblemStatus(1);
429            }
430          }
431          // reset
432          model2->setDisasterHandler(NULL);
433        }
434        // check if clp thought it was in a loop
435        if (model2->status()==3&&!model2->hitMaximumIterations()) {
436          // switch algorithm
437          disasterHandler_->setOsiModel(this);
438          if (inCbcOrOther) {
439            disasterHandler_->setSimplex(model2);
440            disasterHandler_->setWhereFrom(4);
441            model2->setDisasterHandler(disasterHandler_);
442          }
443          model2->dual(0);
444          totalIterations += model2->numberIterations();
445          if (inCbcOrOther) {
446            if(disasterHandler_->inTrouble()) {
447#ifdef COIN_DEVELOP
448              printf("dual trouble a\n");
449#endif
450              if (disasterHandler_->typeOfDisaster()) {
451                // We want to abort
452                abortSearch=true;
453                goto disaster;
454              }
455              // try just going back in
456              disasterHandler_->setPhase(1);
457              model2->dual();
458              totalIterations += model2->numberIterations();
459              if (disasterHandler_->inTrouble()) {
460#ifdef COIN_DEVELOP
461                printf("dual trouble b\n");
462#endif
463                if (disasterHandler_->typeOfDisaster()) {
464                  // We want to abort
465                  abortSearch=true;
466                goto disaster;
467                }
468                // try primal with original basis
469                disasterHandler_->setPhase(2);
470                setBasis(basis_,model2);
471                model2->primal();
472                totalIterations += model2->numberIterations();
473              }
474              if(disasterHandler_->inTrouble()) {
475#ifdef COIN_DEVELOP
476                printf("disaster - treat as infeasible\n");
477#endif
478                if (disasterHandler_->typeOfDisaster()) {
479                  // We want to abort
480                  abortSearch=true;
481                  goto disaster;
482                }
483                model2->setProblemStatus(1);
484              }
485            }
486            // reset
487            model2->setDisasterHandler(NULL);
488          }
489        }
490      }
491      model2->setPerturbation(savePerturbation);
492      if (model2!=solver) {
493        int presolvedStatus = model2->status();
494        pinfo->postsolve(true);
495        delete pinfo;
496        pinfo = NULL;
497       
498        delete model2;
499        int oldStatus=solver->status();
500        solver->setProblemStatus(presolvedStatus);
501        if (solver->logLevel()==63) // for gcc 4.6 bug
502          printf("pstat %d stat %d\n",presolvedStatus,oldStatus);
503        //printf("Resolving from postsolved model\n");
504        // later try without (1) and check duals before solve
505        if (presolvedStatus!=3
506            &&(presolvedStatus||oldStatus==-1)) {
507          if (!inCbcOrOther||presolvedStatus!=1) {
508            disasterHandler_->setOsiModel(this);
509            if (inCbcOrOther) {
510              disasterHandler_->setSimplex(solver); // as "borrowed"
511              disasterHandler_->setWhereFrom(6);
512              solver->setDisasterHandler(disasterHandler_);
513            }
514            solver->primal(1);
515            totalIterations += solver->numberIterations();
516            if (inCbcOrOther) {
517              if(disasterHandler_->inTrouble()) {
518#ifdef COIN_DEVELOP
519                printf("primal trouble a\n");
520#endif
521                if (disasterHandler_->typeOfDisaster()) {
522                  // We want to abort
523                  abortSearch=true;
524                  goto disaster;
525                }
526                // try just going back in (but with dual)
527                disasterHandler_->setPhase(1);
528                solver->dual();
529                totalIterations += solver->numberIterations();
530                if (disasterHandler_->inTrouble()) {
531#ifdef COIN_DEVELOP
532                  printf("primal trouble b\n");
533#endif
534                  if (disasterHandler_->typeOfDisaster()) {
535                    // We want to abort
536                    abortSearch=true;
537                    goto disaster;
538                  }
539                  // try primal with original basis
540                  disasterHandler_->setPhase(2);
541                  setBasis(basis_,solver);
542                  solver->dual();
543                  totalIterations += solver->numberIterations();
544                }
545                if(disasterHandler_->inTrouble()) {
546#ifdef COIN_DEVELOP
547                  printf("disaster - treat as infeasible\n");
548#endif
549                  if (disasterHandler_->typeOfDisaster()) {
550                    // We want to abort
551                    abortSearch=true;
552                    goto disaster;
553                  }
554                  solver->setProblemStatus(1);
555                }
556              }
557              // reset
558              solver->setDisasterHandler(NULL);
559            }
560          }
561        }
562      }
563      lastAlgorithm_=1; // primal
564      //if (solver->numberIterations())
565      //printf("****** iterated %d\n",solver->numberIterations());
566    } else {
567      // do we want crash
568      if (doCrash>0)
569        solver->crash(1000.0,2);
570      else if (doCrash==0)
571        solver->crash(1000.0,0);
572      if (algorithm<0)
573        doPrimal=false;
574      else if (algorithm>0)
575        doPrimal=true;
576      disasterHandler_->setOsiModel(this);
577      disasterHandler_->setSimplex(solver); // as "borrowed"
578      bool inCbcOrOther = (modelPtr_->specialOptions()&0x03000000)!=0;
579      if (!doPrimal) 
580        disasterHandler_->setWhereFrom(4);
581      else
582        disasterHandler_->setWhereFrom(6);
583      if (inCbcOrOther)
584        solver->setDisasterHandler(disasterHandler_);
585      if (!doPrimal) {
586        //printf("doing dual\n");
587        solver->dual(0);
588        totalIterations += solver->numberIterations();
589        if (inCbcOrOther) {
590          if(disasterHandler_->inTrouble()) {
591#ifdef COIN_DEVELOP
592            printf("dual trouble a\n");
593#endif
594            if (disasterHandler_->typeOfDisaster()) {
595              // We want to abort
596              abortSearch=true;
597              goto disaster;
598            }
599            // try just going back in
600            disasterHandler_->setPhase(1);
601            solver->dual();
602            totalIterations += solver->numberIterations();
603            if (disasterHandler_->inTrouble()) {
604#ifdef COIN_DEVELOP
605              printf("dual trouble b\n");
606#endif
607              if (disasterHandler_->typeOfDisaster()) {
608                // We want to abort
609                abortSearch=true;
610                goto disaster;
611              }
612              // try primal with original basis
613              disasterHandler_->setPhase(2);
614              setBasis(basis_,solver);
615              solver->primal();
616              totalIterations += solver->numberIterations();
617            }
618            if(disasterHandler_->inTrouble()) {
619#ifdef COIN_DEVELOP
620              printf("disaster - treat as infeasible\n");
621#endif
622              if (disasterHandler_->typeOfDisaster()) {
623                // We want to abort
624                abortSearch=true;
625                goto disaster;
626              }
627              solver->setProblemStatus(1);
628            }
629          }
630          // reset
631          solver->setDisasterHandler(NULL);
632        }
633        lastAlgorithm_=2; // dual
634        // check if clp thought it was in a loop
635        if (solver->status()==3&&!solver->hitMaximumIterations()) {
636          // switch algorithm
637          solver->primal(0);
638          totalIterations += solver->numberIterations();
639          lastAlgorithm_=1; // primal
640        }
641      } else {
642        //printf("doing primal\n");
643        solver->primal(1);
644        totalIterations += solver->numberIterations();
645        if (inCbcOrOther) {
646          if(disasterHandler_->inTrouble()) {
647#ifdef COIN_DEVELOP
648            printf("primal trouble a\n");
649#endif
650            if (disasterHandler_->typeOfDisaster()) {
651              // We want to abort
652              abortSearch=true;
653              goto disaster;
654            }
655            // try just going back in (but with dual)
656            disasterHandler_->setPhase(1);
657            solver->dual();
658            totalIterations += solver->numberIterations();
659            if (disasterHandler_->inTrouble()) {
660#ifdef COIN_DEVELOP
661              printf("primal trouble b\n");
662#endif
663              if (disasterHandler_->typeOfDisaster()) {
664                // We want to abort
665                abortSearch=true;
666                goto disaster;
667              }
668              // try primal with original basis
669              disasterHandler_->setPhase(2);
670              setBasis(basis_,solver);
671              solver->dual();
672              totalIterations += solver->numberIterations();
673            }
674            if(disasterHandler_->inTrouble()) {
675#ifdef COIN_DEVELOP
676              printf("disaster - treat as infeasible\n");
677#endif
678              if (disasterHandler_->typeOfDisaster()) {
679                // We want to abort
680                abortSearch=true;
681                goto disaster;
682              }
683              solver->setProblemStatus(1);
684            }
685          }
686          // reset
687          solver->setDisasterHandler(NULL);
688        }
689        lastAlgorithm_=1; // primal
690        // check if clp thought it was in a loop
691        if (solver->status()==3&&!solver->hitMaximumIterations()) {
692          // switch algorithm
693          solver->dual(0);
694          totalIterations += solver->numberIterations();
695          lastAlgorithm_=2; // dual
696        }
697      }
698    }
699    // If scaled feasible but unscaled infeasible take action
700    if (!solver->status()&&cleanupScaling_) {
701      solver->cleanup(cleanupScaling_);
702    }
703    basis_ = getBasis(solver);
704    //basis_.print();
705    const double * rowScale2 = solver->rowScale();
706    solver->setSpecialOptions(saveOptions);
707    if (!rowScale1&&rowScale2) {
708      // need to release memory
709      if (!solver->savedRowScale_) {
710        solver->setRowScale(NULL);
711        solver->setColumnScale(NULL);
712      } else {
713        solver->rowScale_=NULL;
714        solver->columnScale_=NULL;
715      }
716    }
717  } else {
718    // User doing nothing and all slack basis
719    ClpSolve options=solveOptions_;
720    bool yesNo;
721    OsiHintStrength strength;
722    getHintParam(OsiDoInBranchAndCut,yesNo,strength);
723    if (yesNo) {
724      solver->setSpecialOptions(solver->specialOptions()|1024);
725    }
726    solver->initialSolve(options);
727    totalIterations += solver->numberIterations();
728    lastAlgorithm_ = 2; // say dual
729    // If scaled feasible but unscaled infeasible take action
730    if (!solver->status()&&cleanupScaling_) {
731      solver->cleanup(cleanupScaling_);
732    }
733    basis_ = getBasis(solver);
734    //basis_.print();
735  }
736  solver->messageHandler()->setLogLevel(saveMessageLevel);
737 disaster:
738  if (deleteSolver) {
739    solver->returnModel(*modelPtr_);
740    delete solver;
741  }
742  if (startFinishOptions) {
743    int save = modelPtr_->logLevel();
744    if (save<2) modelPtr_->setLogLevel(0);
745    modelPtr_->dual(0,startFinishOptions);
746    totalIterations += modelPtr_->numberIterations();
747    modelPtr_->setLogLevel(save);
748  }
749  if (saveSolveType==2) {
750    enableSimplexInterface(doingPrimal);
751  }
752  if (savedObjective) {
753    // fix up
754    modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
755    //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
756    modelPtr_->objective_=savedObjective;
757    if (!modelPtr_->problemStatus_) {
758      CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
759      CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
760      if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
761        CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
762      modelPtr_->computeObjectiveValue();
763    }
764  }
765  modelPtr_->setNumberIterations(totalIterations);
766  handler_->setLogLevel(saveMessageLevel2);
767  if (modelPtr_->problemStatus_==3&&lastAlgorithm_==2)
768    modelPtr_->computeObjectiveValue();
769  // mark so we can pick up objective value quickly
770  modelPtr_->upperIn_=0.0;
771  time1 = CoinCpuTime()-time1;
772  totalTime += time1;
773  assert (!modelPtr_->disasterHandler());
774  if (lastAlgorithm_<1||lastAlgorithm_>2)
775    lastAlgorithm_=1;
776  if (abortSearch) {
777    lastAlgorithm_=-911;
778    modelPtr_->setProblemStatus(4);
779  }
780  modelPtr_->whatsChanged_ |= 0x30000;
781#if 0
782  // delete scaled matrix and rowcopy for safety
783  delete modelPtr_->scaledMatrix_;
784  modelPtr_->scaledMatrix_=NULL;
785  delete modelPtr_->rowCopy_;
786  modelPtr_->rowCopy_=NULL;
787#endif
788  //std::cout<<time1<<" seconds - total "<<totalTime<<std::endl;
789  delete pinfo;
790}
791//-----------------------------------------------------------------------------
792void OsiClpSolverInterface::resolve()
793{
794#ifdef COIN_DEVELOP
795  {
796    int i;
797    int n = getNumCols();
798    const double *lower = getColLower() ;
799    const double *upper = getColUpper() ;
800    for (i=0;i<n;i++) {
801      assert (lower[i]<1.0e12);
802      assert (upper[i]>-1.0e12);
803    }
804    n = getNumRows();
805    lower = getRowLower() ;
806    upper = getRowUpper() ;
807    for (i=0;i<n;i++) {
808      assert (lower[i]<1.0e12);
809      assert (upper[i]>-1.0e12);
810    }
811  }
812#endif
813  if ((stuff_.solverOptions_&65536)!=0) {
814    modelPtr_->fastDual2(&stuff_);
815    return;
816  }
817  if ((specialOptions_&2097152)!=0||(specialOptions_&4194304)!=0) {
818    bool takeHint;
819    OsiHintStrength strength;
820    int algorithm = 0;
821    getHintParam(OsiDoDualInResolve,takeHint,strength);
822    if (strength!=OsiHintIgnore)
823      algorithm = takeHint ? -1 : 1;
824    if (algorithm>0||(specialOptions_&4194304)!=0) {
825      // Gub
826      resolveGub((9*modelPtr_->numberRows())/10);
827      return;
828    }
829  }
830  //void pclp(char *);
831  //pclp("res");
832  bool takeHint;
833  OsiHintStrength strength;
834  bool gotHint = (getHintParam(OsiDoInBranchAndCut,takeHint,strength));
835  assert (gotHint);
836  // mark so we can pick up objective value quickly
837  modelPtr_->upperIn_=0.0;
838  if ((specialOptions_&4096)!=0) {
839    // Quick check to see if optimal
840    modelPtr_->checkSolutionInternal();
841    if (modelPtr_->problemStatus()==0) {
842      modelPtr_->setNumberIterations(0);
843      return;
844    }
845  }
846  int totalIterations=0;
847  bool abortSearch=false;
848  ClpObjective * savedObjective=NULL;
849  double savedDualLimit=modelPtr_->dblParam_[ClpDualObjectiveLimit];
850  if (fakeObjective_) {
851    modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~128));
852    // See if all with costs fixed
853    int numberColumns = modelPtr_->numberColumns_;
854    const double * obj = modelPtr_->objective();
855    const double * lower = modelPtr_->columnLower();
856    const double * upper = modelPtr_->columnUpper();
857    int i;
858    for (i=0;i<numberColumns;i++) {
859      double objValue = obj[i];
860      if (objValue) {
861        if (lower[i]!=upper[i])
862          break;
863      }
864    }
865    if (i==numberColumns) {
866      if ((specialOptions_&524288)==0) {
867        // Set fake
868        savedObjective=modelPtr_->objective_;
869        modelPtr_->objective_=fakeObjective_;
870        modelPtr_->dblParam_[ClpDualObjectiveLimit]=COIN_DBL_MAX;
871      } else {
872        modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()|128);
873      }
874    }
875  }
876  // If using Clp initialSolve and primal - just do here
877  gotHint = (getHintParam(OsiDoDualInResolve,takeHint,strength));
878  assert (gotHint);
879  if (strength!=OsiHintIgnore&&!takeHint&&solveOptions_.getSpecialOption(6)) {
880    ClpSolve options=solveOptions_;
881    // presolve
882    getHintParam(OsiDoPresolveInResolve,takeHint,strength);
883    if (strength!=OsiHintIgnore&&!takeHint)
884      options.setPresolveType(ClpSolve::presolveOff);
885    int saveOptions = modelPtr_->specialOptions();
886    getHintParam(OsiDoInBranchAndCut,takeHint,strength);
887    if (takeHint) {
888      modelPtr_->setSpecialOptions(modelPtr_->specialOptions()|1024);
889    }
890    setBasis(basis_,modelPtr_);
891    modelPtr_->initialSolve(options);
892    lastAlgorithm_ = 1; // say primal
893    // If scaled feasible but unscaled infeasible take action
894    if (!modelPtr_->status()&&cleanupScaling_) {
895      modelPtr_->cleanup(cleanupScaling_);
896    }
897    modelPtr_->setSpecialOptions(saveOptions); // restore
898    basis_ = getBasis(modelPtr_);
899  }
900  int saveSolveType=modelPtr_->solveType();
901  bool doingPrimal = modelPtr_->algorithm()>0;
902  if (saveSolveType==2) {
903    disableSimplexInterface();
904  }
905  int saveOptions = modelPtr_->specialOptions();
906  int startFinishOptions=0;
907  if (specialOptions_!=0x80000000) {
908    if((specialOptions_&1)==0) {
909      startFinishOptions=0;
910      modelPtr_->setSpecialOptions(saveOptions|(64|1024|32768));
911    } else {
912      startFinishOptions=1+4;
913      if ((specialOptions_&8)!=0)
914        startFinishOptions +=2; // allow re-use of factorization
915      if((specialOptions_&4)==0||!takeHint) 
916        modelPtr_->setSpecialOptions(saveOptions|(64|128|512|1024|4096|32768));
917      else
918        modelPtr_->setSpecialOptions(saveOptions|(64|128|512|1024|2048|4096|32768));
919    }
920  } else {
921    modelPtr_->setSpecialOptions(saveOptions|64|32768);
922  }
923  //printf("options %d size %d\n",modelPtr_->specialOptions(),modelPtr_->numberColumns());
924  //modelPtr_->setSolveType(1);
925  // Set message handler to have same levels etc
926  int saveMessageLevel=modelPtr_->logLevel();
927  int messageLevel=messageHandler()->logLevel();
928  bool oldDefault;
929  CoinMessageHandler * saveHandler = NULL;
930  if (!defaultHandler_)
931    saveHandler = modelPtr_->pushMessageHandler(handler_,oldDefault);
932  //printf("basis before dual\n");
933  //basis_.print();
934  setBasis(basis_,modelPtr_);
935#ifdef SAVE_MODEL
936  resolveTry++;
937#if SAVE_MODEL > 1
938  if (resolveTry>=loResolveTry&&
939      resolveTry<=hiResolveTry) {
940    char fileName[20];
941    sprintf(fileName,"save%d.mod",resolveTry);
942    modelPtr_->saveModel(fileName);
943  }
944#endif
945#endif
946  // set reasonable defaults
947  // Switch off printing if asked to
948  gotHint = (getHintParam(OsiDoReducePrint,takeHint,strength));
949  assert (gotHint);
950  if (strength!=OsiHintIgnore&&takeHint) {
951    if (messageLevel>0)
952      messageLevel--;
953  }
954  if (messageLevel<modelPtr_->messageHandler()->logLevel())
955    modelPtr_->messageHandler()->setLogLevel(messageLevel);
956  // See if user set factorization frequency
957  int userFactorizationFrequency = modelPtr_->factorization()->maximumPivots();
958  // scaling
959  if (modelPtr_->solveType()==1) {
960    gotHint = (getHintParam(OsiDoScale,takeHint,strength));
961    assert (gotHint);
962    if (strength==OsiHintIgnore||takeHint) {
963      if (!modelPtr_->scalingFlag())
964        modelPtr_->scaling(3);
965    } else {
966      modelPtr_->scaling(0);
967    }
968  } else {
969    modelPtr_->scaling(0);
970  }
971  // sort out hints;
972  // algorithm -1 force dual, +1 force primal
973  int algorithm = -1;
974  gotHint = (getHintParam(OsiDoDualInResolve,takeHint,strength));
975  assert (gotHint);
976  if (strength!=OsiHintIgnore)
977    algorithm = takeHint ? -1 : 1;
978  //modelPtr_->saveModel("save.bad");
979  // presolve
980  gotHint = (getHintParam(OsiDoPresolveInResolve,takeHint,strength));
981  assert (gotHint);
982  if (strength!=OsiHintIgnore&&takeHint) {
983#ifdef KEEP_SMALL
984    if (smallModel_) {
985      delete [] spareArrays_;
986      spareArrays_ = NULL;
987      delete smallModel_;
988      smallModel_=NULL;
989    }
990#endif
991    ClpPresolve pinfo;
992    if ((specialOptions_&128)!=0) {
993      specialOptions_ &= ~128;
994    }
995    if ((modelPtr_->specialOptions()&1024)!=0) {
996      pinfo.setDoDual(false);
997      pinfo.setDoTripleton(false);
998      pinfo.setDoDupcol(false);
999      pinfo.setDoDuprow(false);
1000      pinfo.setDoSingletonColumn(false);
1001    }
1002    ClpSimplex * model2 = pinfo.presolvedModel(*modelPtr_,1.0e-8);
1003    if (!model2) {
1004      // problem found to be infeasible - whats best?
1005      model2 = modelPtr_;
1006    }
1007    // return number of rows
1008    int * stats = reinterpret_cast<int *> (getApplicationData());
1009    if (stats) {
1010      stats[0]=model2->numberRows();
1011      stats[1]=model2->numberColumns();
1012    }
1013    //printf("rows %d -> %d, columns %d -> %d\n",
1014    //     modelPtr_->numberRows(),model2->numberRows(),
1015    //     modelPtr_->numberColumns(),model2->numberColumns());
1016    // change from 200
1017    if (modelPtr_->factorization()->maximumPivots()==200)
1018      model2->factorization()->maximumPivots(100+model2->numberRows()/50);
1019    else
1020      model2->factorization()->maximumPivots(userFactorizationFrequency);
1021    if (algorithm<0) {
1022      model2->dual();
1023      totalIterations += model2->numberIterations();
1024      // check if clp thought it was in a loop
1025      if (model2->status()==3&&!model2->hitMaximumIterations()) {
1026        // switch algorithm
1027        model2->primal();
1028        totalIterations += model2->numberIterations();
1029      }
1030    } else {
1031      model2->primal(1);
1032      totalIterations += model2->numberIterations();
1033      // check if clp thought it was in a loop
1034      if (model2->status()==3&&!model2->hitMaximumIterations()) {
1035        // switch algorithm
1036        model2->dual();
1037        totalIterations += model2->numberIterations();
1038      }
1039    }
1040    if (model2!=modelPtr_) {
1041      int finalStatus=model2->status();
1042      pinfo.postsolve(true);
1043   
1044      delete model2;
1045      // later try without (1) and check duals before solve
1046      if (finalStatus!=3&&(finalStatus||modelPtr_->status()==-1)) {
1047        modelPtr_->primal(1);
1048        totalIterations += modelPtr_->numberIterations();
1049        lastAlgorithm_=1; // primal
1050        //if (modelPtr_->numberIterations())
1051        //printf("****** iterated %d\n",modelPtr_->numberIterations());
1052      }
1053    }
1054  } else {
1055    //modelPtr_->setLogLevel(63);
1056    //modelPtr_->setDualTolerance(1.0e-7);
1057    if (false&&modelPtr_->scalingFlag_>0&&!modelPtr_->rowScale_&&
1058        !modelPtr_->rowCopy_&&matrixByRow_) {
1059      assert (matrixByRow_->getNumElements()==modelPtr_->clpMatrix()->getNumElements());
1060      modelPtr_->setNewRowCopy(new ClpPackedMatrix(*matrixByRow_));
1061    }
1062    if (algorithm<0) {
1063      //writeMps("try1");
1064      int savePerturbation = modelPtr_->perturbation();
1065      if ((specialOptions_&2)!=0)
1066        modelPtr_->setPerturbation(100);
1067      //modelPtr_->messageHandler()->setLogLevel(1);
1068      //writeMpsNative("bad",NULL,NULL,2,1,1.0);
1069      disasterHandler_->setOsiModel(this);
1070      bool inCbcOrOther = (modelPtr_->specialOptions()&0x03000000)!=0;
1071#if 0
1072      // See how many integers fixed
1073      bool skipCrunch=true;
1074      const char * integerInformation = modelPtr_->integerType_;
1075      if (integerInformation) {
1076        int numberColumns = modelPtr_->numberColumns_;
1077        const double * lower = modelPtr_->columnLower();
1078        const double * upper = modelPtr_->columnUpper();
1079        int target=CoinMax(1,numberColumns/10000);
1080        for (int i=0;i<numberColumns;i++) {
1081          if (integerInformation[i]) {
1082            if (lower[i]==upper[i]) {
1083              target--;
1084              if (!target) {
1085                skipCrunch=false;
1086                break;
1087              }
1088            }
1089          }
1090        }
1091      }
1092#endif
1093      if((specialOptions_&1)==0||(specialOptions_&2048)!=0/*||skipCrunch*/) {
1094        disasterHandler_->setWhereFrom(0); // dual
1095        if (inCbcOrOther)
1096          modelPtr_->setDisasterHandler(disasterHandler_);
1097        bool specialScale;
1098        if ((specialOptions_&131072)!=0&&!modelPtr_->rowScale_) {
1099          modelPtr_->rowScale_ = rowScale_.array();
1100          modelPtr_->columnScale_ = columnScale_.array();
1101          specialScale=true;
1102        } else {
1103          specialScale=false;
1104        }
1105#ifdef KEEP_SMALL
1106        if (smallModel_) {
1107          delete [] spareArrays_;
1108          spareArrays_ = NULL;
1109          delete smallModel_;
1110          smallModel_=NULL;
1111        }
1112#endif
1113#ifdef CBC_STATISTICS
1114        osi_dual++;
1115#endif
1116        modelPtr_->dual(0,startFinishOptions);
1117        totalIterations += modelPtr_->numberIterations();
1118        if (specialScale) {
1119          modelPtr_->rowScale_ = NULL;
1120          modelPtr_->columnScale_ = NULL;
1121        }
1122      } else {
1123#ifdef CBC_STATISTICS
1124        osi_crunch++;
1125#endif
1126        crunch();
1127        totalIterations += modelPtr_->numberIterations();
1128        if (modelPtr_->problemStatus()==4)
1129          goto disaster;
1130        // should have already been fixed if problems
1131        inCbcOrOther=false;
1132      }
1133      if (inCbcOrOther) {
1134        if(disasterHandler_->inTrouble()) {
1135          if (disasterHandler_->typeOfDisaster()) {
1136            // We want to abort
1137            abortSearch=true;
1138            goto disaster;
1139          }
1140          // try just going back in
1141          disasterHandler_->setPhase(1);
1142          modelPtr_->dual();
1143          totalIterations += modelPtr_->numberIterations();
1144          if (disasterHandler_->inTrouble()) {
1145            if (disasterHandler_->typeOfDisaster()) {
1146              // We want to abort
1147              abortSearch=true;
1148              goto disaster;
1149            }
1150            // try primal with original basis
1151            disasterHandler_->setPhase(2);
1152            setBasis(basis_,modelPtr_);
1153            modelPtr_->primal();
1154            totalIterations += modelPtr_->numberIterations();
1155          }
1156          if(disasterHandler_->inTrouble()) {
1157#ifdef COIN_DEVELOP
1158            printf("disaster - treat as infeasible\n");
1159#endif
1160            if (disasterHandler_->typeOfDisaster()) {
1161              // We want to abort
1162              abortSearch=true;
1163              goto disaster;
1164            }
1165            modelPtr_->setProblemStatus(1);
1166          }
1167        }
1168        // reset
1169        modelPtr_->setDisasterHandler(NULL);
1170      }
1171      if (modelPtr_->problemStatus()==4) {
1172        // bad bounds?
1173        modelPtr_->setProblemStatus(1);
1174      }
1175      if (!modelPtr_->problemStatus()&&0) {
1176        int numberColumns = modelPtr_->numberColumns();
1177        const double * columnLower = modelPtr_->columnLower();
1178        const double * columnUpper = modelPtr_->columnUpper();
1179        int nBad=0;
1180        for (int i=0;i<numberColumns;i++) {
1181          if (columnLower[i]==columnUpper[i]&&modelPtr_->getColumnStatus(i)==ClpSimplex::basic) {
1182            nBad++;
1183            modelPtr_->setColumnStatus(i,ClpSimplex::isFixed);
1184          }
1185        }
1186        if (nBad) {
1187          modelPtr_->primal(1);
1188    totalIterations += modelPtr_->numberIterations();
1189          printf("%d fixed basic - %d iterations\n",nBad,modelPtr_->numberIterations());
1190        }
1191      }
1192      assert (modelPtr_->objectiveValue()<1.0e100);
1193      modelPtr_->setPerturbation(savePerturbation);
1194      lastAlgorithm_=2; // dual
1195      // check if clp thought it was in a loop
1196      if (modelPtr_->status()==3&&!modelPtr_->hitMaximumIterations()) {
1197        modelPtr_->setSpecialOptions(saveOptions);
1198        // switch algorithm
1199        //modelPtr_->messageHandler()->setLogLevel(63);
1200        // Allow for catastrophe
1201        int saveMax = modelPtr_->maximumIterations();
1202        int numberIterations = modelPtr_->numberIterations();
1203        int numberRows = modelPtr_->numberRows();
1204        int numberColumns = modelPtr_->numberColumns();
1205        if (modelPtr_->maximumIterations()>100000+numberIterations)
1206          modelPtr_->setMaximumIterations(numberIterations + 1000 + 2*numberRows+numberColumns);
1207        modelPtr_->primal(0,startFinishOptions);
1208        totalIterations += modelPtr_->numberIterations();
1209        modelPtr_->setMaximumIterations(saveMax);
1210        lastAlgorithm_=1; // primal
1211        if (modelPtr_->status()==3&&!modelPtr_->hitMaximumIterations()) {
1212#ifdef COIN_DEVELOP
1213          printf("in trouble - try all slack\n");
1214#endif
1215          CoinWarmStartBasis allSlack;
1216          setBasis(allSlack,modelPtr_);
1217          modelPtr_->dual();
1218          totalIterations += modelPtr_->numberIterations();
1219          if (modelPtr_->status()==3&&!modelPtr_->hitMaximumIterations()) {
1220            if (modelPtr_->numberPrimalInfeasibilities()) {
1221#ifdef COIN_DEVELOP
1222              printf("Real real trouble - treat as infeasible\n");
1223#endif
1224              modelPtr_->setProblemStatus(1);
1225            } else {
1226#ifdef COIN_DEVELOP
1227              printf("Real real trouble - treat as optimal\n");
1228#endif
1229              modelPtr_->setProblemStatus(0);
1230            }
1231          }
1232        }
1233      }
1234      assert (modelPtr_->objectiveValue()<1.0e100);
1235    } else {
1236#ifdef KEEP_SMALL
1237      if (smallModel_) {
1238        delete [] spareArrays_;
1239        spareArrays_ = NULL;
1240        delete smallModel_;
1241        smallModel_=NULL;
1242      }
1243#endif
1244      //printf("doing primal\n");
1245#ifdef CBC_STATISTICS
1246      osi_primal++;
1247#endif
1248      modelPtr_->primal(1,startFinishOptions);
1249      totalIterations += modelPtr_->numberIterations();
1250      lastAlgorithm_=1; // primal
1251      // check if clp thought it was in a loop
1252      if (modelPtr_->status()==3&&!modelPtr_->hitMaximumIterations()) {
1253        // switch algorithm
1254        modelPtr_->dual();
1255        totalIterations += modelPtr_->numberIterations();
1256        lastAlgorithm_=2; // dual
1257      }
1258    }
1259  }
1260  // If scaled feasible but unscaled infeasible take action
1261  //if (!modelPtr_->status()&&cleanupScaling_) {
1262  if (cleanupScaling_) {
1263    modelPtr_->cleanup(cleanupScaling_);
1264  }
1265  basis_ = getBasis(modelPtr_);
1266 disaster:
1267  //printf("basis after dual\n");
1268  //basis_.print();
1269  if (!defaultHandler_)
1270    modelPtr_->popMessageHandler(saveHandler,oldDefault);
1271  modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
1272  if (saveSolveType==2) {
1273    int saveStatus = modelPtr_->problemStatus_;
1274    enableSimplexInterface(doingPrimal);
1275    modelPtr_->problemStatus_=saveStatus;
1276  }
1277#ifdef COIN_DEVELOP_x
1278  extern bool doingDoneBranch;
1279  if (doingDoneBranch) {
1280    if (modelPtr_->numberIterations())
1281      printf("***** done %d iterations after general\n",modelPtr_->numberIterations());
1282    doingDoneBranch=false;
1283  }
1284#endif
1285  modelPtr_->setNumberIterations(totalIterations);
1286  //modelPtr_->setSolveType(saveSolveType);
1287  if (abortSearch) {
1288    lastAlgorithm_=-911;
1289    modelPtr_->setProblemStatus(4);
1290  }
1291  if (savedObjective) {
1292    // fix up
1293    modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
1294    //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
1295    modelPtr_->objective_=savedObjective;
1296    if (!modelPtr_->problemStatus_) {
1297      CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
1298      CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
1299      if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
1300        CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
1301      modelPtr_->computeObjectiveValue();
1302    }
1303  }
1304  modelPtr_->setSpecialOptions(saveOptions); // restore
1305  if (modelPtr_->problemStatus_==3&&lastAlgorithm_==2)
1306    modelPtr_->computeObjectiveValue();
1307  if (lastAlgorithm_<1||lastAlgorithm_>2)
1308    lastAlgorithm_=1;
1309#ifdef SAVE_MODEL
1310  if (resolveTry>=loResolveTry&&
1311      resolveTry<=hiResolveTry) {
1312    printf("resolve %d took %d iterations - algorithm %d\n",resolveTry,modelPtr_->numberIterations(),lastAlgorithm_);
1313  }
1314#endif
1315  // Make sure whatsChanged not out of sync
1316  if (!modelPtr_->columnUpperWork_)
1317    modelPtr_->whatsChanged_ &= ~0xffff;
1318  modelPtr_->whatsChanged_ |= 0x30000;
1319}
1320#include "ClpSimplexOther.hpp"
1321// Resolve an LP relaxation after problem modification (try GUB)
1322void 
1323OsiClpSolverInterface::resolveGub(int needed)
1324{
1325  bool takeHint;
1326  OsiHintStrength strength;
1327  // Switch off printing if asked to
1328  getHintParam(OsiDoReducePrint,takeHint,strength);
1329  int saveMessageLevel=modelPtr_->logLevel();
1330  if (strength!=OsiHintIgnore&&takeHint) {
1331    int messageLevel=messageHandler()->logLevel();
1332    if (messageLevel>0)
1333      modelPtr_->messageHandler()->setLogLevel(messageLevel-1);
1334    else
1335      modelPtr_->messageHandler()->setLogLevel(0);
1336  }
1337  //modelPtr_->messageHandler()->setLogLevel(1);
1338  setBasis(basis_,modelPtr_);
1339  // find gub
1340  int numberRows = modelPtr_->numberRows();
1341  int * which = new int[numberRows];
1342  int numberColumns = modelPtr_->numberColumns();
1343  int * whichC = new int[numberColumns+numberRows];
1344  ClpSimplex * model2 = 
1345    static_cast<ClpSimplexOther *> (modelPtr_)->gubVersion(which,whichC,
1346                                                           needed,100);
1347  if (model2) {
1348    // move in solution
1349    static_cast<ClpSimplexOther *> (model2)->setGubBasis(*modelPtr_,
1350                                                         which,whichC);
1351    model2->setLogLevel(CoinMin(1,model2->logLevel()));
1352    ClpPrimalColumnSteepest steepest(5);
1353    model2->setPrimalColumnPivotAlgorithm(steepest);
1354    //double time1 = CoinCpuTime();
1355    model2->primal();
1356    //printf("Primal took %g seconds\n",CoinCpuTime()-time1);
1357    static_cast<ClpSimplexOther *> (model2)->getGubBasis(*modelPtr_,
1358                                                         which,whichC);
1359    int totalIterations = model2->numberIterations();
1360    delete model2;
1361    //modelPtr_->setLogLevel(63);
1362    modelPtr_->primal(1);
1363    modelPtr_->setNumberIterations(totalIterations+modelPtr_->numberIterations());
1364  } else {
1365    modelPtr_->dual();
1366  }
1367  delete [] which;
1368  delete [] whichC;
1369  basis_ = getBasis(modelPtr_);
1370  modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
1371}
1372// Sort of lexicographic resolve
1373void 
1374OsiClpSolverInterface::lexSolve()
1375{
1376#if 1
1377  abort();
1378#else
1379  ((ClpSimplexPrimal *) modelPtr_)->lexSolve();
1380  printf("itA %d\n",modelPtr_->numberIterations());
1381  modelPtr_->primal();
1382  printf("itB %d\n",modelPtr_->numberIterations());
1383  basis_ = getBasis(modelPtr_);
1384#endif
1385}
1386
1387/* Sets up solver for repeated use by Osi interface.
1388   The normal usage does things like keeping factorization around so can be used.
1389   Will also do things like keep scaling and row copy of matrix if
1390   matrix does not change.
1391   adventure:
1392   0 - safe stuff as above
1393   1 - will take more risks - if it does not work then bug which will be fixed
1394   2 - don't bother doing most extreme termination checks e.g. don't bother
1395       re-factorizing if less than 20 iterations.
1396   3 - Actually safer than 1 (mainly just keeps factorization)
1397
1398   printOut - -1 always skip round common messages instead of doing some work
1399               0 skip if normal defaults
1400               1 leaves
1401  */
1402void 
1403OsiClpSolverInterface::setupForRepeatedUse(int senseOfAdventure, int printOut)
1404{
1405  // First try
1406  switch (senseOfAdventure) {
1407  case 0:
1408    specialOptions_=8;
1409    break;
1410  case 1:
1411    specialOptions_=1+2+8;
1412    break;
1413  case 2:
1414    specialOptions_=1+2+4+8;
1415    break;
1416  case 3:
1417    specialOptions_=1+8;
1418    break;
1419  }
1420  bool stopPrinting=false;
1421  if (printOut<0) {
1422    stopPrinting=true;
1423  } else if (!printOut) {
1424    bool takeHint;
1425    OsiHintStrength strength;
1426    getHintParam(OsiDoReducePrint,takeHint,strength);
1427    int messageLevel=messageHandler()->logLevel();
1428    if (strength!=OsiHintIgnore&&takeHint) 
1429      messageLevel--;
1430    stopPrinting = (messageLevel<=0);
1431  }
1432#ifndef COIN_DEVELOP
1433  if (stopPrinting) {
1434    CoinMessages * messagesPointer = modelPtr_->messagesPointer();
1435    // won't even build messages
1436    messagesPointer->setDetailMessages(100,10000,reinterpret_cast<int *> (NULL));
1437  }
1438#endif
1439}
1440#ifndef NDEBUG
1441// For errors to make sure print to screen
1442// only called in debug mode
1443static void indexError(int index,
1444                        std::string methodName)
1445{
1446  std::cerr<<"Illegal index "<<index<<" in OsiClpSolverInterface::"<<methodName<<std::endl;
1447  throw CoinError("Illegal index",methodName,"OsiClpSolverInterface");
1448}
1449#endif
1450//#############################################################################
1451// Parameter related methods
1452//#############################################################################
1453
1454bool
1455OsiClpSolverInterface::setIntParam(OsiIntParam key, int value)
1456{
1457  return modelPtr_->setIntParam(static_cast<ClpIntParam> (key), value);
1458}
1459
1460//-----------------------------------------------------------------------------
1461
1462bool
1463OsiClpSolverInterface::setDblParam(OsiDblParam key, double value)
1464{
1465  if (key != OsiLastDblParam ) {
1466    if (key==OsiDualObjectiveLimit||key==OsiPrimalObjectiveLimit)
1467      value *= modelPtr_->optimizationDirection();
1468    return modelPtr_->setDblParam(static_cast<ClpDblParam> (key), value);
1469  } else {
1470    return false;
1471  }
1472}
1473
1474//-----------------------------------------------------------------------------
1475
1476bool
1477OsiClpSolverInterface::setStrParam(OsiStrParam key, const std::string & value)
1478{
1479  assert (key!=OsiSolverName);
1480  if (key != OsiLastStrParam ) {
1481    return modelPtr_->setStrParam(static_cast<ClpStrParam> (key), value);
1482  } else {
1483    return false;
1484  }
1485}
1486
1487
1488//-----------------------------------------------------------------------------
1489
1490bool
1491OsiClpSolverInterface::getIntParam(OsiIntParam key, int& value) const 
1492{
1493  return modelPtr_->getIntParam(static_cast<ClpIntParam> (key), value);
1494}
1495
1496//-----------------------------------------------------------------------------
1497
1498bool
1499OsiClpSolverInterface::getDblParam(OsiDblParam key, double& value) const
1500{
1501  if (key != OsiLastDblParam ) {
1502    bool condition =  modelPtr_->getDblParam(static_cast<ClpDblParam> (key), value);
1503    if (key==OsiDualObjectiveLimit||key==OsiPrimalObjectiveLimit)
1504      value *= modelPtr_->optimizationDirection();
1505    return condition;
1506  } else {
1507    return false;
1508  }
1509}
1510
1511//-----------------------------------------------------------------------------
1512
1513bool
1514OsiClpSolverInterface::getStrParam(OsiStrParam key, std::string & value) const
1515{
1516  if ( key==OsiSolverName ) {
1517    value = "clp";
1518    return true;
1519  }
1520  if (key != OsiLastStrParam ) {
1521    return modelPtr_->getStrParam(static_cast<ClpStrParam> (key), value);
1522  } else {
1523    return false;
1524  }
1525}
1526
1527
1528//#############################################################################
1529// Methods returning info on how the solution process terminated
1530//#############################################################################
1531
1532bool OsiClpSolverInterface::isAbandoned() const
1533{
1534  // not sure about -1 (should not happen)
1535  return (modelPtr_->status()==4||modelPtr_->status()==-1||
1536          (modelPtr_->status()==1&&modelPtr_->secondaryStatus()==8));
1537}
1538
1539bool OsiClpSolverInterface::isProvenOptimal() const
1540{
1541
1542  const int stat = modelPtr_->status();
1543  return (stat == 0);
1544}
1545
1546bool OsiClpSolverInterface::isProvenPrimalInfeasible() const
1547{
1548
1549  const int stat = modelPtr_->status();
1550  if (stat != 1)
1551     return false;
1552  return true;
1553}
1554
1555bool OsiClpSolverInterface::isProvenDualInfeasible() const
1556{
1557  const int stat = modelPtr_->status();
1558  return stat == 2;
1559}
1560/*
1561   NOTE - Coding if limit > 1.0e30 says that 1.0e29 is loose bound
1562   so all maximization tests are changed
1563*/
1564bool OsiClpSolverInterface::isPrimalObjectiveLimitReached() const
1565{
1566  double limit = 0.0;
1567  modelPtr_->getDblParam(ClpPrimalObjectiveLimit, limit);
1568  if (fabs(limit) > 1e30) {
1569    // was not ever set
1570    return false;
1571  }
1572   
1573  const double obj = modelPtr_->objectiveValue();
1574  int maxmin = static_cast<int> (modelPtr_->optimizationDirection());
1575
1576  switch (lastAlgorithm_) {
1577   case 0: // no simplex was needed
1578     return maxmin > 0 ? (obj < limit) /*minim*/ : (-obj < limit) /*maxim*/;
1579   case 2: // dual simplex
1580     if (modelPtr_->status() == 0) // optimal
1581        return maxmin > 0 ? (obj < limit) /*minim*/ : (-obj < limit) /*maxim*/;
1582     return false;
1583   case 1: // primal simplex
1584     return maxmin > 0 ? (obj < limit) /*minim*/ : (-obj < limit) /*maxim*/;
1585  }
1586  return false; // fake return
1587}
1588
1589bool OsiClpSolverInterface::isDualObjectiveLimitReached() const
1590{
1591
1592  const int stat = modelPtr_->status();
1593  if (stat == 1)
1594  return true;
1595  double limit = 0.0;
1596  modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
1597  if (fabs(limit) > 1e30) {
1598    // was not ever set
1599    return false;
1600  }
1601   
1602  const double obj = modelPtr_->objectiveValue();
1603  int maxmin = static_cast<int> (modelPtr_->optimizationDirection());
1604
1605  switch (lastAlgorithm_) {
1606   case 0: // no simplex was needed
1607     return maxmin > 0 ? (obj > limit) /*minim*/ : (-obj > limit) /*maxim*/;
1608   case 1: // primal simplex
1609     if (stat == 0) // optimal
1610        return maxmin > 0 ? (obj > limit) /*minim*/ : (-obj > limit) /*maxim*/;
1611     return false;
1612   case 2: // dual simplex
1613     if (stat != 0 && stat != 3)
1614        // over dual limit
1615        return true;
1616     return maxmin > 0 ? (obj > limit) /*minim*/ : (-obj > limit) /*maxim*/;
1617  }
1618  return false; // fake return
1619}
1620
1621bool OsiClpSolverInterface::isIterationLimitReached() const
1622{
1623  const int stat = modelPtr_->status();
1624  return (stat == 3);
1625}
1626
1627//#############################################################################
1628// WarmStart related methods
1629//#############################################################################
1630CoinWarmStart *OsiClpSolverInterface::getEmptyWarmStart () const
1631  { return (static_cast<CoinWarmStart *>(new CoinWarmStartBasis())) ; }
1632
1633CoinWarmStart* OsiClpSolverInterface::getWarmStart() const
1634{
1635
1636  return new CoinWarmStartBasis(basis_);
1637}
1638/* Get warm start information.
1639   Return warm start information for the current state of the solver
1640   interface. If there is no valid warm start information, an empty warm
1641   start object wil be returned.  This does not necessarily create an
1642   object - may just point to one.  must Delete set true if user
1643   should delete returned object.
1644   OsiClp version always returns pointer and false.
1645*/
1646CoinWarmStart* 
1647OsiClpSolverInterface::getPointerToWarmStart(bool & mustDelete) 
1648{
1649  mustDelete = false;
1650  return &basis_;
1651}
1652
1653//-----------------------------------------------------------------------------
1654
1655bool OsiClpSolverInterface::setWarmStart(const CoinWarmStart* warmstart)
1656{
1657  modelPtr_->whatsChanged_ &= 0xffff;
1658  const CoinWarmStartBasis* ws =
1659    dynamic_cast<const CoinWarmStartBasis*>(warmstart);
1660  if (ws) {
1661    basis_ = CoinWarmStartBasis(*ws);
1662    return true;
1663  } else if (!warmstart) {
1664    // create from current basis
1665    basis_ = getBasis(modelPtr_);
1666    return true;
1667  } else {
1668    return false;
1669  }
1670}
1671
1672//#############################################################################
1673// Hotstart related methods (primarily used in strong branching)
1674//#############################################################################
1675void OsiClpSolverInterface::markHotStart()
1676{
1677#ifdef CBC_STATISTICS
1678  osi_hot++;
1679#endif
1680  //printf("HotStart options %x changed %x, small %x spare %x\n",
1681  // specialOptions_,modelPtr_->whatsChanged_,
1682  // smallModel_,spareArrays_);
1683  modelPtr_->setProblemStatus(0);
1684  saveData_.perturbation_=0;
1685  saveData_.specialOptions_ = modelPtr_->specialOptions_;
1686  modelPtr_->specialOptions_ |= 0x1000000;
1687  modelPtr_->specialOptions_ = saveData_.specialOptions_; 
1688  ClpObjective * savedObjective=NULL;
1689  double savedDualLimit=modelPtr_->dblParam_[ClpDualObjectiveLimit];
1690  if (fakeObjective_) {
1691    modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~128));
1692    // See if all with costs fixed
1693    int numberColumns = modelPtr_->numberColumns_;
1694    const double * obj = modelPtr_->objective();
1695    const double * lower = modelPtr_->columnLower();
1696    const double * upper = modelPtr_->columnUpper();
1697    int i;
1698    for (i=0;i<numberColumns;i++) {
1699      double objValue = obj[i];
1700      if (objValue) {
1701        if (lower[i]!=upper[i])
1702          break;
1703      }
1704    }
1705    if (i==numberColumns) {
1706      if ((specialOptions_&524288)==0) {
1707        // Set fake
1708        savedObjective=modelPtr_->objective_;
1709        modelPtr_->objective_=fakeObjective_;
1710        modelPtr_->dblParam_[ClpDualObjectiveLimit]=COIN_DBL_MAX;
1711        saveData_.perturbation_=1;
1712      } else {
1713        modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()|128);
1714      }
1715    }
1716  }
1717#define CLEAN_HOT_START
1718#ifdef CLEAN_HOT_START
1719  if ((specialOptions_&65536)!=0) {
1720    //specialOptions_ |= 65536;
1721    saveData_.scalingFlag_=modelPtr_->logLevel();
1722    if (modelPtr_->logLevel()<2)
1723      modelPtr_->setLogLevel(0);
1724    assert ((specialOptions_&128)==0);
1725    // space for save arrays
1726    int numberColumns = modelPtr_->numberColumns();
1727    int numberRows = modelPtr_->numberRows();
1728    // Get space for strong branching
1729    int size = static_cast<int>((1+4*(numberRows+numberColumns))*sizeof(double));
1730    // and for save of original column bounds
1731    size += static_cast<int>(2*numberColumns*sizeof(double));
1732    size += static_cast<int>((1+4*numberRows+2*numberColumns)*sizeof(int));
1733    size += numberRows+numberColumns;
1734    assert (spareArrays_==NULL);
1735    spareArrays_ = new char[size];
1736    //memset(spareArrays_,0x20,size);
1737    // Setup for strong branching
1738    assert (factorization_==NULL);
1739    if ((specialOptions_&131072)!=0) {
1740      assert (lastNumberRows_>=0);
1741      if (modelPtr_->rowScale_!=rowScale_.array()) {
1742        assert(modelPtr_->columnScale_!=columnScale_.array());
1743        delete [] modelPtr_->rowScale_;
1744        modelPtr_->rowScale_=NULL;
1745        delete [] modelPtr_->columnScale_;
1746        modelPtr_->columnScale_=NULL;
1747        if (lastNumberRows_==modelPtr_->numberRows()) {
1748          // use scaling
1749          modelPtr_->rowScale_ = rowScale_.array();
1750          modelPtr_->columnScale_ = columnScale_.array();
1751        } else {
1752          specialOptions_ &= ~131072;
1753        }
1754      }
1755      lastNumberRows_ = -1 -lastNumberRows_;
1756    }
1757    factorization_ = static_cast<ClpSimplexDual *>(modelPtr_)->setupForStrongBranching(spareArrays_,numberRows,
1758                                                                                       numberColumns,true);
1759    double * arrayD = reinterpret_cast<double *> (spareArrays_);
1760    arrayD[0]=modelPtr_->objectiveValue()* modelPtr_->optimizationDirection();
1761    double * saveSolution = arrayD+1;
1762    double * saveLower = saveSolution + (numberRows+numberColumns);
1763    double * saveUpper = saveLower + (numberRows+numberColumns);
1764    double * saveObjective = saveUpper + (numberRows+numberColumns);
1765    double * saveLowerOriginal = saveObjective + (numberRows+numberColumns);
1766    double * saveUpperOriginal = saveLowerOriginal + numberColumns;
1767    CoinMemcpyN( modelPtr_->columnLower(),numberColumns, saveLowerOriginal);
1768    CoinMemcpyN( modelPtr_->columnUpper(),numberColumns, saveUpperOriginal);
1769#if 0
1770    if (whichRange_&&whichRange_[0]) {
1771      // get ranging information
1772      int numberToDo = whichRange_[0];
1773      int * which = new int [numberToDo];
1774      // Convert column numbers
1775      int * backColumn = whichColumn+numberColumns;
1776      for (int i=0;i<numberToDo;i++) {
1777        int iColumn = whichRange_[i+1];
1778        which[i]=backColumn[iColumn];
1779      }
1780      double * downRange=new double [numberToDo];
1781      double * upRange=new double [numberToDo];
1782      int * whichDown = new int [numberToDo];
1783      int * whichUp = new int [numberToDo];
1784      modelPtr_->gutsOfSolution(NULL,NULL,false);
1785      // Tell code we can increase costs in some cases
1786      modelPtr_->setCurrentDualTolerance(0.0);
1787      ((ClpSimplexOther *) modelPtr_)->dualRanging(numberToDo,which,
1788                                                   upRange, whichUp, downRange, whichDown);
1789      delete [] whichDown;
1790      delete [] whichUp;
1791      delete [] which;
1792      rowActivity_=upRange;
1793      columnActivity_=downRange;
1794    }
1795#endif
1796    if (savedObjective) {
1797      // fix up
1798      modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
1799      //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
1800      modelPtr_->objective_=savedObjective;
1801      if (!modelPtr_->problemStatus_) {
1802        CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
1803        CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
1804        if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
1805          CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
1806        modelPtr_->computeObjectiveValue();
1807      }
1808    }
1809    return;
1810  }
1811#endif
1812  if ((specialOptions_&8192)==0&&false) { // ||(specialOptions_&65536)!=0) {
1813    delete ws_;
1814    ws_ = dynamic_cast<CoinWarmStartBasis*>(getWarmStart());
1815    int numberRows = modelPtr_->numberRows();
1816    rowActivity_= new double[numberRows];
1817    CoinMemcpyN(modelPtr_->primalRowSolution(),numberRows,rowActivity_);
1818    int numberColumns = modelPtr_->numberColumns();
1819    columnActivity_= new double[numberColumns];
1820    CoinMemcpyN(modelPtr_->primalColumnSolution(),numberColumns,columnActivity_);
1821  } else {
1822#if 0
1823    int saveLevel = modelPtr_->logLevel();
1824    modelPtr_->setLogLevel(0);
1825    //modelPtr_->dual();
1826    OsiClpSolverInterface::resolve();
1827    if (modelPtr_->numberIterations()>0)
1828      printf("**** iterated large %d\n",modelPtr_->numberIterations());
1829    //else
1830    //printf("no iterations\n");
1831    modelPtr_->setLogLevel(saveLevel);
1832#endif
1833    // called from CbcNode
1834    int numberColumns = modelPtr_->numberColumns();
1835    int numberRows = modelPtr_->numberRows();
1836    // Get space for crunch and strong branching (too much)
1837    int size = static_cast<int>((1+4*(numberRows+numberColumns))*sizeof(double));
1838    // and for save of original column bounds
1839    size += static_cast<int>(2*numberColumns*sizeof(double));
1840    size += static_cast<int>((1+4*numberRows+2*numberColumns)*sizeof(int));
1841    size += numberRows+numberColumns;
1842#ifdef KEEP_SMALL
1843    if(smallModel_&&(modelPtr_->whatsChanged_&0x30000)!=0x30000) {
1844      //printf("Bounds changed ? %x\n",modelPtr_->whatsChanged_);
1845      delete smallModel_;
1846      smallModel_=NULL;
1847    }
1848    if (!smallModel_) {
1849      delete [] spareArrays_;
1850      spareArrays_ = NULL;
1851    }
1852#endif
1853    if (spareArrays_==NULL) {
1854      //if (smallModel_)
1855      //printf("small model %x\n",smallModel_);
1856      delete smallModel_;
1857      smallModel_=NULL;
1858      spareArrays_ = new char[size];
1859      //memset(spareArrays_,0x20,size);
1860    } else {
1861      double * arrayD = reinterpret_cast<double *> (spareArrays_);
1862      double * saveSolution = arrayD+1;
1863      double * saveLower = saveSolution + (numberRows+numberColumns);
1864      double * saveUpper = saveLower + (numberRows+numberColumns);
1865      double * saveObjective = saveUpper + (numberRows+numberColumns);
1866      double * saveLowerOriginal = saveObjective + (numberRows+numberColumns);
1867      double * saveUpperOriginal = saveLowerOriginal + numberColumns;
1868      double * lowerOriginal = modelPtr_->columnLower();
1869      double * upperOriginal = modelPtr_->columnUpper();
1870      arrayD = saveUpperOriginal + numberColumns;
1871      int * savePivot = reinterpret_cast<int *> (arrayD);
1872      int * whichRow = savePivot+numberRows;
1873      int * whichColumn = whichRow + 3*numberRows;
1874      int nSame=0;
1875      int nSub=0;
1876      for (int i=0;i<numberColumns;i++) {
1877        double lo = lowerOriginal[i];
1878        //char * xx = (char *) (saveLowerOriginal+i);
1879        //assert (xx[0]!=0x20||xx[1]!=0x20);
1880        //xx = (char *) (saveUpperOriginal+i);
1881        //assert (xx[0]!=0x20||xx[1]!=0x20);
1882        double loOld = saveLowerOriginal[i];
1883        //assert (!loOld||fabs(loOld)>1.0e-30);
1884        double up = upperOriginal[i];
1885        double upOld = saveUpperOriginal[i];
1886        if (lo>=loOld&&up<=upOld) {
1887          if (lo==loOld&&up==upOld) {
1888            nSame++;
1889          } else {
1890            nSub++;
1891            //if (!isInteger(i))
1892            //nSub+=10;
1893          }
1894        }
1895      }
1896      //printf("Mark Hot %d bounds same, %d interior, %d bad\n",
1897      //     nSame,nSub,numberColumns-nSame-nSub);
1898      if (nSame<numberColumns) {
1899        if (nSame+nSub<numberColumns) {
1900          delete smallModel_;
1901          smallModel_=NULL;
1902        } else {
1903          // we can fix up (but should we if large number fixed?)
1904          assert (smallModel_);
1905          double * lowerSmall = smallModel_->columnLower();
1906          double * upperSmall = smallModel_->columnUpper();
1907          int numberColumns2 = smallModel_->numberColumns();
1908          for (int i=0;i<numberColumns2;i++) {
1909            int iColumn = whichColumn[i];
1910            lowerSmall[i]=lowerOriginal[iColumn];
1911            upperSmall[i]=upperOriginal[iColumn];
1912          }
1913        }
1914      }
1915    }
1916    double * arrayD = reinterpret_cast<double *> (spareArrays_);
1917    arrayD[0]=modelPtr_->objectiveValue()* modelPtr_->optimizationDirection();
1918    double * saveSolution = arrayD+1;
1919    double * saveLower = saveSolution + (numberRows+numberColumns);
1920    double * saveUpper = saveLower + (numberRows+numberColumns);
1921    double * saveObjective = saveUpper + (numberRows+numberColumns);
1922    double * saveLowerOriginal = saveObjective + (numberRows+numberColumns);
1923    double * saveUpperOriginal = saveLowerOriginal + numberColumns;
1924    arrayD = saveUpperOriginal + numberColumns;
1925    int * savePivot = reinterpret_cast<int *> (arrayD);
1926    int * whichRow = savePivot+numberRows;
1927    int * whichColumn = whichRow + 3*numberRows;
1928    int * arrayI = whichColumn + 2*numberColumns;
1929    //unsigned char * saveStatus = (unsigned char *) (arrayI+1);
1930    // Use dual region
1931    double * rhs = modelPtr_->dualRowSolution();
1932    int nBound=0;
1933    ClpSimplex * small;
1934#ifndef KEEP_SMALL
1935    assert (!smallModel_);
1936    small = static_cast<ClpSimplexOther *> (modelPtr_)->crunch(rhs,whichRow,whichColumn,nBound,true);
1937    bool needSolveInSetupHotStart=true;
1938#else
1939    bool needSolveInSetupHotStart=true;
1940    if (!smallModel_) {
1941#ifndef NDEBUG
1942      CoinFillN(whichRow,3*numberRows+2*numberColumns,-1);
1943#endif
1944      small = static_cast<ClpSimplexOther *> (modelPtr_)->crunch(rhs,whichRow,whichColumn,nBound,true);
1945#ifndef NDEBUG
1946      int nCopy = 3*numberRows+2*numberColumns;
1947      for (int i=0;i<nCopy;i++)
1948        assert (whichRow[i]>=-CoinMax(numberRows,numberColumns)&&whichRow[i]<CoinMax(numberRows,numberColumns));
1949#endif
1950      smallModel_=small;
1951      //int hotIts = small->intParam_[ClpMaxNumIterationHotStart];
1952      //if (5*small->factorization_->maximumPivots()>
1953      //  4*hotIts)
1954      //small->factorization_->maximumPivots(hotIts+1);
1955    } else {
1956      assert((modelPtr_->whatsChanged_&0x30000)==0x30000);
1957      //delete [] spareArrays_;
1958      //spareArrays_ = NULL;
1959      assert (spareArrays_);
1960      int nCopy = 3*numberRows+2*numberColumns;
1961      nBound = whichRow[nCopy];
1962#ifndef NDEBUG
1963      for (int i=0;i<nCopy;i++)
1964        assert (whichRow[i]>=-CoinMax(numberRows,numberColumns)&&whichRow[i]<CoinMax(numberRows,numberColumns));
1965#endif
1966      needSolveInSetupHotStart=false;
1967      small = smallModel_;
1968    }
1969#endif
1970    if (small) {
1971      small->specialOptions_ |= 262144;
1972      small->specialOptions_ &= ~65536;
1973    }
1974    if (small&&(specialOptions_&131072)!=0) {
1975      assert (lastNumberRows_>=0);
1976      int numberRows2 = small->numberRows();
1977      int numberColumns2 = small->numberColumns();
1978      double * rowScale2 = new double [2*numberRows2];
1979      const double * rowScale = rowScale_.array();
1980      double * inverseScale2 = rowScale2+numberRows2;
1981      const double * inverseScale = rowScale+modelPtr_->numberRows_;
1982      int i;
1983      for (i=0;i<numberRows2;i++) {
1984        int iRow = whichRow[i];
1985        rowScale2[i]=rowScale[iRow];
1986        inverseScale2[i]=inverseScale[iRow];
1987      }
1988      small->setRowScale(rowScale2);
1989      double * columnScale2 = new double [2*numberColumns2];
1990      const double * columnScale = columnScale_.array();
1991      inverseScale2 = columnScale2+numberColumns2;
1992      inverseScale = columnScale+modelPtr_->numberColumns_;
1993      for (i=0;i<numberColumns2;i++) {
1994        int iColumn = whichColumn[i];
1995        columnScale2[i]=columnScale[iColumn];
1996        inverseScale2[i]=inverseScale[iColumn];
1997      }
1998      small->setColumnScale(columnScale2);
1999    }
2000    if (!small) {
2001      // should never be infeasible .... but
2002      delete [] spareArrays_;
2003      spareArrays_=NULL;
2004      delete ws_;
2005      ws_ = dynamic_cast<CoinWarmStartBasis*>(getWarmStart());
2006      int numberRows = modelPtr_->numberRows();
2007      rowActivity_= new double[numberRows];
2008      CoinMemcpyN(modelPtr_->primalRowSolution(),numberRows,rowActivity_);
2009      int numberColumns = modelPtr_->numberColumns();
2010      columnActivity_= new double[numberColumns];
2011      CoinMemcpyN(modelPtr_->primalColumnSolution(),numberColumns,columnActivity_);
2012      modelPtr_->setProblemStatus(1);
2013      if (savedObjective) {
2014        // fix up
2015        modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2016        //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2017        modelPtr_->objective_=savedObjective;
2018      }
2019      return;
2020    }
2021    int clpOptions = modelPtr_->specialOptions();
2022    clpOptions &= ~65536;
2023    if((specialOptions_&1)==0) {
2024      small->setSpecialOptions(clpOptions|(64|1024));
2025    } else {
2026      if((specialOptions_&4)==0) 
2027        small->setSpecialOptions(clpOptions|(64|128|512|1024|4096));
2028      else
2029        small->setSpecialOptions(clpOptions|(64|128|512|1024|2048|4096));
2030    }
2031    arrayI[0]=nBound;
2032    assert (smallModel_==NULL||small==smallModel_);
2033    if ((specialOptions_&256)!=0||1) {
2034      // only need to do this on second pass in CbcNode
2035      if (modelPtr_->logLevel()<2) small->setLogLevel(0);
2036      small->specialOptions_ |= 262144;
2037      small->moreSpecialOptions_ = modelPtr_->moreSpecialOptions_;
2038#define SETUP_HOT
2039#ifndef SETUP_HOT
2040      small->dual();
2041#else
2042      assert (factorization_==NULL);
2043      //needSolveInSetupHotStart=true;
2044      ClpFactorization * factorization = static_cast<ClpSimplexDual *>(small)->setupForStrongBranching(spareArrays_,
2045                       numberRows,numberColumns,
2046                                                                                                       needSolveInSetupHotStart);
2047#endif
2048      if (small->numberIterations()>0&&small->logLevel()>2)
2049        printf("**** iterated small %d\n",small->numberIterations());
2050      //small->setLogLevel(0);
2051      // Could be infeasible if forced one way (and other way stopped on iterations)
2052      // could also be stopped on iterations
2053      if (small->status()) {
2054#ifndef KEEP_SMALL
2055        if (small!=modelPtr_)
2056          delete small;
2057        //delete smallModel_;
2058        //smallModel_=NULL;
2059        assert (!smallModel_);
2060#else
2061        assert (small==smallModel_);
2062        if (smallModel_!=modelPtr_) {
2063          delete smallModel_;
2064        }
2065        smallModel_=NULL;
2066#endif
2067        delete [] spareArrays_;
2068        spareArrays_=NULL;
2069        delete ws_;
2070        ws_ = dynamic_cast<CoinWarmStartBasis*>(getWarmStart());
2071        int numberRows = modelPtr_->numberRows();
2072        rowActivity_= new double[numberRows];
2073 CoinMemcpyN(modelPtr_->primalRowSolution(),    numberRows,rowActivity_);
2074        int numberColumns = modelPtr_->numberColumns();
2075        columnActivity_= new double[numberColumns];
2076 CoinMemcpyN(modelPtr_->primalColumnSolution(), numberColumns,columnActivity_);
2077        modelPtr_->setProblemStatus(1);
2078        if (savedObjective) {
2079          // fix up
2080          modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2081          //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2082          modelPtr_->objective_=savedObjective;
2083        }
2084        return;
2085      } else {
2086        // update model
2087        static_cast<ClpSimplexOther *> (modelPtr_)->afterCrunch(*small,whichRow,whichColumn,nBound);
2088      } 
2089#ifndef SETUP_HOT
2090      assert (factorization_==NULL);
2091      factorization_ = static_cast<ClpSimplexDual *>(small)->setupForStrongBranching(spareArrays_,numberRows,
2092                             numberColumns,false);
2093#else
2094      assert (factorization!=NULL || small->problemStatus_ );
2095      factorization_ = factorization;
2096#endif
2097    } else {
2098      assert (factorization_==NULL);
2099      factorization_ = static_cast<ClpSimplexDual *>(small)->setupForStrongBranching(spareArrays_,numberRows,
2100                                                                                     numberColumns,false);
2101    }
2102    smallModel_=small;
2103    if (modelPtr_->logLevel()<2) smallModel_->setLogLevel(0);
2104    // Setup for strong branching
2105    int numberColumns2 = smallModel_->numberColumns();
2106    CoinMemcpyN( modelPtr_->columnLower(),numberColumns, saveLowerOriginal);
2107    CoinMemcpyN( modelPtr_->columnUpper(),numberColumns, saveUpperOriginal);
2108    const double * smallLower = smallModel_->columnLower();
2109    const double * smallUpper = smallModel_->columnUpper();
2110    // But modify if bounds changed in small
2111    for (int i=0;i<numberColumns2;i++) {
2112      int iColumn = whichColumn[i];
2113      saveLowerOriginal[iColumn] = CoinMax(saveLowerOriginal[iColumn],
2114                                           smallLower[i]);
2115      saveUpperOriginal[iColumn] = CoinMin(saveUpperOriginal[iColumn],
2116                                           smallUpper[i]);
2117    }
2118    if (whichRange_&&whichRange_[0]) {
2119      // get ranging information
2120      int numberToDo = whichRange_[0];
2121      int * which = new int [numberToDo];
2122      // Convert column numbers
2123      int * backColumn = whichColumn+numberColumns;
2124      for (int i=0;i<numberToDo;i++) {
2125        int iColumn = whichRange_[i+1];
2126        which[i]=backColumn[iColumn];
2127      }
2128      double * downRange=new double [numberToDo];
2129      double * upRange=new double [numberToDo];
2130      int * whichDown = new int [numberToDo];
2131      int * whichUp = new int [numberToDo];
2132      smallModel_->setFactorization(*factorization_);
2133      smallModel_->gutsOfSolution(NULL,NULL,false);
2134      // Tell code we can increase costs in some cases
2135      smallModel_->setCurrentDualTolerance(0.0);
2136      static_cast<ClpSimplexOther *> (smallModel_)->dualRanging(numberToDo,which,
2137                         upRange, whichUp, downRange, whichDown);
2138      delete [] whichDown;
2139      delete [] whichUp;
2140      delete [] which;
2141      rowActivity_=upRange;
2142      columnActivity_=downRange;
2143    }
2144  }
2145    if (savedObjective) {
2146      // fix up
2147      modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2148      //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2149      modelPtr_->objective_=savedObjective;
2150      if (!modelPtr_->problemStatus_) {
2151        CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
2152        CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
2153        if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
2154          CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
2155        modelPtr_->computeObjectiveValue();
2156      }
2157    }
2158}
2159
2160void OsiClpSolverInterface::solveFromHotStart()
2161{
2162#ifdef KEEP_SMALL
2163  if (!spareArrays_) {
2164    assert (!smallModel_);
2165    assert (modelPtr_->problemStatus_==1);
2166    return;
2167  }
2168#endif
2169  ClpObjective * savedObjective=NULL;
2170  double savedDualLimit=modelPtr_->dblParam_[ClpDualObjectiveLimit];
2171  if (saveData_.perturbation_) {
2172    // Set fake
2173    savedObjective=modelPtr_->objective_;
2174    modelPtr_->objective_=fakeObjective_;
2175    modelPtr_->dblParam_[ClpDualObjectiveLimit]=COIN_DBL_MAX;
2176  }
2177  int numberRows = modelPtr_->numberRows();
2178  int numberColumns = modelPtr_->numberColumns();
2179  modelPtr_->getIntParam(ClpMaxNumIteration,itlimOrig_);
2180  int itlim;
2181  modelPtr_->getIntParam(ClpMaxNumIterationHotStart, itlim);
2182  // Is there an extra copy of scaled bounds
2183  int extraCopy = (modelPtr_->maximumRows_>0) ? modelPtr_->maximumRows_+modelPtr_->maximumColumns_ : 0; 
2184#ifdef CLEAN_HOT_START
2185  if ((specialOptions_&65536)!=0) {
2186    double * arrayD = reinterpret_cast<double *> (spareArrays_);
2187    double saveObjectiveValue = arrayD[0];
2188    double * saveSolution = arrayD+1;
2189    int number = numberRows+numberColumns;
2190    CoinMemcpyN(saveSolution,number,modelPtr_->solutionRegion());
2191    double * saveLower = saveSolution + (numberRows+numberColumns);
2192    CoinMemcpyN(saveLower,number,modelPtr_->lowerRegion());
2193    double * saveUpper = saveLower + (numberRows+numberColumns);
2194    CoinMemcpyN(saveUpper,number,modelPtr_->upperRegion());
2195    double * saveObjective = saveUpper + (numberRows+numberColumns);
2196    CoinMemcpyN(saveObjective,number,modelPtr_->costRegion());
2197    double * saveLowerOriginal = saveObjective + (numberRows+numberColumns);
2198    double * saveUpperOriginal = saveLowerOriginal + numberColumns;
2199    arrayD = saveUpperOriginal + numberColumns;
2200    int * savePivot = reinterpret_cast<int *> (arrayD);
2201    CoinMemcpyN(savePivot,numberRows,modelPtr_->pivotVariable());
2202    int * whichRow = savePivot+numberRows;
2203    int * whichColumn = whichRow + 3*numberRows;
2204    int * arrayI = whichColumn + 2*numberColumns;
2205    unsigned char * saveStatus = reinterpret_cast<unsigned char *> (arrayI+1);
2206    CoinMemcpyN(saveStatus,number,modelPtr_->statusArray());
2207    modelPtr_->setFactorization(*factorization_);
2208    double * columnLower = modelPtr_->columnLower();
2209    double * columnUpper = modelPtr_->columnUpper();
2210    // make sure whatsChanged_ has 1 set
2211    modelPtr_->setWhatsChanged(511);
2212    double * lowerInternal = modelPtr_->lowerRegion();
2213    double * upperInternal = modelPtr_->upperRegion();
2214    double rhsScale = modelPtr_->rhsScale();
2215    const double * columnScale = NULL;
2216    if (modelPtr_->scalingFlag()>0) 
2217      columnScale = modelPtr_->columnScale() ;
2218    // and do bounds in case dual needs them
2219    int iColumn;
2220    for (iColumn=0;iColumn<numberColumns;iColumn++) {
2221      if (columnLower[iColumn]>saveLowerOriginal[iColumn]) {
2222        double value = columnLower[iColumn];
2223        value *= rhsScale;
2224        if (columnScale)
2225          value /= columnScale[iColumn];
2226        lowerInternal[iColumn]=value;
2227        if (extraCopy)
2228          lowerInternal[iColumn+extraCopy]=value;
2229      }
2230      if (columnUpper[iColumn]<saveUpperOriginal[iColumn]) {
2231        double value = columnUpper[iColumn];
2232        value *= rhsScale;
2233        if (columnScale)
2234          value /= columnScale[iColumn];
2235        upperInternal[iColumn]=value;
2236        if (extraCopy)
2237          upperInternal[iColumn+extraCopy]=value;
2238      }
2239    }
2240    //#define REDO_STEEPEST_EDGE
2241#ifdef REDO_STEEPEST_EDGE
2242    if (dynamic_cast<ClpDualRowSteepest *>(modelPtr_->dualRowPivot())) {
2243      // new copy (should think about keeping around)
2244      ClpDualRowSteepest steepest;
2245      modelPtr_->setDualRowPivotAlgorithm(steepest);
2246    }     
2247#endif
2248    // Start of fast iterations
2249    bool alwaysFinish= ((specialOptions_&32)==0) ? true : false;
2250    //modelPtr_->setLogLevel(1);
2251    modelPtr_->setIntParam(ClpMaxNumIteration,itlim);
2252    int saveNumberFake = (static_cast<ClpSimplexDual *>(modelPtr_))->numberFake_;
2253    int status = (static_cast<ClpSimplexDual *>(modelPtr_))->fastDual(alwaysFinish);
2254    (static_cast<ClpSimplexDual *>(modelPtr_))->numberFake_ = saveNumberFake;
2255   
2256    int problemStatus = modelPtr_->problemStatus();
2257    double objectiveValue =modelPtr_->objectiveValue() * modelPtr_->optimizationDirection();
2258    CoinAssert (modelPtr_->problemStatus()||modelPtr_->objectiveValue()<1.0e50);
2259    // make sure plausible
2260    double obj = CoinMax(objectiveValue,saveObjectiveValue);
2261    if (problemStatus==10||problemStatus<0) {
2262      // was trying to clean up or something odd
2263      if (problemStatus==10)
2264        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2265      status=1;
2266    }
2267    if (status) {
2268      // not finished - might be optimal
2269      modelPtr_->checkPrimalSolution(modelPtr_->solutionRegion(0),
2270                                     modelPtr_->solutionRegion(1));
2271      //modelPtr_->gutsOfSolution(NULL,NULL,0);
2272      //if (problemStatus==3)
2273      //modelPtr_->computeObjectiveValue();
2274      objectiveValue =modelPtr_->objectiveValue() *
2275        modelPtr_->optimizationDirection();
2276      obj = CoinMax(objectiveValue,saveObjectiveValue);
2277      if (!modelPtr_->numberDualInfeasibilities()) { 
2278        double limit = 0.0;
2279        modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2280        if (modelPtr_->secondaryStatus()==1&&!problemStatus&&obj<limit) {
2281          obj=limit;
2282          problemStatus=3;
2283        }
2284        if (!modelPtr_->numberPrimalInfeasibilities()&&obj<limit) { 
2285          problemStatus=0;
2286        } else if (problemStatus==10) {
2287          problemStatus=3;
2288        } else if (!modelPtr_->numberPrimalInfeasibilities()) {
2289          problemStatus=1; // infeasible
2290        } 
2291      } else {
2292        // can't say much
2293        //if (problemStatus==3)
2294        //modelPtr_->computeObjectiveValue();
2295        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2296        problemStatus=3;
2297      }
2298    } else if (!problemStatus) {
2299      if (modelPtr_->isDualObjectiveLimitReached()) 
2300        problemStatus=1; // infeasible
2301    }
2302    if (status&&!problemStatus) {
2303      problemStatus=3; // can't be sure
2304      lastAlgorithm_=1;
2305    }
2306    if (problemStatus<0)
2307      problemStatus=3;
2308    modelPtr_->setProblemStatus(problemStatus);
2309    modelPtr_->setObjectiveValue(obj*modelPtr_->optimizationDirection());
2310    double * solution = modelPtr_->primalColumnSolution();
2311    const double * solution2 = modelPtr_->solutionRegion();
2312    // could just do changed bounds - also try double size scale so can use * not /
2313    if (!columnScale) {
2314      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2315        solution[iColumn]= solution2[iColumn];
2316      }
2317    } else {
2318      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2319        solution[iColumn]= solution2[iColumn]*columnScale[iColumn];
2320      }
2321    }
2322    CoinMemcpyN(saveLowerOriginal,numberColumns,columnLower);
2323    CoinMemcpyN(saveUpperOriginal,numberColumns,columnUpper);
2324#if 0
2325    // could combine with loop above
2326    if (modelPtr_==modelPtr_)
2327      modelPtr_->computeObjectiveValue();
2328    if (status&&!problemStatus) {
2329      memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2330      modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2331      modelPtr_->checkSolutionInternal();
2332      //modelPtr_->setLogLevel(1);
2333      //modelPtr_->allSlackBasis();
2334      //modelPtr_->primal(1);
2335      //memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2336      //modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2337      //modelPtr_->checkSolutionInternal();
2338      assert (!modelPtr_->problemStatus());
2339    }
2340#endif
2341    // and back bounds
2342    CoinMemcpyN(saveLower,number,modelPtr_->lowerRegion());
2343    CoinMemcpyN(saveUpper,number,modelPtr_->upperRegion());
2344    if (extraCopy) {
2345      CoinMemcpyN(saveLower,number,modelPtr_->lowerRegion()+extraCopy);
2346      CoinMemcpyN(saveUpper,number,modelPtr_->upperRegion()+extraCopy);
2347    }
2348    modelPtr_->setIntParam(ClpMaxNumIteration,itlimOrig_);
2349    if (savedObjective) {
2350      // fix up
2351      modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2352      //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2353      modelPtr_->objective_=savedObjective;
2354      if (!modelPtr_->problemStatus_) {
2355        CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
2356        CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
2357        if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
2358          CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
2359        modelPtr_->computeObjectiveValue();
2360      }
2361    }
2362    return;
2363  }
2364#endif
2365  if (smallModel_==NULL) {
2366    setWarmStart(ws_);
2367    CoinMemcpyN(           rowActivity_,numberRows,modelPtr_->primalRowSolution());
2368    CoinMemcpyN(columnActivity_,numberColumns,modelPtr_->primalColumnSolution());
2369    modelPtr_->setIntParam(ClpMaxNumIteration,CoinMin(itlim,9999));
2370    resolve();
2371  } else {
2372    assert (spareArrays_);
2373    double * arrayD = reinterpret_cast<double *> (spareArrays_);
2374    double saveObjectiveValue = arrayD[0];
2375    double * saveSolution = arrayD+1;
2376    // double check arrays exist (? for nonlinear)
2377    //if (!smallModel_->solutionRegion())
2378    //smallModel_->createRim(63);
2379    int numberRows2 = smallModel_->numberRows();
2380    int numberColumns2 = smallModel_->numberColumns();
2381    int number = numberRows2+numberColumns2;
2382    CoinMemcpyN(saveSolution,number,smallModel_->solutionRegion());
2383    double * saveLower = saveSolution + (numberRows+numberColumns);
2384    CoinMemcpyN(saveLower,number,smallModel_->lowerRegion());
2385    double * saveUpper = saveLower + (numberRows+numberColumns);
2386    CoinMemcpyN(saveUpper,number,smallModel_->upperRegion());
2387    double * saveObjective = saveUpper + (numberRows+numberColumns);
2388    CoinMemcpyN(saveObjective,number,smallModel_->costRegion());
2389    double * saveLowerOriginal = saveObjective + (numberRows+numberColumns);
2390    double * saveUpperOriginal = saveLowerOriginal + numberColumns;
2391    arrayD = saveUpperOriginal + numberColumns;
2392    int * savePivot = reinterpret_cast<int *> (arrayD);
2393    CoinMemcpyN(savePivot,numberRows2,smallModel_->pivotVariable());
2394    int * whichRow = savePivot+numberRows;
2395    int * whichColumn = whichRow + 3*numberRows;
2396    int * arrayI = whichColumn + 2*numberColumns;
2397    unsigned char * saveStatus = reinterpret_cast<unsigned char *> (arrayI+1);
2398    CoinMemcpyN(saveStatus,number,smallModel_->statusArray());
2399    /* If factorization_ NULL then infeasible
2400       not really sure how could have slipped through.
2401       But this can't make situation worse */
2402    if (factorization_) 
2403      smallModel_->setFactorization(*factorization_);
2404    //int * backColumn = whichColumn+numberColumns;
2405    const double * lowerBig = modelPtr_->columnLower();
2406    const double * upperBig = modelPtr_->columnUpper();
2407    // make sure whatsChanged_ has 1 set
2408    smallModel_->setWhatsChanged(511);
2409    double * lowerSmall = smallModel_->lowerRegion();
2410    double * upperSmall = smallModel_->upperRegion();
2411    double * lowerSmallReal = smallModel_->columnLower();
2412    double * upperSmallReal = smallModel_->columnUpper();
2413    int i;
2414    double rhsScale = smallModel_->rhsScale();
2415    const double * columnScale = NULL;
2416    if (smallModel_->scalingFlag()>0) {
2417      columnScale = smallModel_->columnScale();
2418    }
2419    // and do bounds in case dual needs them
2420    // may be infeasible
2421    for (i=0;i<numberColumns2;i++) {
2422      int iColumn = whichColumn[i];
2423      if (lowerBig[iColumn]>saveLowerOriginal[iColumn]) {
2424        double value = lowerBig[iColumn];
2425        lowerSmallReal[i]=value;
2426        value *= rhsScale;
2427        if (columnScale)
2428          value /= columnScale[i];
2429        lowerSmall[i]=value;
2430      }
2431      if (upperBig[iColumn]<saveUpperOriginal[iColumn]) {
2432        double value = upperBig[iColumn];
2433        upperSmallReal[i]=value;
2434        value *= rhsScale;
2435        if (columnScale)
2436          value /= columnScale[i];
2437        upperSmall[i]=value;
2438      }
2439      if (upperSmall[i]<lowerSmall[i]-1.0e-8)
2440        break;
2441    }
2442    /* If factorization_ NULL then infeasible
2443       not really sure how could have slipped through.
2444       But this can't make situation worse */
2445    bool infeasible = (i<numberColumns2||!factorization_);
2446    // Start of fast iterations
2447    bool alwaysFinish= ((specialOptions_&32)==0) ? true : false;
2448    //smallModel_->setLogLevel(1);
2449    smallModel_->setIntParam(ClpMaxNumIteration,itlim);
2450    int saveNumberFake = (static_cast<ClpSimplexDual *>(smallModel_))->numberFake_;
2451    int status;
2452    if (!infeasible) {
2453#ifdef REDO_STEEPEST_EDGE
2454      if (dynamic_cast<ClpDualRowSteepest *>(smallModel_->dualRowPivot())) {
2455        // new copy
2456        ClpDualRowSteepest steepest;
2457        smallModel_->setDualRowPivotAlgorithm(steepest);
2458      }     
2459#endif
2460      status = static_cast<ClpSimplexDual *>(smallModel_)->fastDual(alwaysFinish);
2461    } else {
2462      status=0;
2463      smallModel_->setProblemStatus(1);
2464    }
2465    (static_cast<ClpSimplexDual *>(smallModel_))->numberFake_ = saveNumberFake;
2466    if (smallModel_->numberIterations()==-98) {
2467      printf("rrrrrrrrrrrr\n");
2468      smallModel_->checkPrimalSolution(smallModel_->solutionRegion(0),
2469                                     smallModel_->solutionRegion(1));
2470      //smallModel_->gutsOfSolution(NULL,NULL,0);
2471      //if (problemStatus==3)
2472      //smallModel_->computeObjectiveValue();
2473      printf("robj %g\n",smallModel_->objectiveValue() *
2474             modelPtr_->optimizationDirection());
2475      writeMps("rr.mps");
2476      smallModel_->writeMps("rr_small.mps");
2477      ClpSimplex temp = *smallModel_;
2478      printf("small\n");
2479      temp.setLogLevel(63);
2480      temp.dual();
2481      double limit = 0.0;
2482      modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2483      if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2484        printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2485               temp.objectiveValue(),
2486               smallModel_->objectiveOffset(),temp.objectiveOffset());
2487      }
2488      printf("big\n");
2489      temp = *modelPtr_;
2490      temp.dual();
2491      if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2492        printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2493               temp.objectiveValue(),
2494               smallModel_->objectiveOffset(),temp.objectiveOffset());
2495      }
2496    }
2497    int problemStatus = smallModel_->problemStatus();
2498    double objectiveValue =smallModel_->objectiveValue() * modelPtr_->optimizationDirection();
2499    CoinAssert (smallModel_->problemStatus()||smallModel_->objectiveValue()<1.0e50);
2500    // make sure plausible
2501    double obj = CoinMax(objectiveValue,saveObjectiveValue);
2502    if (problemStatus==10||problemStatus<0) {
2503      // was trying to clean up or something odd
2504      if (problemStatus==10)
2505        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2506      status=1;
2507    }
2508    if (status) {
2509      // not finished - might be optimal
2510      smallModel_->checkPrimalSolution(smallModel_->solutionRegion(0),
2511                                     smallModel_->solutionRegion(1));
2512      //smallModel_->gutsOfSolution(NULL,NULL,0);
2513      //if (problemStatus==3)
2514      //smallModel_->computeObjectiveValue();
2515      objectiveValue =smallModel_->objectiveValue() *
2516        modelPtr_->optimizationDirection();
2517      if (problemStatus!=10)
2518        obj = CoinMax(objectiveValue,saveObjectiveValue);
2519      if (!smallModel_->numberDualInfeasibilities()) { 
2520        double limit = 0.0;
2521        modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2522        if (smallModel_->secondaryStatus()==1&&!problemStatus&&obj<limit) {
2523#if 0
2524          // switch off
2525          ClpSimplex temp = *smallModel_;
2526          temp.dual();
2527          if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2528            printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2529                   temp.objectiveValue(),
2530                   smallModel_->objectiveOffset(),temp.objectiveOffset());
2531          }
2532          lastAlgorithm_=1;
2533          obj=limit;
2534          problemStatus=10;
2535#else
2536          obj=limit;
2537          problemStatus=3;
2538#endif
2539        }
2540        if (!smallModel_->numberPrimalInfeasibilities()&&obj<limit) { 
2541          problemStatus=0;
2542#if 0
2543          ClpSimplex temp = *smallModel_;
2544          temp.dual();
2545          if (temp.numberIterations())
2546            printf("temp iterated\n");
2547          assert (temp.problemStatus()==0&&temp.objectiveValue()<limit);
2548#endif
2549        } else if (problemStatus==10) {
2550          problemStatus=3;
2551        } else if (!smallModel_->numberPrimalInfeasibilities()) {
2552          problemStatus=1; // infeasible
2553        } 
2554      } else {
2555        // can't say much
2556        //if (problemStatus==3)
2557        //smallModel_->computeObjectiveValue();
2558        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2559        problemStatus=3;
2560      }
2561    } else if (!problemStatus) {
2562      if (smallModel_->isDualObjectiveLimitReached()) 
2563        problemStatus=1; // infeasible
2564    }
2565    if (status&&!problemStatus) {
2566      problemStatus=3; // can't be sure
2567      lastAlgorithm_=1;
2568    }
2569    if (problemStatus<0)
2570      problemStatus=3;
2571    modelPtr_->setProblemStatus(problemStatus);
2572    modelPtr_->setObjectiveValue(obj*modelPtr_->optimizationDirection());
2573    modelPtr_->setSumDualInfeasibilities(smallModel_->sumDualInfeasibilities());
2574    modelPtr_->setNumberDualInfeasibilities(smallModel_->numberDualInfeasibilities());
2575    modelPtr_->setSumPrimalInfeasibilities(smallModel_->sumPrimalInfeasibilities());
2576    modelPtr_->setNumberPrimalInfeasibilities(smallModel_->numberPrimalInfeasibilities());
2577    double * solution = modelPtr_->primalColumnSolution();
2578    const double * solution2 = smallModel_->solutionRegion();
2579    if (!columnScale) {
2580      for (i=0;i<numberColumns2;i++) {
2581        int iColumn = whichColumn[i];
2582        solution[iColumn]= solution2[i];
2583        lowerSmallReal[i]=saveLowerOriginal[iColumn];
2584        upperSmallReal[i]=saveUpperOriginal[iColumn];
2585      }
2586    } else {
2587      for (i=0;i<numberColumns2;i++) {
2588        int iColumn = whichColumn[i];
2589        solution[iColumn]= solution2[i]*columnScale[i];
2590        lowerSmallReal[i]=saveLowerOriginal[iColumn];
2591        upperSmallReal[i]=saveUpperOriginal[iColumn];
2592      }
2593    }
2594    // could combine with loop above
2595    if (modelPtr_==smallModel_)
2596      modelPtr_->computeObjectiveValue();
2597#if 1
2598    if (status&&!problemStatus) {
2599      memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2600      modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2601      modelPtr_->checkSolutionInternal();
2602      //modelPtr_->setLogLevel(1);
2603      //modelPtr_->allSlackBasis();
2604      //modelPtr_->primal(1);
2605      //memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2606      //modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2607      //modelPtr_->checkSolutionInternal();
2608      assert (!modelPtr_->problemStatus());
2609    }
2610#endif
2611    modelPtr_->setNumberIterations(smallModel_->numberIterations());
2612    // and back bounds
2613    CoinMemcpyN(saveLower,number,smallModel_->lowerRegion());
2614    CoinMemcpyN(saveUpper,number,smallModel_->upperRegion());
2615  }
2616  if (savedObjective) {
2617    // fix up
2618    modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2619    //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2620    modelPtr_->objective_=savedObjective;
2621    if (!modelPtr_->problemStatus_) {
2622      CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
2623      CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
2624      if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
2625        CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
2626      modelPtr_->computeObjectiveValue();
2627    }
2628  }
2629  modelPtr_->setIntParam(ClpMaxNumIteration,itlimOrig_);
2630}
2631
2632void OsiClpSolverInterface::unmarkHotStart()
2633{
2634#ifdef CLEAN_HOT_START
2635  if ((specialOptions_&65536)!=0) {
2636    modelPtr_->setLogLevel(saveData_.scalingFlag_);
2637    modelPtr_->deleteRim(0);
2638    if (lastNumberRows_<0) {
2639      specialOptions_ |= 131072;
2640      lastNumberRows_ = -1 -lastNumberRows_;
2641      if (modelPtr_->rowScale_) {
2642        if (modelPtr_->rowScale_!=rowScale_.array()) {
2643          delete [] modelPtr_->rowScale_;
2644          delete [] modelPtr_->columnScale_;
2645        }
2646        modelPtr_->rowScale_=NULL;
2647        modelPtr_->columnScale_=NULL;
2648      }
2649    }
2650    delete factorization_;
2651    delete [] spareArrays_;
2652    smallModel_=NULL;
2653    spareArrays_=NULL;
2654    factorization_=NULL;
2655    delete [] rowActivity_;
2656    delete [] columnActivity_;
2657    rowActivity_=NULL;
2658    columnActivity_=NULL;
2659    return;
2660  }
2661#endif
2662  if (smallModel_==NULL) {
2663    setWarmStart(ws_);
2664    int numberRows = modelPtr_->numberRows();
2665    int numberColumns = modelPtr_->numberColumns();
2666    CoinMemcpyN(           rowActivity_,numberRows,modelPtr_->primalRowSolution());
2667    CoinMemcpyN(columnActivity_,numberColumns,modelPtr_->primalColumnSolution());
2668    delete ws_;
2669    ws_ = NULL;
2670  } else {
2671#ifndef KEEP_SMALL
2672    if (smallModel_!=modelPtr_)
2673      delete smallModel_;
2674    smallModel_=NULL;
2675#else
2676    if (smallModel_==modelPtr_) {
2677      smallModel_=NULL;
2678    } else if (smallModel_) {
2679      if(!spareArrays_) {
2680        delete smallModel_;
2681        smallModel_=NULL;
2682        delete factorization_;
2683        factorization_=NULL;
2684      } else {
2685        static_cast<ClpSimplexDual *> (smallModel_)->cleanupAfterStrongBranching( factorization_);
2686        if ((smallModel_->specialOptions_&4096)==0) {
2687          delete factorization_;
2688        }
2689      }
2690    }
2691#endif
2692    //delete [] spareArrays_;
2693    //spareArrays_=NULL;
2694    factorization_=NULL;
2695  }
2696  delete [] rowActivity_;
2697  delete [] columnActivity_;
2698  rowActivity_=NULL;
2699  columnActivity_=NULL;
2700  // Make sure whatsChanged not out of sync
2701  if (!modelPtr_->columnUpperWork_)
2702    modelPtr_->whatsChanged_ &= ~0xffff;
2703  modelPtr_->specialOptions_ = saveData_.specialOptions_; 
2704}
2705
2706//#############################################################################
2707// Problem information methods (original data)
2708//#############################################################################
2709
2710//------------------------------------------------------------------
2711const char * OsiClpSolverInterface::getRowSense() const
2712{
2713  extractSenseRhsRange();
2714  return rowsense_;
2715}
2716//------------------------------------------------------------------
2717const double * OsiClpSolverInterface::getRightHandSide() const
2718{
2719  extractSenseRhsRange();
2720  return rhs_;
2721}
2722//------------------------------------------------------------------
2723const double * OsiClpSolverInterface::getRowRange() const
2724{
2725  extractSenseRhsRange();
2726  return rowrange_;
2727}
2728//------------------------------------------------------------------
2729// Return information on integrality
2730//------------------------------------------------------------------
2731bool OsiClpSolverInterface::isContinuous(int colNumber) const
2732{
2733  if ( integerInformation_==NULL ) return true;
2734#ifndef NDEBUG
2735  int n = modelPtr_->numberColumns();
2736  if (colNumber<0||colNumber>=n) {
2737    indexError(colNumber,"isContinuous");
2738  }
2739#endif
2740  if ( integerInformation_[colNumber]==0 ) return true;
2741  return false;
2742}
2743bool 
2744OsiClpSolverInterface::isBinary(int colNumber) const
2745{
2746#ifndef NDEBUG
2747  int n = modelPtr_->numberColumns();
2748  if (colNumber<0||colNumber>=n) {
2749    indexError(colNumber,"isBinary");
2750  }
2751#endif
2752  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) {
2753    return false;
2754  } else {
2755    const double * cu = getColUpper();
2756    const double * cl = getColLower();
2757    if ((cu[colNumber]== 1 || cu[colNumber]== 0) && 
2758        (cl[colNumber]== 0 || cl[colNumber]==1))
2759      return true;
2760    else 
2761      return false;
2762  }
2763}
2764//-----------------------------------------------------------------------------
2765bool 
2766OsiClpSolverInterface::isInteger(int colNumber) const
2767{
2768#ifndef NDEBUG
2769  int n = modelPtr_->numberColumns();
2770  if (colNumber<0||colNumber>=n) {
2771    indexError(colNumber,"isInteger");
2772  }
2773#endif
2774  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) 
2775    return false;
2776  else
2777    return true;
2778}
2779//-----------------------------------------------------------------------------
2780bool 
2781OsiClpSolverInterface::isIntegerNonBinary(int colNumber) const
2782{
2783#ifndef NDEBUG
2784  int n = modelPtr_->numberColumns();
2785  if (colNumber<0||colNumber>=n) {
2786    indexError(colNumber,"isIntegerNonBinary");
2787  }
2788#endif
2789  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) {
2790    return false;
2791  } else {
2792    return !isBinary(colNumber);
2793  }
2794}
2795//-----------------------------------------------------------------------------
2796bool 
2797OsiClpSolverInterface::isFreeBinary(int colNumber) const
2798{
2799#ifndef NDEBUG
2800  int n = modelPtr_->numberColumns();
2801  if (colNumber<0||colNumber>=n) {
2802    indexError(colNumber,"isFreeBinary");
2803  }
2804#endif
2805  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) {
2806    return false;
2807  } else {
2808    const double * cu = getColUpper();
2809    const double * cl = getColLower();
2810    if ((cu[colNumber]== 1) && (cl[colNumber]== 0))
2811      return true;
2812    else 
2813      return false;
2814  }
2815}
2816/*  Return array of column length
2817    0 - continuous
2818    1 - binary (may get fixed later)
2819    2 - general integer (may get fixed later)
2820*/
2821const char * 
2822OsiClpSolverInterface::getColType(bool refresh) const
2823{
2824  if (!columnType_||refresh) {
2825    const int numCols = getNumCols() ;
2826    if (!columnType_)
2827      columnType_ = new char [numCols];
2828    if ( integerInformation_==NULL ) {
2829      memset(columnType_,0,numCols);
2830    } else {
2831      const double * cu = getColUpper();
2832      const double * cl = getColLower();
2833      for (int i = 0 ; i < numCols ; ++i) {
2834        if (integerInformation_[i]) {
2835          if ((cu[i]== 1 || cu[i]== 0) && 
2836              (cl[i]== 0 || cl[i]==1))
2837            columnType_[i]=1;
2838          else
2839            columnType_[i]=2;
2840        } else {
2841          columnType_[i]=0;
2842        }
2843      }
2844    }
2845  }
2846  return columnType_;
2847}
2848
2849/* Return true if column is integer but does not have to
2850   be declared as such.
2851   Note: This function returns true if the the column
2852   is binary or a general integer.
2853*/
2854bool 
2855OsiClpSolverInterface::isOptionalInteger(int colNumber) const
2856{
2857#ifndef NDEBUG
2858  int n = modelPtr_->numberColumns();
2859  if (colNumber<0||colNumber>=n) {
2860    indexError(colNumber,"isInteger");
2861  }
2862#endif
2863  if ( integerInformation_==NULL || integerInformation_[colNumber]!=2 ) 
2864    return false;
2865  else
2866    return true;
2867}
2868/* Set the index-th variable to be an optional integer variable */
2869void 
2870OsiClpSolverInterface::setOptionalInteger(int index)
2871{
2872  if (!integerInformation_) {
2873    integerInformation_ = new char[modelPtr_->numberColumns()];
2874    CoinFillN ( integerInformation_, modelPtr_->numberColumns(),static_cast<char> (0));
2875  }
2876#ifndef NDEBUG
2877  int n = modelPtr_->numberColumns();
2878  if (index<0||index>=n) {
2879    indexError(index,"setInteger");
2880  }
2881#endif
2882  integerInformation_[index]=2;
2883  modelPtr_->setInteger(index);
2884}
2885
2886//------------------------------------------------------------------
2887// Row and column copies of the matrix ...
2888//------------------------------------------------------------------
2889const CoinPackedMatrix * OsiClpSolverInterface::getMatrixByRow() const
2890{
2891  if ( matrixByRow_ == NULL ||
2892       matrixByRow_->getNumElements() != 
2893       modelPtr_->clpMatrix()->getNumElements() ) {
2894    delete matrixByRow_;
2895    matrixByRow_ = new CoinPackedMatrix();
2896    matrixByRow_->setExtraGap(0.0);
2897    matrixByRow_->setExtraMajor(0.0);
2898    matrixByRow_->reverseOrderedCopyOf(*modelPtr_->matrix());
2899    //matrixByRow_->removeGaps();
2900#if 0
2901    CoinPackedMatrix back;
2902    std::cout<<"start check"<<std::endl;
2903    back.reverseOrderedCopyOf(*matrixByRow_);
2904    modelPtr_->matrix()->isEquivalent2(back);
2905    std::cout<<"stop check"<<std::endl;
2906#endif
2907  }
2908  assert (matrixByRow_->getNumElements()==modelPtr_->clpMatrix()->getNumElements());
2909  return matrixByRow_;
2910}
2911
2912const CoinPackedMatrix * OsiClpSolverInterface::getMatrixByCol() const
2913{
2914  return modelPtr_->matrix();
2915}
2916 
2917// Get pointer to mutable column-wise copy of matrix (returns NULL if not meaningful)
2918CoinPackedMatrix * 
2919OsiClpSolverInterface::getMutableMatrixByCol() const 
2920{
2921  ClpPackedMatrix * matrix = dynamic_cast<ClpPackedMatrix *>(modelPtr_->matrix_) ;
2922  if (matrix)
2923    return matrix->getPackedMatrix();
2924  else
2925    return NULL;
2926}
2927
2928//------------------------------------------------------------------
2929std::vector<double*> OsiClpSolverInterface::getDualRays(int /*maxNumRays*/,
2930                                                        bool fullRay) const
2931{
2932  return std::vector<double*>(1, modelPtr_->infeasibilityRay(fullRay));
2933}
2934//------------------------------------------------------------------
2935std::vector<double*> OsiClpSolverInterface::getPrimalRays(int /*maxNumRays*/) const
2936{
2937  return std::vector<double*>(1, modelPtr_->unboundedRay());
2938}
2939//#############################################################################
2940void
2941OsiClpSolverInterface::setContinuous(int index)
2942{
2943
2944  if (integerInformation_) {
2945#ifndef NDEBUG
2946    int n = modelPtr_->numberColumns();
2947    if (index<0||index>=n) {
2948      indexError(index,"setContinuous");
2949    }
2950#endif
2951    integerInformation_[index]=0;
2952  }
2953  modelPtr_->setContinuous(index);
2954}
2955//-----------------------------------------------------------------------------
2956void
2957OsiClpSolverInterface::setInteger(int index)
2958{
2959  if (!integerInformation_) {
2960    integerInformation_ = new char[modelPtr_->numberColumns()];
2961    CoinFillN ( integerInformation_, modelPtr_->numberColumns(),static_cast<char> (0));
2962  }
2963#ifndef NDEBUG
2964  int n = modelPtr_->numberColumns();
2965  if (index<0||index>=n) {
2966    indexError(index,"setInteger");
2967  }
2968#endif
2969  integerInformation_[index]=1;
2970  modelPtr_->setInteger(index);
2971}
2972//-----------------------------------------------------------------------------
2973void
2974OsiClpSolverInterface::setContinuous(const int* indices, int len)
2975{
2976  if (integerInformation_) {
2977#ifndef NDEBUG
2978    int n = modelPtr_->numberColumns();
2979#endif
2980    int i;
2981    for (i=0; i<len;i++) {
2982      int colNumber = indices[i];
2983#ifndef NDEBUG
2984      if (colNumber<0||colNumber>=n) {
2985        indexError(colNumber,"setContinuous");
2986      }
2987#endif
2988      integerInformation_[colNumber]=0;
2989      modelPtr_->setContinuous(colNumber);
2990    }
2991  }
2992}
2993//-----------------------------------------------------------------------------
2994void
2995OsiClpSolverInterface::setInteger(const int* indices, int len)
2996{
2997  if (!integerInformation_) {
2998    integerInformation_ = new char[modelPtr_->numberColumns()];
2999    CoinFillN ( integerInformation_, modelPtr_->numberColumns(),static_cast<char> (0));
3000  }
3001#ifndef NDEBUG
3002  int n = modelPtr_->numberColumns();
3003#endif
3004  int i;
3005  for (i=0; i<len;i++) {
3006    int colNumber = indices[i];
3007#ifndef NDEBUG
3008    if (colNumber<0||colNumber>=n) {
3009      indexError(colNumber,"setInteger");
3010    }
3011#endif
3012    integerInformation_[colNumber]=1;
3013    modelPtr_->setInteger(colNumber);
3014  }
3015}
3016/* Set the objective coefficients for all columns
3017    array [getNumCols()] is an array of values for the objective.
3018    This defaults to a series of set operations and is here for speed.
3019*/
3020void 
3021OsiClpSolverInterface::setObjective(const double * array)
3022{
3023  // Say can't gurantee optimal basis etc
3024  lastAlgorithm_=999;
3025  modelPtr_->whatsChanged_ &= (0xffff&~64);
3026  int n = modelPtr_->numberColumns() ;
3027  if (fakeMinInSimplex_) {
3028    std::transform(array,array+n,
3029                   modelPtr_->objective(),std::negate<double>()) ;
3030  } else {
3031    CoinMemcpyN(array,n,modelPtr_->objective());
3032  }
3033}
3034/* Set the lower bounds for all columns
3035    array [getNumCols()] is an array of values for the objective.
3036    This defaults to a series of set operations and is here for speed.
3037*/
3038void 
3039OsiClpSolverInterface::setColLower(const double * array)
3040{
3041  // Say can't gurantee optimal basis etc
3042  lastAlgorithm_=999;
3043  modelPtr_->whatsChanged_ &= (0x1ffff&128);
3044  CoinMemcpyN(array,modelPtr_->numberColumns(),
3045                    modelPtr_->columnLower());
3046}
3047/* Set the upper bounds for all columns
3048    array [getNumCols()] is an array of values for the objective.
3049    This defaults to a series of set operations and is here for speed.
3050*/
3051void 
3052OsiClpSolverInterface::setColUpper(const double * array)
3053{
3054  // Say can't gurantee optimal basis etc
3055  lastAlgorithm_=999;
3056  modelPtr_->whatsChanged_ &= (0x1ffff&256);
3057  CoinMemcpyN(array,modelPtr_->numberColumns(),
3058                    modelPtr_->columnUpper());
3059}
3060//-----------------------------------------------------------------------------
3061void OsiClpSolverInterface::setColSolution(const double * cs) 
3062{
3063  // Say can't gurantee optimal basis etc
3064  lastAlgorithm_=999;
3065  CoinDisjointCopyN(cs,modelPtr_->numberColumns(),
3066                    modelPtr_->primalColumnSolution());
3067  if (modelPtr_->solveType()==2) {
3068    // directly into code as well
3069    CoinDisjointCopyN(cs,modelPtr_->numberColumns(),
3070                      modelPtr_->solutionRegion(1));
3071  }
3072  // compute row activity
3073  memset(modelPtr_->primalRowSolution(),0,modelPtr_->numberRows()*sizeof(double));
3074  modelPtr_->times(1.0,modelPtr_->primalColumnSolution(),modelPtr_->primalRowSolution());
3075}
3076//-----------------------------------------------------------------------------
3077void OsiClpSolverInterface::setRowPrice(const double * rs) 
3078{
3079  CoinDisjointCopyN(rs,modelPtr_->numberRows(),
3080                    modelPtr_->dualRowSolution());
3081  if (modelPtr_->solveType()==2) {
3082    // directly into code as well (? sign )
3083    CoinDisjointCopyN(rs,modelPtr_->numberRows(),
3084                      modelPtr_->djRegion(0));
3085  }
3086  // compute reduced costs
3087  memcpy(modelPtr_->dualColumnSolution(),modelPtr_->objective(),
3088         modelPtr_->numberColumns()*sizeof(double));
3089  modelPtr_->transposeTimes(-1.0,modelPtr_->dualRowSolution(),modelPtr_->dualColumnSolution());
3090}
3091
3092//#############################################################################
3093// Problem modifying methods (matrix)
3094//#############################################################################
3095void 
3096OsiClpSolverInterface::addCol(const CoinPackedVectorBase& vec,
3097                              const double collb, const double colub,   
3098                              const double obj)
3099{
3100  int numberColumns = modelPtr_->numberColumns();
3101  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3102  modelPtr_->resize(modelPtr_->numberRows(),numberColumns+1);
3103  linearObjective_ = modelPtr_->objective();
3104  basis_.resize(modelPtr_->numberRows(),numberColumns+1);
3105  setColBounds(numberColumns,collb,colub);
3106  setObjCoeff(numberColumns,obj);
3107  if (!modelPtr_->clpMatrix())
3108    modelPtr_->createEmptyMatrix();
3109  modelPtr_->matrix()->appendCol(vec);
3110  if (integerInformation_) {
3111    char * temp = new char[numberColumns+1];
3112    CoinMemcpyN(integerInformation_,numberColumns,temp);
3113    delete [] integerInformation_;
3114    integerInformation_ = temp;
3115    integerInformation_[numberColumns]=0;
3116  }
3117  freeCachedResults();
3118}
3119//-----------------------------------------------------------------------------
3120/* Add a column (primal variable) to the problem. */
3121void 
3122OsiClpSolverInterface::addCol(int numberElements, const int * rows, const double * elements,
3123                           const double collb, const double colub,   
3124                           const double obj) 
3125{
3126  CoinPackedVector column(numberElements, rows, elements);
3127  addCol(column,collb,colub,obj);
3128}
3129// Add a named column (primal variable) to the problem.
3130void 
3131OsiClpSolverInterface::addCol(const CoinPackedVectorBase& vec,
3132                      const double collb, const double colub,   
3133                      const double obj, std::string name) 
3134{
3135  int ndx = getNumCols() ;
3136  addCol(vec,collb,colub,obj) ;
3137  setColName(ndx,name) ;
3138}
3139//-----------------------------------------------------------------------------
3140void 
3141OsiClpSolverInterface::addCols(const int numcols,
3142                               const CoinPackedVectorBase * const * cols,
3143                               const double* collb, const double* colub,   
3144                               const double* obj)
3145{
3146  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3147  int numberColumns = modelPtr_->numberColumns();
3148  modelPtr_->resize(modelPtr_->numberRows(),numberColumns+numcols);
3149  linearObjective_ = modelPtr_->objective();
3150  basis_.resize(modelPtr_->numberRows(),numberColumns+numcols);
3151  double * lower = modelPtr_->columnLower()+numberColumns;
3152  double * upper = modelPtr_->columnUpper()+numberColumns;
3153  double * objective = modelPtr_->objective()+numberColumns;
3154  int iCol;
3155  if (collb) {
3156    for (iCol = 0; iCol < numcols; iCol++) {
3157      lower[iCol]= forceIntoRange(collb[iCol], -OsiClpInfinity, OsiClpInfinity);
3158      if (lower[iCol]<-1.0e27)
3159        lower[iCol]=-COIN_DBL_MAX;
3160    }
3161  } else {
3162    CoinFillN ( lower, numcols,0.0);
3163  }
3164  if (colub) {
3165    for (iCol = 0; iCol < numcols; iCol++) {
3166      upper[iCol]= forceIntoRange(colub[iCol], -OsiClpInfinity, OsiClpInfinity);
3167      if (upper[iCol]>1.0e27)
3168        upper[iCol]=COIN_DBL_MAX;
3169    }
3170  } else {
3171    CoinFillN ( upper, numcols,COIN_DBL_MAX);
3172  }
3173  if (obj) {
3174    for (iCol = 0; iCol < numcols; iCol++) {
3175      objective[iCol] = obj[iCol];
3176    }
3177  } else {
3178    CoinFillN ( objective, numcols,0.0);
3179  }
3180  if (!modelPtr_->clpMatrix())
3181    modelPtr_->createEmptyMatrix();
3182  modelPtr_->matrix()->appendCols(numcols,cols);
3183  if (integerInformation_) {
3184    char * temp = new char[numberColumns+numcols];
3185    CoinMemcpyN(integerInformation_,numberColumns,temp);
3186    delete [] integerInformation_;
3187    integerInformation_ = temp;
3188    for (iCol = 0; iCol < numcols; iCol++) 
3189      integerInformation_[numberColumns+iCol]=0;
3190  }
3191  freeCachedResults();
3192}
3193void 
3194OsiClpSolverInterface::addCols(const int numcols,
3195                               const int * columnStarts, const int * rows, const double * elements,
3196                               const double* collb, const double* colub,   
3197                               const double* obj)
3198{
3199  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3200  int numberColumns = modelPtr_->numberColumns();
3201  modelPtr_->resize(modelPtr_->numberRows(),numberColumns+numcols);
3202  linearObjective_ = modelPtr_->objective();
3203  basis_.resize(modelPtr_->numberRows(),numberColumns+numcols);
3204  double * lower = modelPtr_->columnLower()+numberColumns;
3205  double * upper = modelPtr_->columnUpper()+numberColumns;
3206  double * objective = modelPtr_->objective()+numberColumns;
3207  int iCol;
3208  if (collb) {
3209    for (iCol = 0; iCol < numcols; iCol++) {
3210      lower[iCol]= forceIntoRange(collb[iCol], -OsiClpInfinity, OsiClpInfinity);
3211      if (lower[iCol]<-1.0e27)
3212        lower[iCol]=-COIN_DBL_MAX;
3213    }
3214  } else {
3215    CoinFillN ( lower, numcols,0.0);
3216  }
3217  if (colub) {
3218    for (iCol = 0; iCol < numcols; iCol++) {
3219      upper[iCol]= forceIntoRange(colub[iCol], -OsiClpInfinity, OsiClpInfinity);
3220      if (upper[iCol]>1.0e27)
3221        upper[iCol]=COIN_DBL_MAX;
3222    }
3223  } else {
3224    CoinFillN ( upper, numcols,COIN_DBL_MAX);
3225  }
3226  if (obj) {
3227    for (iCol = 0; iCol < numcols; iCol++) {
3228      objective[iCol] = obj[iCol];
3229    }
3230  } else {
3231    CoinFillN ( objective, numcols,0.0);
3232  }
3233  if (!modelPtr_->clpMatrix())
3234    modelPtr_->createEmptyMatrix();
3235  modelPtr_->matrix()->appendCols(numcols,columnStarts,rows,elements);
3236  if (integerInformation_) {
3237    char * temp = new char[numberColumns+numcols];
3238    CoinMemcpyN(integerInformation_,numberColumns,temp);
3239    delete [] integerInformation_;
3240    integerInformation_ = temp;
3241    for (iCol = 0; iCol < numcols; iCol++) 
3242      integerInformation_[numberColumns+iCol]=0;
3243  }
3244  freeCachedResults();
3245}
3246//-----------------------------------------------------------------------------
3247void 
3248OsiClpSolverInterface::deleteCols(const int num, const int * columnIndices)
3249{
3250  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3251  findIntegers(false);
3252  deleteBranchingInfo(num,columnIndices);
3253  modelPtr_->deleteColumns(num,columnIndices);
3254  int nameDiscipline;
3255  getIntParam(OsiNameDiscipline,nameDiscipline) ;
3256  if (num&&nameDiscipline) {
3257    // Very clumsy (and inefficient) - need to sort and then go backwards in ? chunks
3258    int * indices = CoinCopyOfArray(columnIndices,num);
3259    std::sort(indices,indices+num);
3260    int num2=num;
3261    while(num2) {
3262      int next = indices[num2-1];
3263      int firstDelete = num2-1;
3264      int i;
3265      for (i=num2-2;i>=0;i--) {
3266        if (indices[i]+1==next) {
3267          next --;
3268          firstDelete=i;
3269        } else {
3270          break;
3271        }
3272      }
3273      OsiSolverInterface::deleteColNames(indices[firstDelete],num2-firstDelete);
3274      num2 = firstDelete;
3275      assert (num2>=0);
3276    }
3277    delete [] indices;
3278  }
3279  // synchronize integers (again)
3280  if (integerInformation_) {
3281    int numberColumns = modelPtr_->numberColumns();
3282    for (int i=0;i<numberColumns;i++) {
3283      if (modelPtr_->isInteger(i))
3284        integerInformation_[i]=1;
3285      else
3286        integerInformation_[i]=0;
3287    }
3288  }
3289  basis_.deleteColumns(num,columnIndices);
3290  linearObjective_ = modelPtr_->objective();
3291  freeCachedResults();
3292}
3293//-----------------------------------------------------------------------------
3294void 
3295OsiClpSolverInterface::addRow(const CoinPackedVectorBase& vec,
3296                              const double rowlb, const double rowub)
3297{
3298  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3299  freeCachedResults0();
3300  int numberRows = modelPtr_->numberRows();
3301  modelPtr_->resize(numberRows+1,modelPtr_->numberColumns());
3302  basis_.resize(numberRows+1,modelPtr_->numberColumns());
3303  setRowBounds(numberRows,rowlb,rowub);
3304  if (!modelPtr_->clpMatrix())
3305    modelPtr_->createEmptyMatrix();
3306  modelPtr_->matrix()->appendRow(vec);
3307  freeCachedResults1();
3308}
3309//-----------------------------------------------------------------------------
3310void OsiClpSolverInterface::addRow(const CoinPackedVectorBase& vec,
3311                                const double rowlb, const double rowub,
3312                                std::string name)
3313{
3314  int ndx = getNumRows() ;
3315  addRow(vec,rowlb,rowub) ;
3316  setRowName(ndx,name) ;
3317}
3318//-----------------------------------------------------------------------------
3319void 
3320OsiClpSolverInterface::addRow(const CoinPackedVectorBase& vec,
3321                              const char rowsen, const double rowrhs,   
3322                              const double rowrng)
3323{
3324  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3325  freeCachedResults0();
3326  int numberRows = modelPtr_->numberRows();
3327  modelPtr_->resize(numberRows+1,modelPtr_->numberColumns());
3328  basis_.resize(numberRows+1,modelPtr_->numberColumns());
3329  double rowlb = 0, rowub = 0;
3330  convertSenseToBound(rowsen, rowrhs, rowrng, rowlb, rowub);
3331  setRowBounds(numberRows,rowlb,rowub);
3332  if (!modelPtr_->clpMatrix())
3333    modelPtr_->createEmptyMatrix();
3334  modelPtr_->matrix()->appendRow(vec);
3335  freeCachedResults1();
3336}
3337//-----------------------------------------------------------------------------
3338void 
3339OsiClpSolverInterface::addRow(int numberElements, const int * columns, const double * elements,
3340                           const double rowlb, const double rowub) 
3341{
3342  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3343  freeCachedResults0();
3344  int numberRows = modelPtr_->numberRows();
3345  modelPtr_->resize(numberRows+1,modelPtr_->numberColumns());
3346  basis_.resize(numberRows+1,modelPtr_->numberColumns());
3347  setRowBounds(numberRows,rowlb,rowub);
3348  if (!modelPtr_->clpMatrix())
3349    modelPtr_->createEmptyMatrix();
3350  modelPtr_->matrix()->appendRow(numberElements, columns, elements);
3351  CoinBigIndex starts[2];
3352  starts[0]=0;
3353  starts[1]=numberElements;
3354  redoScaleFactors( 1,starts, columns, elements);
3355  freeCachedResults1();
3356}
3357//-----------------------------------------------------------------------------
3358void 
3359OsiClpSolverInterface::addRows(const int numrows,
3360                               const CoinPackedVectorBase * const * rows,
3361                               const double* rowlb, const double* rowub)
3362{
3363  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3364  freeCachedResults0();
3365  int numberRows = modelPtr_->numberRows();
3366  modelPtr_->resize(numberRows+numrows,modelPtr_->numberColumns());
3367  basis_.resize(numberRows+numrows,modelPtr_->numberColumns());
3368  double * lower = modelPtr_->rowLower()+numberRows;
3369  double * upper = modelPtr_->rowUpper()+numberRows;
3370  int iRow;
3371  for (iRow = 0; iRow < numrows; iRow++) {
3372    if (rowlb) 
3373      lower[iRow]= forceIntoRange(rowlb[iRow], -OsiClpInfinity, OsiClpInfinity);
3374    else 
3375      lower[iRow]=-OsiClpInfinity;
3376    if (rowub) 
3377      upper[iRow]= forceIntoRange(rowub[iRow], -OsiClpInfinity, OsiClpInfinity);
3378    else 
3379      upper[iRow]=OsiClpInfinity;
3380    if (lower[iRow]<-1.0e27)
3381      lower[iRow]=-COIN_DBL_MAX;
3382    if (upper[iRow]>1.0e27)
3383      upper[iRow]=COIN_DBL_MAX;
3384  }
3385  if (!modelPtr_->clpMatrix())
3386    modelPtr_->createEmptyMatrix();
3387  modelPtr_->matrix()->appendRows(numrows,rows);
3388  freeCachedResults1();
3389}
3390//-----------------------------------------------------------------------------
3391void 
3392OsiClpSolverInterface::addRows(const int numrows,
3393                               const CoinPackedVectorBase * const * rows,
3394                               const char* rowsen, const double* rowrhs,   
3395                               const double* rowrng)
3396{
3397  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3398  freeCachedResults0();
3399  int numberRows = modelPtr_->numberRows();
3400  modelPtr_->resize(numberRows+numrows,modelPtr_->numberColumns());
3401  basis_.resize(numberRows+numrows,modelPtr_->numberColumns());
3402  double * lower = modelPtr_->rowLower()+numberRows;
3403  double * upper = modelPtr_->rowUpper()+numberRows;
3404  int iRow;
3405  for (iRow = 0; iRow < numrows; iRow++) {
3406    double rowlb = 0, rowub = 0;
3407    convertSenseToBound(rowsen[iRow], rowrhs[iRow], rowrng[iRow], 
3408                        rowlb, rowub);
3409    lower[iRow]= forceIntoRange(rowlb, -OsiClpInfinity, OsiClpInfinity);
3410    upper[iRow]= forceIntoRange(rowub, -OsiClpInfinity, OsiClpInfinity);
3411    if (lower[iRow]<-1.0e27)
3412      lower[iRow]=-COIN_DBL_MAX;
3413    if (upper[iRow]>1.0e27)
3414      upper[iRow]=COIN_DBL_MAX;
3415  }
3416  if (!modelPtr_->clpMatrix())
3417    modelPtr_->createEmptyMatrix();
3418  modelPtr_->matrix()->appendRows(numrows,rows);
3419  freeCachedResults1();
3420}
3421void 
3422OsiClpSolverInterface::addRows(const int numrows,
3423                               const int * rowStarts, const int * columns, const double * element,
3424                               const double* rowlb, const double* rowub)
3425{
3426  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3427  freeCachedResults0();
3428  int numberRows = modelPtr_->numberRows();
3429  modelPtr_->resize(numberRows+numrows,modelPtr_->numberColumns());
3430  basis_.resize(numberRows+numrows,modelPtr_->numberColumns());
3431  double * lower = modelPtr_->rowLower()+numberRows;
3432  double * upper = modelPtr_->rowUpper()+numberRows;
3433  int iRow;
3434  for (iRow = 0; iRow < numrows; iRow++) {
3435    if (rowlb) 
3436      lower[iRow]= forceIntoRange(rowlb[iRow], -OsiClpInfinity, OsiClpInfinity);
3437    else 
3438      lower[iRow]=-OsiClpInfinity;
3439    if (rowub) 
3440      upper[iRow]= forceIntoRange(rowub[iRow], -OsiClpInfinity, OsiClpInfinity);
3441    else 
3442      upper[iRow]=OsiClpInfinity;
3443    if (lower[iRow]<-1.0e27)
3444      lower[iRow]=-COIN_DBL_MAX;
3445    if (upper[iRow]>1.0e27)
3446      upper[iRow]=COIN_DBL_MAX;
3447  }
3448  if (!modelPtr_->clpMatrix())
3449    modelPtr_->createEmptyMatrix();
3450  modelPtr_->matrix()->appendRows(numrows,rowStarts,columns,element);
3451  redoScaleFactors( numrows,rowStarts, columns, element);
3452  freeCachedResults1();
3453}
3454//-----------------------------------------------------------------------------
3455void 
3456OsiClpSolverInterface::deleteRows(const int num, const int * rowIndices)
3457{
3458  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3459  // will still be optimal if all rows basic
3460  bool allBasic=true;
3461  int numBasis = basis_.getNumArtificial();
3462  for (int i=0;i<num;i++) {
3463    int iRow = rowIndices[i];
3464    if (iRow<numBasis) {
3465      if (basis_.getArtifStatus(iRow)!=CoinWarmStartBasis::basic) {
3466        allBasic=false;
3467        break;
3468      }
3469    }
3470  }
3471  int saveAlgorithm = allBasic ? lastAlgorithm_ : 999;
3472  modelPtr_->deleteRows(num,rowIndices);
3473  int nameDiscipline;
3474  getIntParam(OsiNameDiscipline,nameDiscipline) ;
3475  if (num&&nameDiscipline) {
3476    // Very clumsy (and inefficient) - need to sort and then go backwards in ? chunks
3477    int * indices = CoinCopyOfArray(rowIndices,num);
3478    std::sort(indices,indices+num);
3479    int num2=num;
3480    while(num2) {
3481      int next = indices[num2-1];
3482      int firstDelete = num2-1;
3483      int i;
3484      for (i=num2-2;i>=0;i--) {
3485        if (indices[i]+1==next) {
3486          next --;
3487          firstDelete=i;
3488        } else {
3489          break;
3490        }
3491      }
3492      OsiSolverInterface::deleteRowNames(indices[firstDelete],num2-firstDelete);
3493      num2 = firstDelete;
3494      assert (num2>=0);
3495    }
3496    delete [] indices;
3497  }
3498  basis_.deleteRows(num,rowIndices);
3499  CoinPackedMatrix * saveRowCopy = matrixByRow_;
3500  matrixByRow_=NULL;
3501  freeCachedResults();
3502  modelPtr_->setNewRowCopy(NULL);
3503  delete modelPtr_->scaledMatrix_;
3504  modelPtr_->scaledMatrix_=NULL;
3505  if (saveRowCopy) {
3506    matrixByRow_=saveRowCopy;
3507    matrixByRow_->deleteRows(num,rowIndices);
3508    if (matrixByRow_->getNumElements()!=modelPtr_->clpMatrix()->getNumElements()) {
3509      delete matrixByRow_; // odd type matrix
3510      matrixByRow_=NULL;
3511    }
3512  }
3513  lastAlgorithm_ = saveAlgorithm;
3514  if ((specialOptions_&131072)!=0) 
3515    lastNumberRows_=modelPtr_->numberRows();
3516}
3517
3518//#############################################################################
3519// Methods to input a problem
3520//#############################################################################
3521
3522void
3523OsiClpSolverInterface::loadProblem(const CoinPackedMatrix& matrix,
3524                                   const double* collb, const double* colub,   
3525                                   const double* obj,
3526                                   const double* rowlb, const double* rowub)
3527{
3528  modelPtr_->whatsChanged_ = 0;
3529  // Get rid of integer information (modelPtr will get rid of its copy)
3530  delete [] integerInformation_;
3531  integerInformation_=NULL;
3532  modelPtr_->loadProblem(matrix, collb, colub, obj, rowlb, rowub);
3533  linearObjective_ = modelPtr_->objective();
3534  freeCachedResults();
3535  basis_=CoinWarmStartBasis();
3536  if (ws_) {
3537     delete ws_;
3538     ws_ = 0;
3539  }
3540}
3541
3542//-----------------------------------------------------------------------------
3543
3544/*
3545  Expose the method that takes ClpMatrixBase. User request. Can't hurt, given
3546  the number of non-OSI methods already here.
3547*/
3548void OsiClpSolverInterface::loadProblem (const ClpMatrixBase& matrix,
3549                                         const double* collb,
3550                                         const double* colub,   
3551                                         const double* obj,
3552                                         const double* rowlb,
3553                                         const double* rowub)
3554{
3555  modelPtr_->whatsChanged_ = 0;
3556  // Get rid of integer information (modelPtr will get rid of its copy)
3557  delete [] integerInformation_;
3558  integerInformation_=NULL;
3559  modelPtr_->loadProblem(matrix,collb,colub,obj,rowlb,rowub);
3560  linearObjective_ = modelPtr_->objective();
3561  freeCachedResults();
3562  basis_=CoinWarmStartBasis();
3563  if (ws_) {
3564     delete ws_;
3565     ws_ = 0;
3566  }
3567}
3568
3569//-----------------------------------------------------------------------------
3570
3571void
3572OsiClpSolverInterface::assignProblem(CoinPackedMatrix*& matrix,
3573                                     double*& collb, double*& colub,
3574                                     double*& obj,
3575                                     double*& rowlb, double*& rowub)
3576{
3577  modelPtr_->whatsChanged_ = 0;
3578  // Get rid of integer information (modelPtr will get rid of its copy)
3579  loadProblem(*matrix, collb, colub, obj, rowlb, rowub);
3580  delete matrix;   matrix = NULL;
3581  delete[] collb;  collb = NULL;
3582  delete[] colub;  colub = NULL;
3583  delete[] obj;    obj = NULL;
3584  delete[] rowlb;  rowlb = NULL;
3585  delete[] rowub;  rowub = NULL;
3586}
3587
3588//-----------------------------------------------------------------------------
3589
3590void
3591OsiClpSolverInterface::loadProblem(const CoinPackedMatrix& matrix,
3592                                   const double* collb, const double* colub,
3593                                   const double* obj,
3594                                   const char* rowsen, const double* rowrhs,   
3595                                   const double* rowrng)
3596{
3597  modelPtr_->whatsChanged_ = 0;
3598  // Get rid of integer information (modelPtr will get rid of its copy)
3599  // assert( rowsen != NULL );
3600  // assert( rowrhs != NULL );
3601  // If any of Rhs NULLs then create arrays
3602  int numrows = matrix.getNumRows();
3603  const char * rowsenUse = rowsen;
3604  if (!rowsen) {
3605    char * rowsen = new char [numrows];
3606    for (int i=0;i<numrows;i++)
3607      rowsen[i]='G';
3608    rowsenUse = rowsen;
3609  } 
3610  const double * rowrhsUse = rowrhs;
3611  if (!rowrhs) {
3612    double * rowrhs = new double [numrows];
3613    for (int i=0;i<numrows;i++)
3614      rowrhs[i]=0.0;
3615    rowrhsUse = rowrhs;
3616  }
3617  const double * rowrngUse = rowrng;
3618  if (!rowrng) {
3619    double * rowrng = new double [numrows];
3620    for (int i=0;i<numrows;i++)
3621      rowrng[i]=0.0;
3622    rowrngUse = rowrng;
3623  }
3624  double * rowlb = new double[numrows];
3625  double * rowub = new double[numrows];
3626  for (int i = numrows-1; i >= 0; --i) {   
3627    convertSenseToBound(rowsenUse[i],rowrhsUse[i],rowrngUse[i],rowlb[i],rowub[i]);
3628  }
3629  if (rowsen!=rowsenUse)
3630    delete [] rowsenUse;
3631  if (rowrhs!=rowrhsUse)
3632    delete [] rowrhsUse;
3633  if (rowrng!=rowrngUse)
3634    delete [] rowrngUse;
3635  loadProblem(matrix, collb, colub, obj, rowlb, rowub);
3636  delete [] rowlb;
3637  delete [] rowub;
3638}
3639
3640//-----------------------------------------------------------------------------
3641
3642void
3643OsiClpSolverInterface::assignProblem(CoinPackedMatrix*& matrix,
3644                                     double*& collb, double*& colub,
3645                                     double*& obj,
3646                                     char*& rowsen, double*& rowrhs,
3647                                     double*& rowrng)
3648{
3649  modelPtr_->whatsChanged_ = 0;
3650  // Get rid of integer information (modelPtr will get rid of its copy)
3651  loadProblem(*matrix, collb, colub, obj, rowsen, rowrhs, rowrng);
3652  delete matrix;   matrix = NULL;
3653  delete[] collb;  collb = NULL;
3654  delete[] colub;  colub = NULL;
3655  delete[] obj;    obj = NULL;
3656  delete[] rowsen; rowsen = NULL;
3657  delete[] rowrhs; rowrhs = NULL;
3658  delete[] rowrng; rowrng = NULL;
3659}
3660
3661//-----------------------------------------------------------------------------
3662
3663void
3664OsiClpSolverInterface::loadProblem(const int numcols, const int numrows,
3665                                   const CoinBigIndex * start, const int* index,
3666                                   const double* value,
3667                                   const double* collb, const double* colub,
3668                                   const double* obj,
3669                                   const double* rowlb, const double* rowub)
3670{
3671  modelPtr_->whatsChanged_ = 0;
3672  // Get rid of integer information (modelPtr will get rid of its copy)
3673  delete [] integerInformation_;
3674  integerInformation_=NULL;
3675  modelPtr_->loadProblem(numcols, numrows, start,  index,
3676            value, collb, colub, obj,
3677            rowlb,  rowub);
3678  linearObjective_ = modelPtr_->objective();
3679  freeCachedResults();
3680  basis_=CoinWarmStartBasis();
3681  if (ws_) {
3682     delete ws_;
3683     ws_ = 0;
3684  }
3685}
3686//-----------------------------------------------------------------------------
3687
3688void
3689OsiClpSolverInterface::loadProblem(const int numcols, const int numrows,
3690                                   const CoinBigIndex * start, const int* index,
3691                                   const double* value,
3692                                   const double* collb, const double* colub,
3693                                   const double* obj,
3694                                   const char* rowsen, const double* rowrhs,   
3695                                   const double* rowrng)
3696{
3697  modelPtr_->whatsChanged_ = 0;
3698  // Get rid of integer information (modelPtr will get rid of its copy)
3699  // If any of Rhs NULLs then create arrays
3700  const char * rowsenUse = rowsen;
3701  if (!rowsen) {
3702    char * rowsen = new char [numrows];
3703    for (int i=0;i<numrows;i++)
3704      rowsen[i]='G';
3705    rowsenUse = rowsen;
3706  } 
3707  const double * rowrhsUse = rowrhs;
3708  if (!rowrhs) {
3709    double * rowrhs = new double [numrows];
3710    for (int i=0;i<numrows;i++)
3711      rowrhs[i]=0.0;
3712    rowrhsUse = rowrhs;
3713  }
3714  const double * rowrngUse = rowrng;
3715  if (!rowrng) {
3716    double * rowrng = new double [numrows];
3717    for (int i=0;i<numrows;i++)
3718      rowrng[i]=0.0;
3719    rowrngUse = rowrng;
3720  }
3721  double * rowlb = new double[numrows];
3722  double * rowub = new double[numrows];
3723  for (int i = numrows-1; i >= 0; --i) {   
3724    convertSenseToBound(rowsenUse[i],rowrhsUse[i],rowrngUse[i],rowlb[i],rowub[i]);
3725  }
3726  if (rowsen!=rowsenUse)
3727    delete [] rowsenUse;
3728  if (rowrhs!=rowrhsUse)
3729    delete [] rowrhsUse;
3730  if (rowrng!=rowrngUse)
3731    delete [] rowrngUse;
3732  loadProblem(numcols, numrows, start,  index, value, collb, colub, obj,
3733              rowlb,  rowub);
3734  delete[] rowlb;
3735  delete[] rowub;
3736}
3737// This loads a model from a coinModel object - returns number of errors
3738int 
3739OsiClpSolverInterface::loadFromCoinModel (  CoinModel & modelObject, bool keepSolution)
3740{
3741  modelPtr_->whatsChanged_ = 0;
3742  int numberErrors = 0;
3743  // Set arrays for normal use
3744  double * rowLower = modelObject.rowLowerArray();
3745  double * rowUpper = modelObject.rowUpperArray();
3746  double * columnLower = modelObject.columnLowerArray();
3747  double * columnUpper = modelObject.columnUpperArray();
3748  double * objective = modelObject.objectiveArray();
3749  int * integerType = modelObject.integerTypeArray();
3750  double * associated = modelObject.associatedArray();
3751  // If strings then do copies
3752  if (modelObject.stringsExist()) {
3753    numberErrors = modelObject.createArrays(rowLower, rowUpper, columnLower, columnUpper,
3754                                            objective, integerType,associated);
3755  }
3756  CoinPackedMatrix matrix;
3757  modelObject.createPackedMatrix(matrix,associated);
3758  int numberRows = modelObject.numberRows();
3759  int numberColumns = modelObject.numberColumns();
3760  CoinWarmStart * ws = getWarmStart();
3761  bool restoreBasis = keepSolution && numberRows&&numberRows==getNumRows()&&
3762    numberColumns==getNumCols();
3763  loadProblem(matrix, 
3764              columnLower, columnUpper, objective, rowLower, rowUpper);
3765  if (restoreBasis)
3766    setWarmStart(ws);
3767  delete ws;
3768  // Do names if wanted
3769  int numberItems;
3770  numberItems = modelObject.rowNames()->numberItems();
3771  if (numberItems) {
3772    const char *const * rowNames=modelObject.rowNames()->names();
3773    modelPtr_->copyRowNames(rowNames,0,numberItems);
3774  }
3775  numberItems = modelObject.columnNames()->numberItems();
3776  if (numberItems) {
3777    const char *const * columnNames=modelObject.columnNames()->names();
3778    modelPtr_->copyColumnNames(columnNames,0,numberItems);
3779  }
3780  // Do integers if wanted
3781  assert(integerType);
3782  for (int iColumn=0;iColumn<numberColumns;iColumn++) {
3783    if (integerType[iColumn])
3784      setInteger(iColumn);
3785  }
3786  if (rowLower!=modelObject.rowLowerArray()||
3787      columnLower!=modelObject.columnLowerArray()) {
3788    delete [] rowLower;
3789    delete [] rowUpper;
3790    delete [] columnLower;
3791    delete [] columnUpper;
3792    delete [] objective;
3793    delete [] integerType;
3794    delete [] associated;
3795    //if (numberErrors)
3796    //  handler_->message(CLP_BAD_STRING_VALUES,messages_)
3797    //    <<numberErrors
3798    //    <<CoinMessageEol;
3799  }
3800  modelPtr_->optimizationDirection_ = modelObject.optimizationDirection(); 
3801  return numberErrors;
3802}
3803
3804//-----------------------------------------------------------------------------
3805// Write mps files
3806//-----------------------------------------------------------------------------
3807
3808void OsiClpSolverInterface::writeMps(const char * filename,
3809                                     const char * extension,
3810                                     double objSense) const
3811{
3812  std::string f(filename);
3813  std::string e(extension);
3814  std::string fullname;
3815  if (e!="") {
3816    fullname = f + "." + e;
3817  } else {
3818    // no extension so no trailing period
3819    fullname = f;
3820  }
3821  // get names
3822  const char * const * const rowNames = modelPtr_->rowNamesAsChar();
3823  const char * const * const columnNames = modelPtr_->columnNamesAsChar();
3824  // Fall back on Osi version - possibly with names
3825  OsiSolverInterface::writeMpsNative(fullname.c_str(), 
3826                                     const_cast<const char **>(rowNames),
3827                                     const_cast<const char **>(columnNames),0,2,objSense,
3828                                     numberSOS_,setInfo_);
3829  if (rowNames) {
3830    modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_+1);
3831    modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
3832  }
3833}
3834
3835int 
3836OsiClpSolverInterface::writeMpsNative(const char *filename, 
3837                  const char ** rowNames, const char ** columnNames,
3838                  int formatType,int numberAcross,double objSense) const 
3839{
3840  return OsiSolverInterface::writeMpsNative(filename, rowNames, columnNames,
3841                               formatType, numberAcross,objSense,
3842                                            numberSOS_,setInfo_);
3843}
3844
3845//#############################################################################
3846// CLP specific public interfaces
3847//#############################################################################
3848
3849ClpSimplex * OsiClpSolverInterface::getModelPtr() const
3850{
3851  int saveAlgorithm = lastAlgorithm_;
3852  //freeCachedResults();
3853  lastAlgorithm_ = saveAlgorithm;
3854  //bool inCbcOrOther = (modelPtr_->specialOptions()&0x03000000)!=0;
3855  return modelPtr_;
3856}
3857
3858//-------------------------------------------------------------------
3859
3860//#############################################################################
3861// Constructors, destructors clone and assignment
3862//#############################################################################
3863//-------------------------------------------------------------------
3864// Default Constructor
3865//-------------------------------------------------------------------
3866OsiClpSolverInterface::OsiClpSolverInterface ()
3867:
3868OsiSolverInterface(),
3869rowsense_(NULL),
3870rhs_(NULL),
3871rowrange_(NULL),
3872ws_(NULL),
3873rowActivity_(NULL),
3874columnActivity_(NULL),
3875numberSOS_(0),
3876setInfo_(NULL),
3877smallModel_(NULL),
3878factorization_(NULL),
3879smallestElementInCut_(1.0e-15),
3880smallestChangeInCut_(1.0e-10),
3881largestAway_(-1.0),
3882spareArrays_(NULL),
3883matrixByRow_(NULL),
3884matrixByRowAtContinuous_(NULL), 
3885integerInformation_(NULL),
3886whichRange_(NULL),
3887fakeMinInSimplex_(false),
3888linearObjective_(NULL),
3889cleanupScaling_(0),
3890specialOptions_(0x80000000),
3891baseModel_(NULL),
3892lastNumberRows_(0),
3893continuousModel_(NULL),
3894fakeObjective_(NULL)
3895{
3896  //printf("in default %x\n",this);
3897  modelPtr_=NULL;
3898  notOwned_=false;
3899  disasterHandler_ = new OsiClpDisasterHandler();
3900  reset();
3901}
3902
3903//-------------------------------------------------------------------
3904// Clone
3905//-------------------------------------------------------------------
3906OsiSolverInterface * OsiClpSolverInterface::clone(bool CopyData) const
3907{
3908  //printf("in clone %x\n",this);
3909  OsiClpSolverInterface * newSolver;
3910  if (CopyData) {
3911    newSolver = new OsiClpSolverInterface(*this);
3912  } else {
3913    newSolver = new OsiClpSolverInterface();
3914  }
3915#if 0
3916  const double * obj = newSolver->getObjCoefficients();
3917  const double * oldObj = getObjCoefficients();
3918  if(newSolver->getNumCols()>3787)
3919    printf("%x - obj %x (from %x) val %g\n",newSolver,obj,oldObj,obj[3787]);
3920#endif
3921  return newSolver;
3922}
3923
3924
3925//-------------------------------------------------------------------
3926// Copy constructor
3927//-------------------------------------------------------------------
3928OsiClpSolverInterface::OsiClpSolverInterface (
3929                  const OsiClpSolverInterface & rhs)
3930: OsiSolverInterface(rhs),
3931rowsense_(NULL),
3932rhs_(NULL),
3933rowrange_(NULL),
3934ws_(NULL),
3935rowActivity_(NULL),
3936columnActivity_(NULL),
3937  stuff_(rhs.stuff_),
3938numberSOS_(rhs.numberSOS_),
3939setInfo_(NULL),
3940smallModel_(NULL),
3941factorization_(NULL),
3942smallestElementInCut_(rhs.smallestElementInCut_),
3943smallestChangeInCut_(rhs.smallestChangeInCut_),
3944largestAway_(-1.0),
3945spareArrays_(NULL),
3946basis_(),
3947itlimOrig_(9999999),
3948lastAlgorithm_(0),
3949notOwned_(false),
3950matrixByRow_(NULL),
3951matrixByRowAtContinuous_(NULL), 
3952integerInformation_(NULL),
3953whichRange_(NULL),
3954fakeMinInSimplex_(rhs.fakeMinInSimplex_)
3955{
3956  //printf("in copy %x - > %x\n",&rhs,this);
3957  if ( rhs.modelPtr_  ) 
3958    modelPtr_ = new ClpSimplex(*rhs.modelPtr_);
3959  else
3960    modelPtr_ = new ClpSimplex();
3961  if ( rhs.baseModel_  ) 
3962    baseModel_ = new ClpSimplex(*rhs.baseModel_);
3963  else
3964    baseModel_ = NULL;
3965  if ( rhs.continuousModel_  ) 
3966    continuousModel_ = new ClpSimplex(*rhs.continuousModel_);
3967  else
3968    continuousModel_ = NULL;
3969  if (rhs.matrixByRowAtContinuous_)
3970    matrixByRowAtContinuous_ = new CoinPackedMatrix(*rhs.matrixByRowAtContinuous_);
3971  if ( rhs.disasterHandler_  ) 
3972    disasterHandler_ = dynamic_cast<OsiClpDisasterHandler *>(rhs.disasterHandler_->clone());
3973  else
3974    disasterHandler_ = NULL;
3975  if (rhs.fakeObjective_)
3976    fakeObjective_ = new ClpLinearObjective(*rhs.fakeObjective_);
3977  else
3978    fakeObjective_ = NULL;
3979  linearObjective_ = modelPtr_->objective();
3980  if ( rhs.ws_ ) 
3981    ws_ = new CoinWarmStartBasis(*rhs.ws_);
3982  basis_ = rhs.basis_;
3983  if (rhs.integerInformation_) {
3984    int numberColumns = modelPtr_->numberColumns();
3985    integerInformation_ = new char[numberColumns];
3986    CoinMemcpyN(rhs.integerInformation_,        numberColumns,integerInformation_);
3987  }
3988  saveData_ = rhs.saveData_;
3989  solveOptions_  = rhs.solveOptions_;
3990  cleanupScaling_ = rhs.cleanupScaling_;
3991  specialOptions_ = rhs.specialOptions_;
3992  lastNumberRows_ = rhs.lastNumberRows_;
3993  rowScale_ = rhs.rowScale_;
3994  columnScale_ = rhs.columnScale_;
3995  fillParamMaps();
3996  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
3997  if (numberSOS_) {
3998    setInfo_ = new CoinSet[numberSOS_];
3999    for (int i=0;i<numberSOS_;i++)
4000      setInfo_[i]=rhs.setInfo_[i];
4001  }
4002}
4003
4004// Borrow constructor - only delete one copy
4005OsiClpSolverInterface::OsiClpSolverInterface (ClpSimplex * rhs,
4006                                              bool reallyOwn)
4007:
4008OsiSolverInterface(),
4009rowsense_(NULL),
4010rhs_(NULL),
4011rowrange_(NULL),
4012ws_(NULL),
4013rowActivity_(NULL),
4014columnActivity_(NULL),
4015numberSOS_(0),
4016setInfo_(NULL),
4017smallModel_(NULL),
4018factorization_(NULL),
4019smallestElementInCut_(1.0e-15),
4020smallestChangeInCut_(1.0e-10),
4021largestAway_(-1.0),
4022spareArrays_(NULL),
4023basis_(),
4024itlimOrig_(9999999),
4025lastAlgorithm_(0),
4026notOwned_(false),
4027matrixByRow_(NULL),
4028matrixByRowAtContinuous_(NULL), 
4029integerInformation_(NULL),
4030whichRange_(NULL),
4031fakeMinInSimplex_(false),
4032cleanupScaling_(0),
4033specialOptions_(0x80000000),
4034baseModel_(NULL),
4035lastNumberRows_(0),
4036continuousModel_(NULL),
4037fakeObjective_(NULL)
4038{
4039  disasterHandler_ = new OsiClpDisasterHandler();
4040  //printf("in borrow %x - > %x\n",&rhs,this);
4041  modelPtr_ = rhs;
4042  basis_.resize(modelPtr_->numberRows(),modelPtr_->numberColumns());
4043  linearObjective_ = modelPtr_->objective();
4044  if (rhs) {
4045    notOwned_=!reallyOwn;
4046
4047    if (rhs->integerInformation()) {
4048      int numberColumns = modelPtr_->numberColumns();
4049      integerInformation_ = new char[numberColumns];
4050      CoinMemcpyN(rhs->integerInformation(),    numberColumns,integerInformation_);
4051    }
4052  }
4053  fillParamMaps();
4054}
4055   
4056// Releases so won't error
4057void 
4058OsiClpSolverInterface::releaseClp()
4059{
4060  modelPtr_=NULL;
4061  notOwned_=false;
4062}
4063   
4064//-------------------------------------------------------------------
4065// Destructor
4066//-------------------------------------------------------------------
4067OsiClpSolverInterface::~OsiClpSolverInterface ()
4068{
4069  //printf("in destructor %x\n",this);
4070  freeCachedResults();
4071  if (!notOwned_)
4072    delete modelPtr_;
4073  delete baseModel_;
4074  delete continuousModel_;
4075  delete disasterHandler_;
4076  delete fakeObjective_;
4077  delete ws_;
4078  delete [] rowActivity_;
4079  delete [] columnActivity_;
4080  delete [] setInfo_;
4081#ifdef KEEP_SMALL
4082  if (smallModel_) {
4083    delete [] spareArrays_;
4084    spareArrays_ = NULL;
4085    delete smallModel_;
4086    smallModel_=NULL;
4087  }
4088#endif
4089  assert(smallModel_==NULL);
4090  assert(factorization_==NULL);
4091  assert(spareArrays_==NULL);
4092  delete [] integerInformation_;
4093  delete matrixByRowAtContinuous_;
4094  delete matrixByRow_;
4095}
4096
4097//-------------------------------------------------------------------
4098// Assignment operator
4099//-------------------------------------------------------------------
4100OsiClpSolverInterface &
4101OsiClpSolverInterface::operator=(const OsiClpSolverInterface& rhs)
4102{
4103  if (this != &rhs) {   
4104    //printf("in = %x - > %x\n",&rhs,this);
4105    OsiSolverInterface::operator=(rhs);
4106    freeCachedResults();
4107    if (!notOwned_)
4108      delete modelPtr_;
4109    delete ws_;
4110    if ( rhs.modelPtr_  ) 
4111      modelPtr_ = new ClpSimplex(*rhs.modelPtr_);
4112    delete baseModel_;
4113    if ( rhs.baseModel_  ) 
4114      baseModel_ = new ClpSimplex(*rhs.baseModel_);
4115    else
4116      baseModel_ = NULL;
4117    delete continuousModel_;
4118    if ( rhs.continuousModel_  ) 
4119      continuousModel_ = new ClpSimplex(*rhs.continuousModel_);
4120    else
4121      continuousModel_ = NULL;
4122    delete matrixByRowAtContinuous_;
4123    delete matrixByRow_;
4124    matrixByRow_=NULL;
4125    if (rhs.matrixByRowAtContinuous_)
4126      matrixByRowAtContinuous_ = new CoinPackedMatrix(*rhs.matrixByRowAtContinuous_);
4127    else
4128      matrixByRowAtContinuous_=NULL;
4129    delete disasterHandler_;
4130    if ( rhs.disasterHandler_  ) 
4131      disasterHandler_ = dynamic_cast<OsiClpDisasterHandler *>(rhs.disasterHandler_->clone());
4132    else
4133      disasterHandler_ = NULL;
4134    delete fakeObjective_;
4135    if (rhs.fakeObjective_)
4136      fakeObjective_ = new ClpLinearObjective(*rhs.fakeObjective_);
4137    else
4138      fakeObjective_ = NULL;
4139    notOwned_=false;
4140    linearObjective_ = modelPtr_->objective();
4141    saveData_ = rhs.saveData_;
4142    solveOptions_  = rhs.solveOptions_;
4143    cleanupScaling_ = rhs.cleanupScaling_;
4144    specialOptions_ = rhs.specialOptions_;
4145    lastNumberRows_ = rhs.lastNumberRows_;
4146    rowScale_ = rhs.rowScale_;
4147    columnScale_ = rhs.columnScale_;
4148    basis_ = rhs.basis_;
4149    stuff_ = rhs.stuff_;
4150    if (rhs.integerInformation_) {
4151      int numberColumns = modelPtr_->numberColumns();
4152      integerInformation_ = new char[numberColumns];
4153      CoinMemcpyN(rhs.integerInformation_,      numberColumns,integerInformation_);
4154    }
4155    if ( rhs.ws_ ) 
4156      ws_ = new CoinWarmStartBasis(*rhs.ws_);
4157    else
4158      ws_=NULL;
4159    delete [] rowActivity_;
4160    delete [] columnActivity_;
4161    rowActivity_=NULL;
4162    columnActivity_=NULL;
4163    delete [] setInfo_;
4164    numberSOS_ = rhs.numberSOS_;
4165    setInfo_=NULL;
4166    if (numberSOS_) {
4167      setInfo_ = new CoinSet[numberSOS_];
4168      for (int i=0;i<numberSOS_;i++)
4169        setInfo_[i]=rhs.setInfo_[i];
4170    }
4171    assert(smallModel_==NULL);
4172    assert(factorization_==NULL);
4173    smallestElementInCut_ = rhs.smallestElementInCut_;
4174    smallestChangeInCut_ = rhs.smallestChangeInCut_;
4175    largestAway_ = -1.0;
4176    assert(spareArrays_==NULL);
4177    basis_ = rhs.basis_;
4178    fillParamMaps();
4179    messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4180  }
4181  return *this;
4182}
4183
4184//#############################################################################
4185// Applying cuts
4186//#############################################################################
4187
4188void OsiClpSolverInterface::applyRowCut( const OsiRowCut & rowCut )
4189{
4190  applyRowCuts(1, &rowCut);
4191}
4192/* Apply a collection of row cuts which are all effective.
4193   applyCuts seems to do one at a time which seems inefficient.
4194*/
4195void 
4196OsiClpSolverInterface::applyRowCuts(int numberCuts, const OsiRowCut * cuts)
4197{
4198  if (numberCuts) {
4199    // Say can't gurantee optimal basis etc
4200    lastAlgorithm_=999;
4201
4202    // Thanks to js
4203    const OsiRowCut * * cutsp = new const OsiRowCut * [numberCuts];
4204    for (int i=0;i<numberCuts;i++) 
4205      cutsp[i] = &cuts[i];
4206   
4207    applyRowCuts(numberCuts, cutsp);
4208   
4209    delete [] cutsp;
4210  }
4211}
4212/* Apply a collection of row cuts which are all effective.
4213   applyCuts seems to do one at a time which seems inefficient.
4214*/
4215void 
4216OsiClpSolverInterface::applyRowCuts(int numberCuts, const OsiRowCut ** cuts)
4217{
4218  int i;
4219  if (!numberCuts)
4220    return;
4221  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
4222  CoinPackedMatrix * saveRowCopy = matrixByRow_;
4223  matrixByRow_=NULL;
4224#if 0 // was #ifndef NDEBUG
4225  int nameDiscipline;
4226  getIntParam(OsiNameDiscipline,nameDiscipline) ;
4227  assert (!nameDiscipline);
4228#endif
4229  freeCachedResults0();
4230  // Say can't gurantee optimal basis etc
4231  lastAlgorithm_=999;
4232  int numberRows = modelPtr_->numberRows();
4233  modelPtr_->resize(numberRows+numberCuts,modelPtr_->numberColumns());
4234  basis_.resize(numberRows+numberCuts,modelPtr_->numberColumns());
4235  // redo as relaxed - use modelPtr_-> addRows with starts etc
4236  int size = 0;
4237  for (i=0;i<numberCuts;i++) 
4238    size += cuts[i]->row().getNumElements();
4239  CoinBigIndex * starts = new CoinBigIndex [numberCuts+1];
4240  int * indices = new int[size];
4241  double * elements = new double[size];
4242  double * lower = modelPtr_->rowLower()+numberRows;
4243  double * upper = modelPtr_->rowUpper()+numberRows;
4244  const double * columnLower = modelPtr_->columnLower();
4245  const double * columnUpper = modelPtr_->columnUpper();
4246  size=0;
4247  for (i=0;i<numberCuts;i++) {
4248    double rowLb = cuts[i]->lb();
4249    double rowUb = cuts[i]->ub();
4250    int n=cuts[i]->row().getNumElements();
4251    const int * index = cuts[i]->row().getIndices();
4252    const double * elem = cuts[i]->row().getElements();
4253    starts[i]=size;
4254    for (int j=0;j<n;j++) {
4255      double value = elem[j];
4256      int column = index[j];
4257      if (fabs(value)>=smallestChangeInCut_) {
4258        // always take
4259        indices[size]=column;
4260        elements[size++]=value;
4261      } else if (fabs(value)>=smallestElementInCut_) {
4262        double lowerValue = columnLower[column];
4263        double upperValue = columnUpper[column];
4264        double difference = upperValue-lowerValue;
4265        if (difference<1.0e20&&difference*fabs(value)<smallestChangeInCut_&&
4266            (rowLb<-1.0e20||rowUb>1.0e20)) {
4267          // Take out and adjust to relax
4268          //printf("small el %g adjusted\n",value);
4269          if (rowLb>-1.0e20) {
4270            // just lower bound on row
4271            if (value>0.0) {
4272              // pretend at upper
4273              rowLb -= value*upperValue;
4274            } else {
4275              // pretend at lower
4276              rowLb -= value*lowerValue;
4277            }
4278          } else {
4279            // just upper bound on row
4280            if (value>0.0) {
4281              // pretend at lower
4282              rowUb -= value*lowerValue;
4283            } else {
4284              // pretend at upper
4285              rowUb -= value*upperValue;
4286            }
4287          }
4288        } else {
4289          // take (unwillingly)
4290          indices[size]=column;
4291          elements[size++]=value;
4292        }
4293      } else {
4294        //printf("small el %g ignored\n",value);
4295      }
4296    }
4297    lower[i]= forceIntoRange(rowLb, -OsiClpInfinity, OsiClpInfinity);
4298    upper[i]= forceIntoRange(rowUb, -OsiClpInfinity, OsiClpInfinity);
4299    if (lower[i]<-1.0e27)
4300      lower[i]=-COIN_DBL_MAX;
4301    if (upper[i]>1.0e27)
4302      upper[i]=COIN_DBL_MAX;
4303  }
4304  starts[numberCuts]=size;
4305 if (!modelPtr_->clpMatrix())
4306    modelPtr_->createEmptyMatrix();
4307  //modelPtr_->matrix()->appendRows(numberCuts,rows);
4308  modelPtr_->clpMatrix()->appendMatrix(numberCuts,0,starts,indices,elements);
4309  modelPtr_->setNewRowCopy(NULL);
4310  modelPtr_->setClpScaledMatrix(NULL);
4311  freeCachedResults1();
4312  redoScaleFactors( numberCuts,starts, indices, elements);
4313  if (saveRowCopy) {
4314#if 1
4315    matrixByRow_=saveRowCopy;
4316    matrixByRow_->appendRows(numberCuts,starts,indices,elements,0);
4317    if (matrixByRow_->getNumElements()!=modelPtr_->clpMatrix()->getNumElements()) {
4318      delete matrixByRow_; // odd type matrix
4319      matrixByRow_=NULL;
4320    }
4321#else
4322    delete saveRowCopy;
4323#endif
4324  }
4325  delete [] starts;
4326  delete [] indices;
4327  delete [] elements;
4328
4329}
4330//#############################################################################
4331// Apply Cuts
4332//#############################################################################
4333
4334OsiSolverInterface::ApplyCutsReturnCode
4335OsiClpSolverInterface::applyCuts( const OsiCuts & cs, double effectivenessLb ) 
4336{
4337  OsiSolverInterface::ApplyCutsReturnCode retVal;
4338  int i;
4339
4340  // Loop once for each column cut
4341  for ( i=0; i<cs.sizeColCuts(); i ++ ) {
4342    if ( cs.colCut(i).effectiveness() < effectivenessLb ) {
4343      retVal.incrementIneffective();
4344      continue;
4345    }
4346    if ( !cs.colCut(i).consistent() ) {
4347      retVal.incrementInternallyInconsistent();
4348      continue;
4349    }
4350    if ( !cs.colCut(i).consistent(*this) ) {
4351      retVal.incrementExternallyInconsistent();
4352      continue;
4353    }
4354    if ( cs.colCut(i).infeasible(*this) ) {
4355      retVal.incrementInfeasible();
4356      continue;
4357    }
4358    applyColCut( cs.colCut(i) );
4359    retVal.incrementApplied();
4360  }
4361
4362  // Loop once for each row cut
4363  const OsiRowCut ** addCuts = new const OsiRowCut* [cs.sizeRowCuts()];
4364  int nAdd=0;
4365  for ( i=0; i<cs.sizeRowCuts(); i ++ ) {
4366    if ( cs.rowCut(i).effectiveness() < effectivenessLb ) {
4367      retVal.incrementIneffective();
4368      continue;
4369    }
4370    if ( !cs.rowCut(i).consistent() ) {
4371      retVal.incrementInternallyInconsistent();
4372      continue;
4373    }
4374    if ( !cs.rowCut(i).consistent(*this) ) {
4375      retVal.incrementExternallyInconsistent();
4376      continue;
4377    }
4378    if ( cs.rowCut(i).infeasible(*this) ) {
4379      retVal.incrementInfeasible();
4380      continue;
4381    }
4382    addCuts[nAdd++] = cs.rowCutPtr(i);
4383    retVal.incrementApplied();
4384  }
4385  // now apply
4386  applyRowCuts(nAdd,addCuts);
4387  delete [] addCuts;
4388 
4389  return retVal;
4390}
4391// Extend scale factors
4392void 
4393OsiClpSolverInterface::redoScaleFactors(int numberAdd,const CoinBigIndex * starts,
4394                                        const int * indices, const double * elements)
4395{
4396  if ((specialOptions_&131072)!=0) {
4397    int numberRows = modelPtr_->numberRows()-numberAdd;
4398    assert (lastNumberRows_==numberRows); // ???
4399    int iRow;
4400    int newNumberRows = numberRows + numberAdd;
4401    rowScale_.extend(static_cast<int>(2*newNumberRows*sizeof(double)));
4402    double * rowScale = rowScale_.array();
4403    double * oldInverseScale = rowScale + lastNumberRows_;
4404    double * inverseRowScale = rowScale + newNumberRows;
4405    for (iRow=lastNumberRows_-1;iRow>=0;iRow--)
4406      inverseRowScale[iRow] = oldInverseScale[iRow] ;
4407    //int numberColumns = baseModel_->numberColumns();
4408    const double * columnScale = columnScale_.array();
4409    //const double * inverseColumnScale = columnScale + numberColumns;
4410    // Geometric mean on row scales
4411    // adjust arrays
4412    rowScale += lastNumberRows_;
4413    inverseRowScale += lastNumberRows_;
4414    for (iRow=0;iRow<numberAdd;iRow++) {
4415      CoinBigIndex j;
4416      double largest=1.0e-20;
4417      double smallest=1.0e50;
4418      for (j=starts[iRow];j<starts[iRow+1];j++) {
4419        int iColumn = indices[j];
4420        double value = fabs(elements[j]);
4421        // Don't bother with tiny elements
4422        if (value>1.0e-20) {
4423          value *= columnScale[iColumn];
4424          largest = CoinMax(largest,value);
4425          smallest = CoinMin(smallest,value);
4426        }
4427      }
4428      double scale=sqrt(smallest*largest);
4429      scale=CoinMax(1.0e-10,CoinMin(1.0e10,scale));
4430      inverseRowScale[iRow]=scale;
4431      rowScale[iRow]=1.0/scale;
4432    }
4433    lastNumberRows_=newNumberRows;
4434  }
4435}
4436// Delete all scale factor stuff and reset option
4437void OsiClpSolverInterface::deleteScaleFactors()
4438{
4439  delete baseModel_;
4440  baseModel_=NULL;
4441  lastNumberRows_=0;
4442  specialOptions_ &= ~131072;
4443}
4444//-----------------------------------------------------------------------------
4445
4446void OsiClpSolverInterface::applyColCut( const OsiColCut & cc )
4447{
4448  modelPtr_->whatsChanged_ &= (0x1ffff&~(128|256));
4449  // Say can't gurantee optimal basis etc
4450  lastAlgorithm_=999;
4451  double * lower = modelPtr_->columnLower();
4452  double * upper = modelPtr_->columnUpper();
4453  const CoinPackedVector & lbs = cc.lbs();
4454  const CoinPackedVector & ubs = cc.ubs();
4455  int i;
4456
4457  for ( i=0; i<lbs.getNumElements(); i++ ) {
4458    int iCol = lbs.getIndices()[i];
4459    double value = lbs.getElements()[i];
4460    if ( value > lower[iCol] )
4461      lower[iCol]= value;
4462  }
4463  for ( i=0; i<ubs.getNumElements(); i++ ) {
4464    int iCol = ubs.getIndices()[i];
4465    double value = ubs.getElements()[i];
4466    if ( value < upper[iCol] )
4467      upper[iCol]= value;
4468  }
4469}
4470//#############################################################################
4471// Private methods
4472//#############################################################################
4473
4474
4475//-------------------------------------------------------------------
4476
4477void OsiClpSolverInterface::freeCachedResults() const
4478{ 
4479  // Say can't gurantee optimal basis etc
4480  lastAlgorithm_=999;
4481  delete [] rowsense_;
4482  delete [] rhs_;
4483  delete [] rowrange_;
4484  delete matrixByRow_;
4485  //delete ws_;
4486  rowsense_=NULL;
4487  rhs_=NULL;
4488  rowrange_=NULL;
4489  matrixByRow_=NULL;
4490  //ws_ = NULL;
4491  if (!notOwned_&&modelPtr_) {
4492    if(modelPtr_->scaledMatrix_) {
4493      delete modelPtr_->scaledMatrix_;
4494      modelPtr_->scaledMatrix_=NULL;
4495    }
4496    if (modelPtr_->clpMatrix()) {
4497      modelPtr_->clpMatrix()->refresh(modelPtr_); // make sure all clean
4498#ifndef NDEBUG
4499      ClpPackedMatrix * clpMatrix = dynamic_cast<ClpPackedMatrix *> (modelPtr_->clpMatrix());
4500      if (clpMatrix) {
4501        if (clpMatrix->getNumCols())
4502          assert (clpMatrix->getNumRows()==modelPtr_->getNumRows());
4503        if (clpMatrix->getNumRows())
4504          assert (clpMatrix->getNumCols()==modelPtr_->getNumCols());
4505      }
4506#endif
4507    }
4508  }
4509}
4510
4511//-------------------------------------------------------------------
4512
4513void OsiClpSolverInterface::freeCachedResults0() const
4514{ 
4515  delete [] rowsense_;
4516  delete [] rhs_;
4517  delete [] rowrange_;
4518  rowsense_=NULL;
4519  rhs_=NULL;
4520  rowrange_=NULL;
4521}
4522
4523//-------------------------------------------------------------------
4524
4525void OsiClpSolverInterface::freeCachedResults1() const
4526{ 
4527  // Say can't gurantee optimal basis etc
4528  lastAlgorithm_=999;
4529  delete matrixByRow_;
4530  matrixByRow_=NULL;
4531  //ws_ = NULL;
4532  if (modelPtr_&&modelPtr_->clpMatrix()) {
4533    delete modelPtr_->scaledMatrix_;
4534    modelPtr_->scaledMatrix_=NULL;
4535    modelPtr_->clpMatrix()->refresh(modelPtr_); // make sure all clean
4536#ifndef NDEBUG
4537    ClpPackedMatrix * clpMatrix = dynamic_cast<ClpPackedMatrix *> (modelPtr_->clpMatrix());
4538    if (clpMatrix) {
4539      assert (clpMatrix->getNumRows()==modelPtr_->getNumRows());
4540      assert (clpMatrix->getNumCols()==modelPtr_->getNumCols());
4541    }
4542#endif
4543  }
4544}
4545
4546//------------------------------------------------------------------
4547void OsiClpSolverInterface::extractSenseRhsRange() const
4548{
4549  if (rowsense_ == NULL) {
4550    // all three must be NULL
4551    assert ((rhs_ == NULL) && (rowrange_ == NULL));
4552   
4553    int nr=modelPtr_->numberRows();
4554    if ( nr!=0 ) {
4555      rowsense_ = new char[nr];
4556      rhs_ = new double[nr];
4557      rowrange_ = new double[nr];
4558      std::fill(rowrange_,rowrange_+nr,0.0);
4559     
4560      const double * lb = modelPtr_->rowLower();
4561      const double * ub = modelPtr_->rowUpper();
4562     
4563      int i;
4564      for ( i=0; i<nr; i++ ) {
4565        convertBoundToSense(lb[i], ub[i], rowsense_[i], rhs_[i], rowrange_[i]);
4566      }
4567    }
4568  }
4569}
4570// Set language
4571void 
4572OsiClpSolverInterface::newLanguage(CoinMessages::Language language)
4573{
4574  modelPtr_->newLanguage(language);
4575  OsiSolverInterface::newLanguage(language);
4576}
4577//#############################################################################
4578
4579void
4580OsiClpSolverInterface::fillParamMaps()
4581{
4582  assert (static_cast<int> (OsiMaxNumIteration)==        static_cast<int>(ClpMaxNumIteration));
4583  assert (static_cast<int> (OsiMaxNumIterationHotStart)==static_cast<int>(ClpMaxNumIterationHotStart));
4584  //assert (static_cast<int> (OsiLastIntParam)==           static_cast<int>(ClpLastIntParam));
4585 
4586  assert (static_cast<int> (OsiDualObjectiveLimit)==  static_cast<int>(ClpDualObjectiveLimit));
4587  assert (static_cast<int> (OsiPrimalObjectiveLimit)==static_cast<int>(ClpPrimalObjectiveLimit));
4588  assert (static_cast<int> (OsiDualTolerance)==       static_cast<int>(ClpDualTolerance));
4589  assert (static_cast<int> (OsiPrimalTolerance)==     static_cast<int>(ClpPrimalTolerance));
4590  assert (static_cast<int> (OsiObjOffset)==           static_cast<int>(ClpObjOffset));
4591  //assert (static_cast<int> (OsiLastDblParam)==        static_cast<int>(ClpLastDblParam));
4592 
4593  assert (static_cast<int> (OsiProbName)==    static_cast<int> (ClpProbName));
4594  //strParamMap_[OsiLastStrParam] = ClpLastStrParam;
4595}
4596// Sets up basis
4597void 
4598OsiClpSolverInterface::setBasis ( const CoinWarmStartBasis & basis)
4599{
4600  setBasis(basis,modelPtr_);
4601  setWarmStart(&basis); 
4602}
4603// Warm start
4604CoinWarmStartBasis
4605OsiClpSolverInterface::getBasis(ClpSimplex * model) const
4606{
4607  int iRow,iColumn;
4608  int numberRows = model->numberRows();
4609  int numberColumns = model->numberColumns();
4610  CoinWarmStartBasis basis;
4611  basis.setSize(numberColumns,numberRows);
4612  if (model->statusExists()) {
4613    // Flip slacks
4614    int lookupA[]={0,1,3,2,0,2};
4615    for (iRow=0;iRow<numberRows;iRow++) {
4616      int iStatus = model->getRowStatus(iRow);
4617      iStatus = lookupA[iStatus];
4618      basis.setArtifStatus(iRow,static_cast<CoinWarmStartBasis::Status> (iStatus));
4619    }
4620    int lookupS[]={0,1,2,3,0,3};
4621    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4622      int iStatus = model->getColumnStatus(iColumn);
4623      iStatus = lookupS[iStatus];
4624      basis.setStructStatus(iColumn,static_cast<CoinWarmStartBasis::Status> (iStatus));
4625    }
4626  }
4627  //basis.print();
4628  return basis;
4629}
4630// Warm start from statusArray
4631CoinWarmStartBasis * 
4632OsiClpSolverInterface::getBasis(const unsigned char * statusArray) const 
4633{
4634  int iRow,iColumn;
4635  int numberRows = modelPtr_->numberRows();
4636  int numberColumns = modelPtr_->numberColumns();
4637  CoinWarmStartBasis * basis = new CoinWarmStartBasis();
4638  basis->setSize(numberColumns,numberRows);
4639  // Flip slacks
4640  int lookupA[]={0,1,3,2,0,2};
4641  for (iRow=0;iRow<numberRows;iRow++) {
4642    int iStatus = statusArray[numberColumns+iRow]&7;
4643    iStatus = lookupA[iStatus];
4644    basis->setArtifStatus(iRow,static_cast<CoinWarmStartBasis::Status> (iStatus));
4645  }
4646  int lookupS[]={0,1,2,3,0,3};
4647  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4648    int iStatus = statusArray[iColumn]&7;
4649    iStatus = lookupS[iStatus];
4650    basis->setStructStatus(iColumn,static_cast<CoinWarmStartBasis::Status> (iStatus));
4651  }
4652  //basis->print();
4653  return basis;
4654}
4655// Sets up basis
4656void 
4657OsiClpSolverInterface::setBasis ( const CoinWarmStartBasis & basis,
4658                                  ClpSimplex * model)
4659{
4660  // Say can't gurantee optimal basis etc
4661  lastAlgorithm_=999;
4662  // transform basis to status arrays
4663  int iRow,iColumn;
4664  int numberRows = model->numberRows();
4665  int numberColumns = model->numberColumns();
4666  if (!model->statusExists()) {
4667    /*
4668      get status arrays
4669      ClpBasis would seem to have overheads and we will need
4670      extra bits anyway.
4671    */
4672    model->createStatus();
4673  }
4674  if (basis.getNumArtificial()!=numberRows||
4675      basis.getNumStructural()!=numberColumns) {
4676    CoinWarmStartBasis basis2 = basis;
4677    // resize
4678    basis2.resize(numberRows,numberColumns);
4679    // move status
4680    model->createStatus();
4681    // For rows lower and upper are flipped
4682    for (iRow=0;iRow<numberRows;iRow++) {
4683      int stat = basis2.getArtifStatus(iRow);
4684      if (stat>1)
4685        stat = 5 - stat; // so 2->3 and 3->2
4686      model->setRowStatus(iRow, static_cast<ClpSimplex::Status> (stat));
4687    }
4688    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4689      model->setColumnStatus(iColumn,
4690                             static_cast<ClpSimplex::Status> (basis2.getStructStatus(iColumn)));
4691    }
4692  } else {
4693    // move status
4694    model->createStatus();
4695    // For rows lower and upper are flipped
4696    for (iRow=0;iRow<numberRows;iRow++) {
4697      int stat = basis.getArtifStatus(iRow);
4698      if (stat>1)
4699        stat = 5 - stat; // so 2->3 and 3->2
4700      model->setRowStatus(iRow, static_cast<ClpSimplex::Status> (stat));
4701    }
4702    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4703      model->setColumnStatus(iColumn,
4704                             static_cast<ClpSimplex::Status> (basis.getStructStatus(iColumn)));
4705    }
4706  }
4707}
4708// Warm start difference from basis_ to statusArray
4709CoinWarmStartDiff * 
4710OsiClpSolverInterface::getBasisDiff(const unsigned char * statusArray) const
4711{
4712  int iRow,iColumn;
4713  int numberRows = modelPtr_->numberRows();
4714  int numberColumns = modelPtr_->numberColumns();
4715  CoinWarmStartBasis basis;
4716  basis.setSize(numberColumns,numberRows);
4717  assert (modelPtr_->statusExists());
4718  int lookupS[]={0,1,2,3,0,3};
4719  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4720    int iStatus = statusArray[iColumn]&7;
4721    iStatus = lookupS[iStatus];
4722    basis.setStructStatus(iColumn,static_cast<CoinWarmStartBasis::Status> (iStatus));
4723  }
4724  statusArray += numberColumns;
4725  // Flip slacks
4726  int lookupA[]={0,1,3,2,0,2};
4727  for (iRow=0;iRow<numberRows;iRow++) {
4728    int iStatus = statusArray[iRow]&7;
4729    iStatus = lookupA[iStatus];
4730    basis.setArtifStatus(iRow,static_cast<CoinWarmStartBasis::Status> (iStatus));
4731  }
4732  // Now basis is what we want while basis_ is old
4733  CoinWarmStartDiff * difference = basis.generateDiff(&basis_);
4734  return difference;
4735}
4736/*
4737  Read an mps file from the given filename - returns number of errors
4738  (see CoinMpsIO class)
4739*/
4740int 
4741OsiClpSolverInterface::readMps(const char *filename,
4742                               const char *extension ) 
4743{
4744  // Get rid of integer stuff
4745  delete [] integerInformation_;
4746  integerInformation_=NULL;
4747  freeCachedResults();
4748 
4749  CoinMpsIO m;
4750  m.setInfinity(getInfinity());
4751  m.passInMessageHandler(modelPtr_->messageHandler());
4752  *m.messagesPointer()=modelPtr_->coinMessages();
4753
4754  delete [] setInfo_;
4755  setInfo_=NULL;
4756  numberSOS_=0;
4757  CoinSet ** sets=NULL;
4758  // Temporarily reduce log level to get CoinMpsIO to shut up.
4759  int saveLogLevel = modelPtr_->messageHandler()->logLevel() ;
4760  modelPtr_->messageHandler()->setLogLevel(0) ;
4761  int numberErrors = m.readMps(filename,extension,numberSOS_,sets);
4762  modelPtr_->messageHandler()->setLogLevel(saveLogLevel) ;
4763  if (numberSOS_) {
4764    setInfo_ = new CoinSet[numberSOS_];
4765    for (int i=0;i<numberSOS_;i++) {
4766      setInfo_[i]=*sets[i];
4767      delete sets[i];
4768    }
4769    delete [] sets;
4770  }
4771  handler_->message(COIN_SOLVER_MPS,messages_)
4772    <<m.getProblemName()<< numberErrors <<CoinMessageEol;
4773  if (!numberErrors) {
4774
4775    // set objective function offest
4776    setDblParam(OsiObjOffset,m.objectiveOffset());
4777
4778    // set problem name
4779    setStrParam(OsiProbName,m.getProblemName());
4780   
4781    // no errors
4782    loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
4783                m.getObjCoefficients(),m.getRowSense(),m.getRightHandSide(),
4784                m.getRowRange());
4785    const char * integer = m.integerColumns();
4786    int nCols=m.getNumCols();
4787    int nRows=m.getNumRows();
4788    if (integer) {
4789      int i,n=0;
4790      int * index = new int [nCols];
4791      for (i=0;i<nCols;i++) {
4792        if (integer[i]) {
4793          index[n++]=i;
4794        }
4795      }
4796      setInteger(index,n);
4797      delete [] index;
4798      if (n) 
4799        modelPtr_->copyInIntegerInformation(integer);
4800    }
4801
4802    // set objective name
4803    setObjName(m.getObjectiveName());
4804
4805    // Always keep names
4806    int nameDiscipline;
4807    getIntParam(OsiNameDiscipline,nameDiscipline) ;
4808    int iRow;
4809    std::vector<std::string> rowNames = std::vector<std::string> ();
4810    std::vector<std::string> columnNames = std::vector<std::string> ();
4811    rowNames.reserve(nRows);
4812    for (iRow=0;iRow<nRows;iRow++) {
4813      const char * name = m.rowName(iRow);
4814      rowNames.push_back(name);
4815      if (nameDiscipline) 
4816        OsiSolverInterface::setRowName(iRow,name) ;
4817    }
4818   
4819    int iColumn;
4820    columnNames.reserve(nCols);
4821    for (iColumn=0;iColumn<nCols;iColumn++) {
4822      const char * name = m.columnName(iColumn);
4823      columnNames.push_back(name);
4824      if (nameDiscipline) 
4825        OsiSolverInterface::setColName(iColumn,name) ;
4826    }
4827    modelPtr_->copyNames(rowNames,columnNames);
4828  }
4829  return numberErrors;
4830}
4831int 
4832OsiClpSolverInterface::readMps(const char *filename, const char*extension,
4833                            int & numberSets, CoinSet ** & sets)
4834{
4835  int numberErrors = readMps(filename,extension);
4836  numberSets= numberSOS_;
4837  sets = &setInfo_;
4838  return numberErrors;
4839}
4840/* Read an mps file from the given filename returns
4841   number of errors (see OsiMpsReader class) */
4842int 
4843OsiClpSolverInterface::readMps(const char *filename,bool keepNames,bool allowErrors)
4844{
4845  // Get rid of integer stuff
4846  delete [] integerInformation_;
4847  integerInformation_=NULL;
4848  freeCachedResults();
4849 
4850  CoinMpsIO m;
4851  m.setInfinity(getInfinity());
4852  m.passInMessageHandler(modelPtr_->messageHandler());
4853  *m.messagesPointer()=modelPtr_->coinMessages();
4854  m.setSmallElementValue(CoinMax(modelPtr_->getSmallElementValue(), 
4855                                 m.getSmallElementValue()));
4856
4857  delete [] setInfo_;
4858  setInfo_=NULL;
4859  numberSOS_=0;
4860  CoinSet ** sets=NULL;
4861  int numberErrors = m.readMps(filename,"",numberSOS_,sets);
4862  if (numberSOS_) {
4863    setInfo_ = new CoinSet[numberSOS_];
4864    for (int i=0;i<numberSOS_;i++) {
4865      setInfo_[i]=*sets[i];
4866      delete sets[i];
4867    }
4868    delete [] sets;
4869  }
4870  handler_->message(COIN_SOLVER_MPS,messages_)
4871    <<m.getProblemName()<< numberErrors <<CoinMessageEol;
4872  if (!numberErrors||((numberErrors>0&&numberErrors<100000)&&allowErrors)) {
4873
4874    // set objective function offest
4875    setDblParam(OsiObjOffset,m.objectiveOffset());
4876
4877    // set problem name
4878    setStrParam(OsiProbName,m.getProblemName());
4879
4880    // set objective name
4881    setObjName(m.getObjectiveName());
4882
4883    // no errors
4884    loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
4885                m.getObjCoefficients(),m.getRowSense(),m.getRightHandSide(),
4886                m.getRowRange());
4887    int nCols=m.getNumCols();
4888    // get quadratic part
4889    if (m.reader()->whichSection (  ) == COIN_QUAD_SECTION ) {
4890      int * start=NULL;
4891      int * column = NULL;
4892      double * element = NULL;
4893      int status=m.readQuadraticMps(NULL,start,column,element,2);
4894      if (!status) 
4895        modelPtr_->loadQuadraticObjective(nCols,start,column,element);
4896      delete [] start;
4897      delete [] column;
4898      delete [] element;
4899    }
4900    const char * integer = m.integerColumns();
4901    int nRows=m.getNumRows();
4902    if (integer) {
4903      int i,n=0;
4904      int * index = new int [nCols];
4905      for (i=0;i<nCols;i++) {
4906        if (integer[i]) {
4907          index[n++]=i;
4908        }
4909      }
4910      setInteger(index,n);
4911      delete [] index;
4912      if (n) 
4913        modelPtr_->copyInIntegerInformation(integer);
4914    }
4915    if (keepNames) {
4916      // keep names
4917      int nameDiscipline;
4918      getIntParam(OsiNameDiscipline,nameDiscipline) ;
4919      int iRow;
4920      std::vector<std::string> rowNames = std::vector<std::string> ();
4921      std::vector<std::string> columnNames = std::vector<std::string> ();
4922      rowNames.reserve(nRows);
4923      for (iRow=0;iRow<nRows;iRow++) {
4924        const char * name = m.rowName(iRow);
4925        rowNames.push_back(name);
4926        if (nameDiscipline) 
4927          OsiSolverInterface::setRowName(iRow,name) ;
4928      }
4929     
4930      int iColumn;
4931      columnNames.reserve(nCols);
4932      for (iColumn=0;iColumn<nCols;iColumn++) {
4933        const char * name = m.columnName(iColumn);
4934        columnNames.push_back(name);
4935        if (nameDiscipline) 
4936          OsiSolverInterface::setColName(iColumn,name) ;
4937      }
4938      modelPtr_->copyNames(rowNames,columnNames);
4939    }
4940  }
4941  return numberErrors;
4942}
4943// Read file in LP format (with names)
4944int 
4945OsiClpSolverInterface::readLp(const char *filename, const double epsilon )
4946{
4947  CoinLpIO m;
4948  m.passInMessageHandler(modelPtr_->messageHandler());
4949  *m.messagesPointer()=modelPtr_->coinMessages();
4950  m.readLp(filename, epsilon);
4951  freeCachedResults();
4952
4953  // set objective function offest
4954  setDblParam(OsiObjOffset, 0);
4955
4956  // set problem name
4957  setStrParam(OsiProbName, m.getProblemName());
4958
4959  // set objective name
4960  setObjName(m.getObjName());
4961
4962  // no errors
4963  loadProblem(*m.getMatrixByRow(), m.getColLower(), m.getColUpper(),
4964              m.getObjCoefficients(), m.getRowLower(), m.getRowUpper());
4965
4966  const char *integer = m.integerColumns();
4967  int nCols = m.getNumCols();
4968  int nRows = m.getNumRows();
4969  if (integer) {
4970    int i, n = 0;
4971    int *index = new int [nCols];
4972    for (i=0; i<nCols; i++) {
4973      if (integer[i]) {
4974        index[n++] = i;
4975      }
4976    }
4977    setInteger(index,n);
4978    delete [] index;
4979  }
4980  // Always keep names
4981  int nameDiscipline;
4982  getIntParam(OsiNameDiscipline,nameDiscipline) ;
4983  int iRow;
4984  std::vector<std::string> rowNames = std::vector<std::string> ();
4985  std::vector<std::string> columnNames = std::vector<std::string> ();
4986  rowNames.reserve(nRows);
4987  for (iRow=0;iRow<nRows;iRow++) {
4988    const char * name = m.rowName(iRow);
4989    rowNames.push_back(name);
4990    if (nameDiscipline) 
4991      OsiSolverInterface::setRowName(iRow,name) ;
4992  }
4993 
4994  int iColumn;
4995  columnNames.reserve(nCols);
4996  for (iColumn=0;iColumn<nCols;iColumn++) {
4997    const char * name = m.columnName(iColumn);
4998    columnNames.push_back(name);
4999    if (nameDiscipline) 
5000      OsiSolverInterface::setColName(iColumn,name) ;
5001  }
5002  modelPtr_->copyNames(rowNames,columnNames);
5003  return(0);
5004}
5005/* Write the problem into an Lp file of the given filename.
5006   If objSense is non zero then -1.0 forces the code to write a
5007   maximization objective and +1.0 to write a minimization one.
5008   If 0.0 then solver can do what it wants.
5009   This version calls writeLpNative with names */
5010void 
5011OsiClpSolverInterface::writeLp(const char *filename,
5012                               const char *extension ,
5013                               double epsilon ,
5014                               int numberAcross ,
5015                               int decimals ,
5016                               double objSense ,
5017                               bool changeNameOnRange) const
5018{
5019  std::string f(filename);
5020  std::string e(extension);
5021  std::string fullname;
5022  if (e!="") {
5023    fullname = f + "." + e;
5024  } else {
5025    // no extension so no trailing period
5026    fullname = f;
5027  }
5028  // get names
5029  const char * const * const rowNames = modelPtr_->rowNamesAsChar();
5030  const char * const * const columnNames = modelPtr_->columnNamesAsChar();
5031  // Fall back on Osi version - possibly with names
5032  OsiSolverInterface::writeLpNative(fullname.c_str(), 
5033                                    rowNames,columnNames, epsilon, numberAcross,
5034                                    decimals, objSense,changeNameOnRange);
5035  if (rowNames) {
5036    modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_+1);
5037    modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
5038  }
5039}
5040void 
5041OsiClpSolverInterface::writeLp(FILE * fp,
5042                               double epsilon ,
5043                               int numberAcross ,
5044                               int decimals ,
5045                               double objSense ,
5046                               bool changeNameOnRange) const
5047{
5048  // get names
5049  const char * const * const rowNames = modelPtr_->rowNamesAsChar();
5050  const char * const * const columnNames = modelPtr_->columnNamesAsChar();
5051  // Fall back on Osi version - possibly with names
5052  OsiSolverInterface::writeLpNative(fp,
5053                                    rowNames,columnNames, epsilon, numberAcross,
5054                                    decimals, objSense,changeNameOnRange);
5055  if (rowNames) {
5056    modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_+1);
5057    modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
5058  }
5059}
5060/*
5061  I (JJF) am getting incredibly annoyed because I can't just replace a matrix.
5062  The default behavior of this is do nothing so only use where that would not matter
5063  e.g. strengthening a matrix for MIP
5064*/
5065void 
5066OsiClpSolverInterface::replaceMatrixOptional(const CoinPackedMatrix & matrix)
5067{
5068  modelPtr_->whatsChanged_ &= (0xffff&~(2|4|8));
5069  replaceMatrix(matrix);
5070}
5071// And if it does matter (not used at present)
5072void 
5073OsiClpSolverInterface::replaceMatrix(const CoinPackedMatrix & matrix)
5074{
5075  modelPtr_->whatsChanged_ &= (0xffff&~(2|4|8));
5076  delete modelPtr_->matrix_;
5077  delete modelPtr_->rowCopy_;
5078  modelPtr_->rowCopy_=NULL;
5079  if (matrix.isColOrdered()) {
5080    modelPtr_->matrix_=new ClpPackedMatrix(matrix);
5081  } else {
5082    CoinPackedMatrix matrix2;
5083    matrix2.setExtraGap(0.0);
5084    matrix2.setExtraMajor(0.0);
5085    matrix2.reverseOrderedCopyOf(matrix);
5086    modelPtr_->matrix_=new ClpPackedMatrix(matrix2);
5087  }   
5088  modelPtr_->matrix_->setDimensions(modelPtr_->numberRows_,modelPtr_->numberColumns_);
5089  freeCachedResults();
5090}
5091// Get pointer to array[getNumCols()] of primal solution vector
5092const double * 
5093OsiClpSolverInterface::getColSolution() const 
5094{ 
5095  if (modelPtr_->solveType()!=2) {
5096    return modelPtr_->primalColumnSolution();
5097  } else {
5098    // simplex interface
5099    return modelPtr_->solutionRegion(1);
5100  }
5101}
5102 
5103// Get pointer to array[getNumRows()] of dual prices
5104const double * 
5105OsiClpSolverInterface::getRowPrice() const
5106{ 
5107  if (modelPtr_->solveType()!=2) {
5108    return modelPtr_->dualRowSolution();
5109  } else {
5110    // simplex interface
5111    //return modelPtr_->djRegion(0);
5112    return modelPtr_->dualRowSolution();
5113  }
5114}
5115 
5116// Get a pointer to array[getNumCols()] of reduced costs
5117const double * 
5118OsiClpSolverInterface::getReducedCost() const 
5119{ 
5120  if (modelPtr_->solveType()!=2) {
5121    return modelPtr_->dualColumnSolution();
5122  } else {
5123    // simplex interface
5124    return modelPtr_->djRegion(1);
5125  }
5126}
5127
5128/* Get pointer to array[getNumRows()] of row activity levels (constraint
5129   matrix times the solution vector */
5130const double * 
5131OsiClpSolverInterface::getRowActivity() const 
5132{ 
5133  if (modelPtr_->solveType()!=2) {
5134    return modelPtr_->primalRowSolution();
5135  } else {
5136    // simplex interface
5137    return modelPtr_->solutionRegion(0);
5138  }
5139}
5140double 
5141OsiClpSolverInterface::getObjValue() const 
5142{
5143  if (modelPtr_->numberIterations()||modelPtr_->upperIn_!=-COIN_DBL_MAX) {
5144    // This does not pass unitTest when getObjValue is called before solve.
5145    //printf("obj a %g %g\n",modelPtr_->objectiveValue(),
5146    //     OsiSolverInterface::getObjValue());
5147    if (fakeMinInSimplex_)
5148      return -modelPtr_->objectiveValue() ;
5149    else
5150      return modelPtr_->objectiveValue();
5151  } else {
5152    return OsiSolverInterface::getObjValue();
5153  }
5154}
5155
5156/* Set an objective function coefficient */
5157void 
5158OsiClpSolverInterface::setObjCoeff( int elementIndex, double elementValue )
5159{
5160  modelPtr_->whatsChanged_ &= 0xffff;
5161  // Say can't gurantee optimal basis etc
5162  lastAlgorithm_=999;
5163#ifndef NDEBUG
5164  int n = modelPtr_->numberColumns();
5165  if (elementIndex<0||elementIndex>=n) {
5166    indexError(elementIndex,"setObjCoeff");
5167  }
5168#endif
5169  modelPtr_->setObjectiveCoefficient(elementIndex,
5170    ((fakeMinInSimplex_)?-elementValue:elementValue));
5171}
5172
5173/* Set a single column lower bound<br>
5174   Use -DBL_MAX for -infinity. */
5175void 
5176OsiClpSolverInterface::setColLower( int index, double elementValue )
5177{
5178  modelPtr_->whatsChanged_ &= 0x1ffff;
5179#ifndef NDEBUG
5180  int n = modelPtr_->numberColumns();
5181  if (index<0||index>=n) {
5182    indexError(index,"setColLower");
5183  }
5184#endif
5185  double currentValue = modelPtr_->columnActivity_[index];
5186  bool changed=(currentValue<elementValue-modelPtr_->primalTolerance()||
5187                index>=basis_.getNumStructural()||
5188                basis_.getStructStatus(index)==CoinWarmStartBasis::atLowerBound);
5189  // Say can't gurantee optimal basis etc
5190  if (changed)
5191    lastAlgorithm_=999;
5192  if (!modelPtr_->lower_)
5193    modelPtr_->whatsChanged_ &= ~0xffff; // switch off
5194  modelPtr_->setColumnLower(index,elementValue);
5195}
5196     
5197/* Set a single column upper bound<br>
5198   Use DBL_MAX for infinity. */
5199void 
5200OsiClpSolverInterface::setColUpper( int index, double elementValue )
5201{
5202  modelPtr_->whatsChanged_ &= 0x1ffff;
5203#ifndef NDEBUG
5204  int n = modelPtr_->numberColumns();
5205  if (index<0||index>=n) {
5206    indexError(index,"setColUpper");
5207  }
5208#endif
5209  double currentValue = modelPtr_->columnActivity_[index];
5210  bool changed=(currentValue>elementValue+modelPtr_->primalTolerance()||
5211                index>=basis_.getNumStructural()||
5212                basis_.getStructStatus(index)==CoinWarmStartBasis::atUpperBound);
5213  // Say can't gurantee optimal basis etc
5214  if (changed)
5215    lastAlgorithm_=999;
5216  if (!modelPtr_->upper_)
5217    modelPtr_->whatsChanged_ &= ~0xffff; // switch off
5218  modelPtr_->setColumnUpper(index,elementValue);
5219}
5220
5221/* Set a single column lower and upper bound */
5222void 
5223OsiClpSolverInterface::setColBounds( int elementIndex,
5224                                     double lower, double upper )
5225{
5226  modelPtr_->whatsChanged_ &= 0x1ffff;
5227  // Say can't gurantee optimal basis etc
5228  lastAlgorithm_=999;
5229#ifndef NDEBUG
5230  int n = modelPtr_->numberColumns();
5231  if (elementIndex<0||elementIndex>=n) {
5232    indexError(elementIndex,"setColBounds");
5233  }
5234#endif
5235  if (!modelPtr_->lower_)
5236    modelPtr_->whatsChanged_ &= ~0xffff; // switch off
5237  modelPtr_->setColumnBounds(elementIndex,lower,upper);
5238}
5239void OsiClpSolverInterface::setColSetBounds(const int* indexFirst,
5240                                            const int* indexLast,
5241                                            const double* boundList)
5242{
5243  modelPtr_->whatsChanged_ &= 0x1ffff;
5244  // Say can't gurantee optimal basis etc
5245  lastAlgorithm_=999;
5246#ifndef NDEBUG
5247  int n = modelPtr_->numberColumns();
5248  const int * indexFirst2=indexFirst;
5249  while (indexFirst2 != indexLast) {
5250    const int iColumn=*indexFirst2++;
5251    if (iColumn<0||iColumn>=n) {
5252      indexError(iColumn,"setColSetBounds");
5253    }
5254  }
5255#endif
5256  modelPtr_->setColSetBounds(indexFirst,indexLast,boundList);
5257}
5258//------------------------------------------------------------------
5259/* Set a single row lower bound<br>
5260   Use -DBL_MAX for -infinity. */
5261void 
5262OsiClpSolverInterface::setRowLower( int elementIndex, double elementValue ) {
5263  // Say can't gurantee optimal basis etc
5264  lastAlgorithm_=999;
5265  modelPtr_->whatsChanged_ &= 0xffff;
5266#ifndef NDEBUG
5267  int n = modelPtr_->numberRows();
5268  if (elementIndex<0||elementIndex>=n) {
5269    indexError(elementIndex,"setRowLower");
5270  }
5271#endif
5272  modelPtr_->setRowLower(elementIndex , elementValue);
5273  if (rowsense_!=NULL) {
5274    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5275    convertBoundToSense(modelPtr_->rowLower_[elementIndex], 
5276                        modelPtr_->rowUpper_[elementIndex],
5277                        rowsense_[elementIndex], rhs_[elementIndex], rowrange_[elementIndex]);
5278  }
5279}
5280     
5281/* Set a single row upper bound<br>
5282   Use DBL_MAX for infinity. */
5283void 
5284OsiClpSolverInterface::setRowUpper( int elementIndex, double elementValue ) {
5285  modelPtr_->whatsChanged_ &= 0xffff;
5286  // Say can't guarantee optimal basis etc
5287  lastAlgorithm_=999;
5288#ifndef NDEBUG
5289  int n = modelPtr_->numberRows();
5290  if (elementIndex<0||elementIndex>=n) {
5291    indexError(elementIndex,"setRowUpper");
5292  }
5293#endif
5294  modelPtr_->setRowUpper(elementIndex , elementValue);
5295  if (rowsense_!=NULL) {
5296    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5297    convertBoundToSense(modelPtr_->rowLower_[elementIndex], 
5298                        modelPtr_->rowUpper_[elementIndex],
5299                        rowsense_[elementIndex], rhs_[elementIndex], rowrange_[elementIndex]);
5300  }
5301}
5302   
5303/* Set a single row lower and upper bound */
5304void 
5305OsiClpSolverInterface::setRowBounds( int elementIndex,
5306              double lower, double upper ) {
5307  modelPtr_->whatsChanged_ &= 0xffff;
5308  // Say can't gurantee optimal basis etc
5309  lastAlgorithm_=999;
5310#ifndef NDEBUG
5311  int n = modelPtr_->numberRows();
5312  if (elementIndex<0||elementIndex>=n) {
5313    indexError(elementIndex,"setRowBounds");
5314  }
5315#endif
5316  modelPtr_->setRowBounds(elementIndex,lower,upper);
5317  if (rowsense_!=NULL) {
5318    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5319    convertBoundToSense(modelPtr_->rowLower_[elementIndex], 
5320                        modelPtr_->rowUpper_[elementIndex],
5321                        rowsense_[elementIndex], rhs_[elementIndex], rowrange_[elementIndex]);
5322  }
5323}
5324//-----------------------------------------------------------------------------
5325void
5326OsiClpSolverInterface::setRowType(int i, char sense, double rightHandSide,
5327                                  double range)
5328{
5329  modelPtr_->whatsChanged_ &= 0xffff;
5330  // Say can't gurantee optimal basis etc
5331  lastAlgorithm_=999;
5332#ifndef NDEBUG
5333  int n = modelPtr_->numberRows();
5334  if (i<0||i>=n) {
5335    indexError(i,"setRowType");
5336  }
5337#endif
5338  double lower = 0, upper = 0;
5339  convertSenseToBound(sense, rightHandSide, range, lower, upper);
5340  setRowBounds(i, lower, upper);
5341  // If user is using sense then set
5342  if (rowsense_) {
5343    rowsense_[i] = sense;
5344    rhs_[i] = rightHandSide;
5345    rowrange_[i] = range;
5346  }
5347}
5348// Set name of row
5349void 
5350//OsiClpSolverInterface::setRowName(int rowIndex, std::string & name)
5351OsiClpSolverInterface::setRowName(int rowIndex, std::string name) 
5352{
5353  if (rowIndex>=0&&rowIndex<modelPtr_->numberRows()) {
5354    int nameDiscipline;
5355    getIntParam(OsiNameDiscipline,nameDiscipline) ;
5356    if (nameDiscipline) {
5357      modelPtr_->setRowName(rowIndex,name);
5358      OsiSolverInterface::setRowName(rowIndex,name) ;
5359    }
5360  }
5361}
5362// Return name of row if one exists or Rnnnnnnn
5363// we ignore maxLen
5364std::string
5365OsiClpSolverInterface::getRowName(int rowIndex, unsigned int /*maxLen*/) const
5366{ 
5367  if (rowIndex == getNumRows())
5368    return getObjName();
5369  int useNames;
5370  getIntParam (OsiNameDiscipline,useNames);
5371  if (useNames)
5372    return modelPtr_->getRowName(rowIndex);
5373  else
5374    return dfltRowColName('r',rowIndex);
5375}
5376   
5377// Set name of col
5378void 
5379//OsiClpSolverInterface::setColName(int colIndex, std::string & name)
5380OsiClpSolverInterface::setColName(int colIndex, std::string name) 
5381{
5382  if (colIndex>=0&&colIndex<modelPtr_->numberColumns()) {
5383    int nameDiscipline;
5384    getIntParam(OsiNameDiscipline,nameDiscipline) ;
5385    if (nameDiscipline) {
5386      modelPtr_->setColumnName(colIndex,name);
5387      OsiSolverInterface::setColName(colIndex,name) ;
5388    }
5389  }
5390}
5391// Return name of col if one exists or Rnnnnnnn
5392std::string
5393OsiClpSolverInterface::getColName(int colIndex, unsigned int /*maxLen*/) const
5394{
5395  int useNames;
5396  getIntParam (OsiNameDiscipline,useNames);
5397  if (useNames)
5398    return modelPtr_->getColumnName(colIndex);
5399  else
5400    return dfltRowColName('c',colIndex);
5401}
5402   
5403   
5404//-----------------------------------------------------------------------------
5405void OsiClpSolverInterface::setRowSetBounds(const int* indexFirst,
5406                                            const int* indexLast,
5407                                            const double* boundList)
5408{
5409  modelPtr_->whatsChanged_ &= 0xffff;
5410  // Say can't gurantee optimal basis etc
5411  lastAlgorithm_=999;
5412#ifndef NDEBUG
5413  int n = modelPtr_->numberRows();
5414  const int * indexFirst2=indexFirst;
5415  while (indexFirst2 != indexLast) {
5416    const int iColumn=*indexFirst2++;
5417    if (iColumn<0||iColumn>=n) {
5418      indexError(iColumn,"setColumnSetBounds");
5419    }
5420  }
5421#endif
5422  modelPtr_->setRowSetBounds(indexFirst,indexLast,boundList);
5423  if (rowsense_ != NULL) {
5424    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5425    double * lower = modelPtr_->rowLower();
5426    double * upper = modelPtr_->rowUpper();
5427    while (indexFirst != indexLast) {
5428      const int iRow=*indexFirst++;
5429      convertBoundToSense(lower[iRow], upper[iRow],
5430                          rowsense_[iRow], rhs_[iRow], rowrange_[iRow]);
5431    }
5432  }
5433}
5434//-----------------------------------------------------------------------------
5435void
5436OsiClpSolverInterface::setRowSetTypes(const int* indexFirst,
5437                                      const int* indexLast,
5438                                      const char* senseList,
5439                                      const double* rhsList,
5440                                      const double* rangeList)
5441{
5442  modelPtr_->whatsChanged_ &= 0xffff;
5443  // Say can't gurantee optimal basis etc
5444  lastAlgorithm_=999;
5445#ifndef NDEBUG
5446  int n = modelPtr_->numberRows();
5447#endif
5448  const int len = static_cast<int>(indexLast - indexFirst);
5449  while (indexFirst != indexLast) {
5450    const int iRow= *indexFirst++;
5451#ifndef NDEBUG
5452    if (iRow<0||iRow>=n) {
5453      indexError(iRow,"isContinuous");
5454    }
5455#endif
5456    double lowerValue = 0;
5457    double upperValue = 0;
5458    if (rangeList){
5459      convertSenseToBound(*senseList++, *rhsList++, *rangeList++,
5460                          lowerValue, upperValue);
5461    } else {
5462      convertSenseToBound(*senseList++, *rhsList++, 0,
5463                          lowerValue, upperValue);
5464    }
5465    modelPtr_->setRowBounds(iRow,lowerValue,upperValue);
5466  }
5467  if (rowsense_ != NULL) {
5468    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5469    indexFirst -= len;
5470    senseList -= len;
5471    rhsList -= len;
5472    if (rangeList)
5473       rangeList -= len;
5474    while (indexFirst != indexLast) {
5475      const int iRow=*indexFirst++;
5476      rowsense_[iRow] = *senseList++;
5477      rhs_[iRow] = *rhsList++;
5478      if (rangeList)
5479         rowrange_[iRow] = *rangeList++;
5480    }
5481  }
5482}
5483
5484/*
5485  Clp's copy-in/copy-out design paradigm is a challenge for the simplex modes.
5486  Normal operation goes like this:
5487    * startup() loads clp's work arrays, performing scaling for numerical
5488      stability and compensating for max.
5489    * clp solves the problem
5490    * finish() unloads the work arrays into answer arrays, undoing scaling
5491      and max compensation.
5492  There are two solutions: undo scaling and max on demand, or make them
5493  into noops. The various getBInv* methods undo scaling on demand (but
5494  see special option 512) and do not need to worry about max. Other get
5495  solution methods are not coded to do this, so the second approach is
5496  used. For simplex modes, turn off scaling (necessary for both primal and
5497  dual solutions) and temporarily convert max to min (necessary for dual
5498  solution). This makes the unscaling in getBInv* superfluous, but don't
5499  remove it. Arguably the better solution here would be to go through and
5500  add unscaling and max compensation to the get solution methods. Look for
5501  fakeMinInSimplex to see the places this propagates to.
5502
5503  TODO: setRowPrice never has worked properly, and I didn't try to fix it in
5504        this go-round.
5505
5506  As of 100907, change applied to [enable|disable]Factorization (mode 1).
5507  Limitation of [enable|disable]SimplexInterface (mode 2) noted in
5508  documentation.  -- lh, 100907 --
5509*/
5510/*
5511  Enables normal operation of subsequent functions.  This method is supposed
5512  to ensure that all typical things (like reduced costs, etc.) are updated
5513  when individual pivots are executed and can be queried by other methods
5514*/
5515void 
5516OsiClpSolverInterface::enableSimplexInterface(bool doingPrimal)
5517{
5518  modelPtr_->whatsChanged_ &= 0xffff;
5519  if (modelPtr_->solveType()==2)
5520    return;
5521  assert (modelPtr_->solveType()==1);
5522  int saveIts = modelPtr_->numberIterations_;
5523  modelPtr_->setSolveType(2);
5524  if (doingPrimal)
5525    modelPtr_->setAlgorithm(1);
5526  else
5527    modelPtr_->setAlgorithm(-1);
5528  // Do initialization
5529  saveData_ = modelPtr_->saveData();
5530  saveData_.scalingFlag_=modelPtr_->scalingFlag();
5531  modelPtr_->scaling(0);
5532  specialOptions_ = 0x80000000;
5533  // set infeasibility cost up
5534  modelPtr_->setInfeasibilityCost(1.0e12);
5535  ClpDualRowDantzig dantzig;
5536  modelPtr_->setDualRowPivotAlgorithm(dantzig);
5537  ClpPrimalColumnDantzig dantzigP;
5538  dantzigP.saveWeights(modelPtr_,0); // set modelPtr
5539  modelPtr_->setPrimalColumnPivotAlgorithm(dantzigP);
5540  int saveOptions = modelPtr_->specialOptions_;
5541  modelPtr_->specialOptions_ &= ~262144;
5542  delete modelPtr_->scaledMatrix_;
5543  modelPtr_->scaledMatrix_=NULL;
5544  // make sure using standard factorization
5545  modelPtr_->factorization()->forceOtherFactorization(4);
5546#ifdef NDEBUG
5547  modelPtr_->startup(0);
5548#else
5549  int returnCode=modelPtr_->startup(0);
5550  assert (!returnCode||returnCode==2);
5551#endif
5552  modelPtr_->specialOptions_=saveOptions;
5553  modelPtr_->numberIterations_=saveIts;
5554}
5555
5556//Undo whatever setting changes the above method had to make
5557void 
5558OsiClpSolverInterface::disableSimplexInterface()
5559{
5560  modelPtr_->whatsChanged_ &= 0xffff;
5561  assert (modelPtr_->solveType()==2);
5562  // declare optimality anyway  (for message handler)
5563  modelPtr_->setProblemStatus(0);
5564  modelPtr_->setSolveType(1);
5565  // message will not appear anyway
5566  int saveMessageLevel=modelPtr_->messageHandler()->logLevel();
5567  modelPtr_->messageHandler()->setLogLevel(0);
5568  modelPtr_->finish();
5569  modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
5570  modelPtr_->restoreData(saveData_);
5571  modelPtr_->scaling(saveData_.scalingFlag_);
5572  ClpDualRowSteepest steepest;
5573  modelPtr_->setDualRowPivotAlgorithm(steepest);
5574  ClpPrimalColumnSteepest steepestP;
5575  modelPtr_->setPrimalColumnPivotAlgorithm(steepestP);
5576  basis_ = getBasis(modelPtr_);
5577  modelPtr_->setSolveType(1);
5578}
5579
5580/*
5581  Force scaling off. If the client thinks we're maximising, arrange it so
5582  that clp sees minimisation while the client still sees maximisation. In
5583  keeping with the spirit of the getBInv methods, special option 512 will
5584  leave all work to the client.
5585*/
5586void 
5587OsiClpSolverInterface::enableFactorization() const
5588{
5589  saveData_.specialOptions_=specialOptions_;
5590  // Try to preserve work regions, reuse factorization
5591  if ((specialOptions_&(1+8))!=1+8)
5592    setSpecialOptionsMutable((1+8)|specialOptions_);
5593  // Are we allowed to make the output sensible to mere mortals?
5594  if ((specialOptions_&512)==0) {
5595    // Force scaling to off
5596    saveData_.scalingFlag_ = modelPtr_->scalingFlag() ;
5597    modelPtr_->scaling(0) ;
5598    // Temporarily force to min but keep a copy of original objective.
5599    if (getObjSense() < 0.0) {
5600      fakeMinInSimplex_ = true ;
5601      modelPtr_->setOptimizationDirection(1.0) ;
5602      double *c = modelPtr_->objective() ;
5603      int n = getNumCols() ;
5604      linearObjective_ = new double[n] ;
5605      CoinMemcpyN(c,n,linearObjective_) ;
5606      std::transform(c,c+n,c,std::negate<double>()) ;
5607    }
5608  }
5609  int saveStatus = modelPtr_->problemStatus_;
5610#ifdef NDEBUG
5611  modelPtr_->startup(0);
5612#else
5613  int returnCode=modelPtr_->startup(0);
5614  assert (!returnCode||returnCode==2);
5615#endif
5616  modelPtr_->problemStatus_=saveStatus;
5617}
5618
5619/*
5620  Undo enableFactorization. Retrieve the special options and scaling and
5621  remove the temporary objective used to fake minimisation in clp.
5622*/
5623void 
5624OsiClpSolverInterface::disableFactorization() const
5625{
5626  specialOptions_=saveData_.specialOptions_;
5627  // declare optimality anyway  (for message handler)
5628  modelPtr_->setProblemStatus(0);
5629  // message will not appear anyway
5630  int saveMessageLevel=modelPtr_->messageHandler()->logLevel();
5631  modelPtr_->messageHandler()->setLogLevel(0);
5632  modelPtr_->finish();
5633  modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
5634  // Client asked for transforms on the way in, so back out.
5635  if ((specialOptions_&512)==0) {
5636    modelPtr_->scaling(saveData_.scalingFlag_) ;
5637    if (fakeMinInSimplex_ == true) {
5638      fakeMinInSimplex_ = false ;
5639      modelPtr_->setOptimizationDirection(-1.0) ;
5640      double *c = modelPtr_->objective() ;
5641      int n = getNumCols() ;
5642      std::transform(c,c+n,c,std::negate<double>()) ;
5643      delete[] linearObjective_ ;
5644    }
5645  }
5646}
5647
5648
5649
5650/* The following two methods may be replaced by the
5651   methods of OsiSolverInterface using OsiWarmStartBasis if:
5652   1. OsiWarmStartBasis resize operation is implemented
5653   more efficiently and
5654   2. It is ensured that effects on the solver are the same
5655   
5656   Returns a basis status of the structural/artificial variables
5657*/
5658void 
5659OsiClpSolverInterface::getBasisStatus(int* cstat, int* rstat) const
5660{
5661  int iRow,iColumn;
5662  int numberRows = modelPtr_->numberRows();
5663  int numberColumns = modelPtr_->numberColumns();
5664  const double * pi = modelPtr_->dualRowSolution();
5665  const double * dj = modelPtr_->dualColumnSolution();