source: stable/1.15/Clp/src/OsiClp/OsiClpSolverInterface.cpp @ 1931

Last change on this file since 1931 was 1928, checked in by stefan, 7 years ago

revert r1925 (so we can merge r1926) and sync with rev 1927

File size: 301.2 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    // Start of fast iterations
2241    bool alwaysFinish= ((specialOptions_&32)==0) ? true : false;
2242    //modelPtr_->setLogLevel(1);
2243    modelPtr_->setIntParam(ClpMaxNumIteration,itlim);
2244    int saveNumberFake = (static_cast<ClpSimplexDual *>(modelPtr_))->numberFake_;
2245    int status = (static_cast<ClpSimplexDual *>(modelPtr_))->fastDual(alwaysFinish);
2246    (static_cast<ClpSimplexDual *>(modelPtr_))->numberFake_ = saveNumberFake;
2247   
2248    int problemStatus = modelPtr_->problemStatus();
2249    double objectiveValue =modelPtr_->objectiveValue() * modelPtr_->optimizationDirection();
2250    CoinAssert (modelPtr_->problemStatus()||modelPtr_->objectiveValue()<1.0e50);
2251    // make sure plausible
2252    double obj = CoinMax(objectiveValue,saveObjectiveValue);
2253    if (problemStatus==10||problemStatus<0) {
2254      // was trying to clean up or something odd
2255      if (problemStatus==10)
2256        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2257      status=1;
2258    }
2259    if (status) {
2260      // not finished - might be optimal
2261      modelPtr_->checkPrimalSolution(modelPtr_->solutionRegion(0),
2262                                     modelPtr_->solutionRegion(1));
2263      //modelPtr_->gutsOfSolution(NULL,NULL,0);
2264      //if (problemStatus==3)
2265      //modelPtr_->computeObjectiveValue();
2266      objectiveValue =modelPtr_->objectiveValue() *
2267        modelPtr_->optimizationDirection();
2268      obj = CoinMax(objectiveValue,saveObjectiveValue);
2269      if (!modelPtr_->numberDualInfeasibilities()) { 
2270        double limit = 0.0;
2271        modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2272        if (modelPtr_->secondaryStatus()==1&&!problemStatus&&obj<limit) {
2273          obj=limit;
2274          problemStatus=3;
2275        }
2276        if (!modelPtr_->numberPrimalInfeasibilities()&&obj<limit) { 
2277          problemStatus=0;
2278        } else if (problemStatus==10) {
2279          problemStatus=3;
2280        } else if (!modelPtr_->numberPrimalInfeasibilities()) {
2281          problemStatus=1; // infeasible
2282        } 
2283      } else {
2284        // can't say much
2285        //if (problemStatus==3)
2286        //modelPtr_->computeObjectiveValue();
2287        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2288        problemStatus=3;
2289      }
2290    } else if (!problemStatus) {
2291      if (modelPtr_->isDualObjectiveLimitReached()) 
2292        problemStatus=1; // infeasible
2293    }
2294    if (status&&!problemStatus) {
2295      problemStatus=3; // can't be sure
2296      lastAlgorithm_=1;
2297    }
2298    if (problemStatus<0)
2299      problemStatus=3;
2300    modelPtr_->setProblemStatus(problemStatus);
2301    modelPtr_->setObjectiveValue(obj*modelPtr_->optimizationDirection());
2302    double * solution = modelPtr_->primalColumnSolution();
2303    const double * solution2 = modelPtr_->solutionRegion();
2304    // could just do changed bounds - also try double size scale so can use * not /
2305    if (!columnScale) {
2306      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2307        solution[iColumn]= solution2[iColumn];
2308      }
2309    } else {
2310      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2311        solution[iColumn]= solution2[iColumn]*columnScale[iColumn];
2312      }
2313    }
2314    CoinMemcpyN(saveLowerOriginal,numberColumns,columnLower);
2315    CoinMemcpyN(saveUpperOriginal,numberColumns,columnUpper);
2316#if 0
2317    // could combine with loop above
2318    if (modelPtr_==modelPtr_)
2319      modelPtr_->computeObjectiveValue();
2320    if (status&&!problemStatus) {
2321      memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2322      modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2323      modelPtr_->checkSolutionInternal();
2324      //modelPtr_->setLogLevel(1);
2325      //modelPtr_->allSlackBasis();
2326      //modelPtr_->primal(1);
2327      //memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2328      //modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2329      //modelPtr_->checkSolutionInternal();
2330      assert (!modelPtr_->problemStatus());
2331    }
2332#endif
2333    // and back bounds
2334    CoinMemcpyN(saveLower,number,modelPtr_->lowerRegion());
2335    CoinMemcpyN(saveUpper,number,modelPtr_->upperRegion());
2336    if (extraCopy) {
2337      CoinMemcpyN(saveLower,number,modelPtr_->lowerRegion()+extraCopy);
2338      CoinMemcpyN(saveUpper,number,modelPtr_->upperRegion()+extraCopy);
2339    }
2340    modelPtr_->setIntParam(ClpMaxNumIteration,itlimOrig_);
2341    if (savedObjective) {
2342      // fix up
2343      modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2344      //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2345      modelPtr_->objective_=savedObjective;
2346      if (!modelPtr_->problemStatus_) {
2347        CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
2348        CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
2349        if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
2350          CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
2351        modelPtr_->computeObjectiveValue();
2352      }
2353    }
2354    return;
2355  }
2356#endif
2357  if (smallModel_==NULL) {
2358    setWarmStart(ws_);
2359    CoinMemcpyN(           rowActivity_,numberRows,modelPtr_->primalRowSolution());
2360    CoinMemcpyN(columnActivity_,numberColumns,modelPtr_->primalColumnSolution());
2361    modelPtr_->setIntParam(ClpMaxNumIteration,CoinMin(itlim,9999));
2362    resolve();
2363  } else {
2364    assert (spareArrays_);
2365    double * arrayD = reinterpret_cast<double *> (spareArrays_);
2366    double saveObjectiveValue = arrayD[0];
2367    double * saveSolution = arrayD+1;
2368    // double check arrays exist (? for nonlinear)
2369    //if (!smallModel_->solutionRegion())
2370    //smallModel_->createRim(63);
2371    int numberRows2 = smallModel_->numberRows();
2372    int numberColumns2 = smallModel_->numberColumns();
2373    int number = numberRows2+numberColumns2;
2374    CoinMemcpyN(saveSolution,number,smallModel_->solutionRegion());
2375    double * saveLower = saveSolution + (numberRows+numberColumns);
2376    CoinMemcpyN(saveLower,number,smallModel_->lowerRegion());
2377    double * saveUpper = saveLower + (numberRows+numberColumns);
2378    CoinMemcpyN(saveUpper,number,smallModel_->upperRegion());
2379    double * saveObjective = saveUpper + (numberRows+numberColumns);
2380    CoinMemcpyN(saveObjective,number,smallModel_->costRegion());
2381    double * saveLowerOriginal = saveObjective + (numberRows+numberColumns);
2382    double * saveUpperOriginal = saveLowerOriginal + numberColumns;
2383    arrayD = saveUpperOriginal + numberColumns;
2384    int * savePivot = reinterpret_cast<int *> (arrayD);
2385    CoinMemcpyN(savePivot,numberRows2,smallModel_->pivotVariable());
2386    int * whichRow = savePivot+numberRows;
2387    int * whichColumn = whichRow + 3*numberRows;
2388    int * arrayI = whichColumn + 2*numberColumns;
2389    unsigned char * saveStatus = reinterpret_cast<unsigned char *> (arrayI+1);
2390    CoinMemcpyN(saveStatus,number,smallModel_->statusArray());
2391    /* If factorization_ NULL then infeasible
2392       not really sure how could have slipped through.
2393       But this can't make situation worse */
2394    if (factorization_) 
2395      smallModel_->setFactorization(*factorization_);
2396    //int * backColumn = whichColumn+numberColumns;
2397    const double * lowerBig = modelPtr_->columnLower();
2398    const double * upperBig = modelPtr_->columnUpper();
2399    // make sure whatsChanged_ has 1 set
2400    smallModel_->setWhatsChanged(511);
2401    double * lowerSmall = smallModel_->lowerRegion();
2402    double * upperSmall = smallModel_->upperRegion();
2403    double * lowerSmallReal = smallModel_->columnLower();
2404    double * upperSmallReal = smallModel_->columnUpper();
2405    int i;
2406    double rhsScale = smallModel_->rhsScale();
2407    const double * columnScale = NULL;
2408    if (smallModel_->scalingFlag()>0) {
2409      columnScale = smallModel_->columnScale();
2410    }
2411    // and do bounds in case dual needs them
2412    // may be infeasible
2413    for (i=0;i<numberColumns2;i++) {
2414      int iColumn = whichColumn[i];
2415      if (lowerBig[iColumn]>saveLowerOriginal[iColumn]) {
2416        double value = lowerBig[iColumn];
2417        lowerSmallReal[i]=value;
2418        value *= rhsScale;
2419        if (columnScale)
2420          value /= columnScale[i];
2421        lowerSmall[i]=value;
2422      }
2423      if (upperBig[iColumn]<saveUpperOriginal[iColumn]) {
2424        double value = upperBig[iColumn];
2425        upperSmallReal[i]=value;
2426        value *= rhsScale;
2427        if (columnScale)
2428          value /= columnScale[i];
2429        upperSmall[i]=value;
2430      }
2431      if (upperSmall[i]<lowerSmall[i]-1.0e-8)
2432        break;
2433    }
2434    /* If factorization_ NULL then infeasible
2435       not really sure how could have slipped through.
2436       But this can't make situation worse */
2437    bool infeasible = (i<numberColumns2||!factorization_);
2438    // Start of fast iterations
2439    bool alwaysFinish= ((specialOptions_&32)==0) ? true : false;
2440    //smallModel_->setLogLevel(1);
2441    smallModel_->setIntParam(ClpMaxNumIteration,itlim);
2442    int saveNumberFake = (static_cast<ClpSimplexDual *>(smallModel_))->numberFake_;
2443    int status;
2444    if (!infeasible) {
2445      status = static_cast<ClpSimplexDual *>(smallModel_)->fastDual(alwaysFinish);
2446    } else {
2447      status=0;
2448      smallModel_->setProblemStatus(1);
2449    }
2450    (static_cast<ClpSimplexDual *>(smallModel_))->numberFake_ = saveNumberFake;
2451    if (smallModel_->numberIterations()==-98) {
2452      printf("rrrrrrrrrrrr\n");
2453      smallModel_->checkPrimalSolution(smallModel_->solutionRegion(0),
2454                                     smallModel_->solutionRegion(1));
2455      //smallModel_->gutsOfSolution(NULL,NULL,0);
2456      //if (problemStatus==3)
2457      //smallModel_->computeObjectiveValue();
2458      printf("robj %g\n",smallModel_->objectiveValue() *
2459             modelPtr_->optimizationDirection());
2460      writeMps("rr.mps");
2461      smallModel_->writeMps("rr_small.mps");
2462      ClpSimplex temp = *smallModel_;
2463      printf("small\n");
2464      temp.setLogLevel(63);
2465      temp.dual();
2466      double limit = 0.0;
2467      modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2468      if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2469        printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2470               temp.objectiveValue(),
2471               smallModel_->objectiveOffset(),temp.objectiveOffset());
2472      }
2473      printf("big\n");
2474      temp = *modelPtr_;
2475      temp.dual();
2476      if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2477        printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2478               temp.objectiveValue(),
2479               smallModel_->objectiveOffset(),temp.objectiveOffset());
2480      }
2481    }
2482    int problemStatus = smallModel_->problemStatus();
2483    double objectiveValue =smallModel_->objectiveValue() * modelPtr_->optimizationDirection();
2484    CoinAssert (smallModel_->problemStatus()||smallModel_->objectiveValue()<1.0e50);
2485    // make sure plausible
2486    double obj = CoinMax(objectiveValue,saveObjectiveValue);
2487    if (problemStatus==10||problemStatus<0) {
2488      // was trying to clean up or something odd
2489      if (problemStatus==10)
2490        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2491      status=1;
2492    }
2493    if (status) {
2494      // not finished - might be optimal
2495      smallModel_->checkPrimalSolution(smallModel_->solutionRegion(0),
2496                                     smallModel_->solutionRegion(1));
2497      //smallModel_->gutsOfSolution(NULL,NULL,0);
2498      //if (problemStatus==3)
2499      //smallModel_->computeObjectiveValue();
2500      objectiveValue =smallModel_->objectiveValue() *
2501        modelPtr_->optimizationDirection();
2502      if (problemStatus!=10)
2503        obj = CoinMax(objectiveValue,saveObjectiveValue);
2504      if (!smallModel_->numberDualInfeasibilities()) { 
2505        double limit = 0.0;
2506        modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2507        if (smallModel_->secondaryStatus()==1&&!problemStatus&&obj<limit) {
2508#if 0
2509          // switch off
2510          ClpSimplex temp = *smallModel_;
2511          temp.dual();
2512          if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2513            printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2514                   temp.objectiveValue(),
2515                   smallModel_->objectiveOffset(),temp.objectiveOffset());
2516          }
2517          lastAlgorithm_=1;
2518          obj=limit;
2519          problemStatus=10;
2520#else
2521          obj=limit;
2522          problemStatus=3;
2523#endif
2524        }
2525        if (!smallModel_->numberPrimalInfeasibilities()&&obj<limit) { 
2526          problemStatus=0;
2527#if 0
2528          ClpSimplex temp = *smallModel_;
2529          temp.dual();
2530          if (temp.numberIterations())
2531            printf("temp iterated\n");
2532          assert (temp.problemStatus()==0&&temp.objectiveValue()<limit);
2533#endif
2534        } else if (problemStatus==10) {
2535          problemStatus=3;
2536        } else if (!smallModel_->numberPrimalInfeasibilities()) {
2537          problemStatus=1; // infeasible
2538        } 
2539      } else {
2540        // can't say much
2541        //if (problemStatus==3)
2542        //smallModel_->computeObjectiveValue();
2543        lastAlgorithm_=1; // so won't fail on cutoff (in CbcNode)
2544        problemStatus=3;
2545      }
2546    } else if (!problemStatus) {
2547      if (smallModel_->isDualObjectiveLimitReached()) 
2548        problemStatus=1; // infeasible
2549    }
2550    if (status&&!problemStatus) {
2551      problemStatus=3; // can't be sure
2552      lastAlgorithm_=1;
2553    }
2554    if (problemStatus<0)
2555      problemStatus=3;
2556    modelPtr_->setProblemStatus(problemStatus);
2557    modelPtr_->setObjectiveValue(obj*modelPtr_->optimizationDirection());
2558    modelPtr_->setSumDualInfeasibilities(smallModel_->sumDualInfeasibilities());
2559    modelPtr_->setNumberDualInfeasibilities(smallModel_->numberDualInfeasibilities());
2560    modelPtr_->setSumPrimalInfeasibilities(smallModel_->sumPrimalInfeasibilities());
2561    modelPtr_->setNumberPrimalInfeasibilities(smallModel_->numberPrimalInfeasibilities());
2562    double * solution = modelPtr_->primalColumnSolution();
2563    const double * solution2 = smallModel_->solutionRegion();
2564    if (!columnScale) {
2565      for (i=0;i<numberColumns2;i++) {
2566        int iColumn = whichColumn[i];
2567        solution[iColumn]= solution2[i];
2568        lowerSmallReal[i]=saveLowerOriginal[iColumn];
2569        upperSmallReal[i]=saveUpperOriginal[iColumn];
2570      }
2571    } else {
2572      for (i=0;i<numberColumns2;i++) {
2573        int iColumn = whichColumn[i];
2574        solution[iColumn]= solution2[i]*columnScale[i];
2575        lowerSmallReal[i]=saveLowerOriginal[iColumn];
2576        upperSmallReal[i]=saveUpperOriginal[iColumn];
2577      }
2578    }
2579    // could combine with loop above
2580    if (modelPtr_==smallModel_)
2581      modelPtr_->computeObjectiveValue();
2582#if 1
2583    if (status&&!problemStatus) {
2584      memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2585      modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2586      modelPtr_->checkSolutionInternal();
2587      //modelPtr_->setLogLevel(1);
2588      //modelPtr_->allSlackBasis();
2589      //modelPtr_->primal(1);
2590      //memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2591      //modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2592      //modelPtr_->checkSolutionInternal();
2593      assert (!modelPtr_->problemStatus());
2594    }
2595#endif
2596    modelPtr_->setNumberIterations(smallModel_->numberIterations());
2597    // and back bounds
2598    CoinMemcpyN(saveLower,number,smallModel_->lowerRegion());
2599    CoinMemcpyN(saveUpper,number,smallModel_->upperRegion());
2600  }
2601  if (savedObjective) {
2602    // fix up
2603    modelPtr_->dblParam_[ClpDualObjectiveLimit]=savedDualLimit;
2604    //modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2605    modelPtr_->objective_=savedObjective;
2606    if (!modelPtr_->problemStatus_) {
2607      CoinZeroN(modelPtr_->dual_,modelPtr_->numberRows_);
2608      CoinZeroN(modelPtr_->reducedCost_,modelPtr_->numberColumns_);
2609      if (modelPtr_->dj_&&(modelPtr_->whatsChanged_&1)!=0)
2610        CoinZeroN(modelPtr_->dj_,modelPtr_->numberColumns_+modelPtr_->numberRows_);
2611      modelPtr_->computeObjectiveValue();
2612    }
2613  }
2614  modelPtr_->setIntParam(ClpMaxNumIteration,itlimOrig_);
2615}
2616
2617void OsiClpSolverInterface::unmarkHotStart()
2618{
2619#ifdef CLEAN_HOT_START
2620  if ((specialOptions_&65536)!=0) {
2621    modelPtr_->setLogLevel(saveData_.scalingFlag_);
2622    modelPtr_->deleteRim(0);
2623    if (lastNumberRows_<0) {
2624      specialOptions_ |= 131072;
2625      lastNumberRows_ = -1 -lastNumberRows_;
2626      if (modelPtr_->rowScale_) {
2627        if (modelPtr_->rowScale_!=rowScale_.array()) {
2628          delete [] modelPtr_->rowScale_;
2629          delete [] modelPtr_->columnScale_;
2630        }
2631        modelPtr_->rowScale_=NULL;
2632        modelPtr_->columnScale_=NULL;
2633      }
2634    }
2635    delete factorization_;
2636    delete [] spareArrays_;
2637    smallModel_=NULL;
2638    spareArrays_=NULL;
2639    factorization_=NULL;
2640    delete [] rowActivity_;
2641    delete [] columnActivity_;
2642    rowActivity_=NULL;
2643    columnActivity_=NULL;
2644    return;
2645  }
2646#endif
2647  if (smallModel_==NULL) {
2648    setWarmStart(ws_);
2649    int numberRows = modelPtr_->numberRows();
2650    int numberColumns = modelPtr_->numberColumns();
2651    CoinMemcpyN(           rowActivity_,numberRows,modelPtr_->primalRowSolution());
2652    CoinMemcpyN(columnActivity_,numberColumns,modelPtr_->primalColumnSolution());
2653    delete ws_;
2654    ws_ = NULL;
2655  } else {
2656#ifndef KEEP_SMALL
2657    if (smallModel_!=modelPtr_)
2658      delete smallModel_;
2659    smallModel_=NULL;
2660#else
2661    if (smallModel_==modelPtr_) {
2662      smallModel_=NULL;
2663    } else if (smallModel_) {
2664      if(!spareArrays_) {
2665        delete smallModel_;
2666        smallModel_=NULL;
2667        delete factorization_;
2668        factorization_=NULL;
2669      } else {
2670        static_cast<ClpSimplexDual *> (smallModel_)->cleanupAfterStrongBranching( factorization_);
2671        if ((smallModel_->specialOptions_&4096)==0) {
2672          delete factorization_;
2673        }
2674      }
2675    }
2676#endif
2677    //delete [] spareArrays_;
2678    //spareArrays_=NULL;
2679    factorization_=NULL;
2680  }
2681  delete [] rowActivity_;
2682  delete [] columnActivity_;
2683  rowActivity_=NULL;
2684  columnActivity_=NULL;
2685  // Make sure whatsChanged not out of sync
2686  if (!modelPtr_->columnUpperWork_)
2687    modelPtr_->whatsChanged_ &= ~0xffff;
2688  modelPtr_->specialOptions_ = saveData_.specialOptions_; 
2689}
2690
2691//#############################################################################
2692// Problem information methods (original data)
2693//#############################################################################
2694
2695//------------------------------------------------------------------
2696const char * OsiClpSolverInterface::getRowSense() const
2697{
2698  extractSenseRhsRange();
2699  return rowsense_;
2700}
2701//------------------------------------------------------------------
2702const double * OsiClpSolverInterface::getRightHandSide() const
2703{
2704  extractSenseRhsRange();
2705  return rhs_;
2706}
2707//------------------------------------------------------------------
2708const double * OsiClpSolverInterface::getRowRange() const
2709{
2710  extractSenseRhsRange();
2711  return rowrange_;
2712}
2713//------------------------------------------------------------------
2714// Return information on integrality
2715//------------------------------------------------------------------
2716bool OsiClpSolverInterface::isContinuous(int colNumber) const
2717{
2718  if ( integerInformation_==NULL ) return true;
2719#ifndef NDEBUG
2720  int n = modelPtr_->numberColumns();
2721  if (colNumber<0||colNumber>=n) {
2722    indexError(colNumber,"isContinuous");
2723  }
2724#endif
2725  if ( integerInformation_[colNumber]==0 ) return true;
2726  return false;
2727}
2728bool 
2729OsiClpSolverInterface::isBinary(int colNumber) const
2730{
2731#ifndef NDEBUG
2732  int n = modelPtr_->numberColumns();
2733  if (colNumber<0||colNumber>=n) {
2734    indexError(colNumber,"isBinary");
2735  }
2736#endif
2737  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) {
2738    return false;
2739  } else {
2740    const double * cu = getColUpper();
2741    const double * cl = getColLower();
2742    if ((cu[colNumber]== 1 || cu[colNumber]== 0) && 
2743        (cl[colNumber]== 0 || cl[colNumber]==1))
2744      return true;
2745    else 
2746      return false;
2747  }
2748}
2749//-----------------------------------------------------------------------------
2750bool 
2751OsiClpSolverInterface::isInteger(int colNumber) const
2752{
2753#ifndef NDEBUG
2754  int n = modelPtr_->numberColumns();
2755  if (colNumber<0||colNumber>=n) {
2756    indexError(colNumber,"isInteger");
2757  }
2758#endif
2759  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) 
2760    return false;
2761  else
2762    return true;
2763}
2764//-----------------------------------------------------------------------------
2765bool 
2766OsiClpSolverInterface::isIntegerNonBinary(int colNumber) const
2767{
2768#ifndef NDEBUG
2769  int n = modelPtr_->numberColumns();
2770  if (colNumber<0||colNumber>=n) {
2771    indexError(colNumber,"isIntegerNonBinary");
2772  }
2773#endif
2774  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) {
2775    return false;
2776  } else {
2777    return !isBinary(colNumber);
2778  }
2779}
2780//-----------------------------------------------------------------------------
2781bool 
2782OsiClpSolverInterface::isFreeBinary(int colNumber) const
2783{
2784#ifndef NDEBUG
2785  int n = modelPtr_->numberColumns();
2786  if (colNumber<0||colNumber>=n) {
2787    indexError(colNumber,"isFreeBinary");
2788  }
2789#endif
2790  if ( integerInformation_==NULL || integerInformation_[colNumber]==0 ) {
2791    return false;
2792  } else {
2793    const double * cu = getColUpper();
2794    const double * cl = getColLower();
2795    if ((cu[colNumber]== 1) && (cl[colNumber]== 0))
2796      return true;
2797    else 
2798      return false;
2799  }
2800}
2801/*  Return array of column length
2802    0 - continuous
2803    1 - binary (may get fixed later)
2804    2 - general integer (may get fixed later)
2805*/
2806const char * 
2807OsiClpSolverInterface::getColType(bool refresh) const
2808{
2809  if (!columnType_||refresh) {
2810    const int numCols = getNumCols() ;
2811    if (!columnType_)
2812      columnType_ = new char [numCols];
2813    if ( integerInformation_==NULL ) {
2814      memset(columnType_,0,numCols);
2815    } else {
2816      const double * cu = getColUpper();
2817      const double * cl = getColLower();
2818      for (int i = 0 ; i < numCols ; ++i) {
2819        if (integerInformation_[i]) {
2820          if ((cu[i]== 1 || cu[i]== 0) && 
2821              (cl[i]== 0 || cl[i]==1))
2822            columnType_[i]=1;
2823          else
2824            columnType_[i]=2;
2825        } else {
2826          columnType_[i]=0;
2827        }
2828      }
2829    }
2830  }
2831  return columnType_;
2832}
2833
2834/* Return true if column is integer but does not have to
2835   be declared as such.
2836   Note: This function returns true if the the column
2837   is binary or a general integer.
2838*/
2839bool 
2840OsiClpSolverInterface::isOptionalInteger(int colNumber) const
2841{
2842#ifndef NDEBUG
2843  int n = modelPtr_->numberColumns();
2844  if (colNumber<0||colNumber>=n) {
2845    indexError(colNumber,"isInteger");
2846  }
2847#endif
2848  if ( integerInformation_==NULL || integerInformation_[colNumber]!=2 ) 
2849    return false;
2850  else
2851    return true;
2852}
2853/* Set the index-th variable to be an optional integer variable */
2854void 
2855OsiClpSolverInterface::setOptionalInteger(int index)
2856{
2857  if (!integerInformation_) {
2858    integerInformation_ = new char[modelPtr_->numberColumns()];
2859    CoinFillN ( integerInformation_, modelPtr_->numberColumns(),static_cast<char> (0));
2860  }
2861#ifndef NDEBUG
2862  int n = modelPtr_->numberColumns();
2863  if (index<0||index>=n) {
2864    indexError(index,"setInteger");
2865  }
2866#endif
2867  integerInformation_[index]=2;
2868  modelPtr_->setInteger(index);
2869}
2870
2871//------------------------------------------------------------------
2872// Row and column copies of the matrix ...
2873//------------------------------------------------------------------
2874const CoinPackedMatrix * OsiClpSolverInterface::getMatrixByRow() const
2875{
2876  if ( matrixByRow_ == NULL ||
2877       matrixByRow_->getNumElements() != 
2878       modelPtr_->clpMatrix()->getNumElements() ) {
2879    delete matrixByRow_;
2880    matrixByRow_ = new CoinPackedMatrix();
2881    matrixByRow_->setExtraGap(0.0);
2882    matrixByRow_->setExtraMajor(0.0);
2883    matrixByRow_->reverseOrderedCopyOf(*modelPtr_->matrix());
2884    //matrixByRow_->removeGaps();
2885#if 0
2886    CoinPackedMatrix back;
2887    std::cout<<"start check"<<std::endl;
2888    back.reverseOrderedCopyOf(*matrixByRow_);
2889    modelPtr_->matrix()->isEquivalent2(back);
2890    std::cout<<"stop check"<<std::endl;
2891#endif
2892  }
2893  assert (matrixByRow_->getNumElements()==modelPtr_->clpMatrix()->getNumElements());
2894  return matrixByRow_;
2895}
2896
2897const CoinPackedMatrix * OsiClpSolverInterface::getMatrixByCol() const
2898{
2899  return modelPtr_->matrix();
2900}
2901 
2902// Get pointer to mutable column-wise copy of matrix (returns NULL if not meaningful)
2903CoinPackedMatrix * 
2904OsiClpSolverInterface::getMutableMatrixByCol() const 
2905{
2906  ClpPackedMatrix * matrix = dynamic_cast<ClpPackedMatrix *>(modelPtr_->matrix_) ;
2907  if (matrix)
2908    return matrix->getPackedMatrix();
2909  else
2910    return NULL;
2911}
2912
2913//------------------------------------------------------------------
2914std::vector<double*> OsiClpSolverInterface::getDualRays(int /*maxNumRays*/,
2915                                                        bool fullRay) const
2916{
2917  return std::vector<double*>(1, modelPtr_->infeasibilityRay(fullRay));
2918}
2919//------------------------------------------------------------------
2920std::vector<double*> OsiClpSolverInterface::getPrimalRays(int /*maxNumRays*/) const
2921{
2922  return std::vector<double*>(1, modelPtr_->unboundedRay());
2923}
2924//#############################################################################
2925void
2926OsiClpSolverInterface::setContinuous(int index)
2927{
2928
2929  if (integerInformation_) {
2930#ifndef NDEBUG
2931    int n = modelPtr_->numberColumns();
2932    if (index<0||index>=n) {
2933      indexError(index,"setContinuous");
2934    }
2935#endif
2936    integerInformation_[index]=0;
2937  }
2938  modelPtr_->setContinuous(index);
2939}
2940//-----------------------------------------------------------------------------
2941void
2942OsiClpSolverInterface::setInteger(int index)
2943{
2944  if (!integerInformation_) {
2945    integerInformation_ = new char[modelPtr_->numberColumns()];
2946    CoinFillN ( integerInformation_, modelPtr_->numberColumns(),static_cast<char> (0));
2947  }
2948#ifndef NDEBUG
2949  int n = modelPtr_->numberColumns();
2950  if (index<0||index>=n) {
2951    indexError(index,"setInteger");
2952  }
2953#endif
2954  integerInformation_[index]=1;
2955  modelPtr_->setInteger(index);
2956}
2957//-----------------------------------------------------------------------------
2958void
2959OsiClpSolverInterface::setContinuous(const int* indices, int len)
2960{
2961  if (integerInformation_) {
2962#ifndef NDEBUG
2963    int n = modelPtr_->numberColumns();
2964#endif
2965    int i;
2966    for (i=0; i<len;i++) {
2967      int colNumber = indices[i];
2968#ifndef NDEBUG
2969      if (colNumber<0||colNumber>=n) {
2970        indexError(colNumber,"setContinuous");
2971      }
2972#endif
2973      integerInformation_[colNumber]=0;
2974      modelPtr_->setContinuous(colNumber);
2975    }
2976  }
2977}
2978//-----------------------------------------------------------------------------
2979void
2980OsiClpSolverInterface::setInteger(const int* indices, int len)
2981{
2982  if (!integerInformation_) {
2983    integerInformation_ = new char[modelPtr_->numberColumns()];
2984    CoinFillN ( integerInformation_, modelPtr_->numberColumns(),static_cast<char> (0));
2985  }
2986#ifndef NDEBUG
2987  int n = modelPtr_->numberColumns();
2988#endif
2989  int i;
2990  for (i=0; i<len;i++) {
2991    int colNumber = indices[i];
2992#ifndef NDEBUG
2993    if (colNumber<0||colNumber>=n) {
2994      indexError(colNumber,"setInteger");
2995    }
2996#endif
2997    integerInformation_[colNumber]=1;
2998    modelPtr_->setInteger(colNumber);
2999  }
3000}
3001/* Set the objective coefficients for all columns
3002    array [getNumCols()] is an array of values for the objective.
3003    This defaults to a series of set operations and is here for speed.
3004*/
3005void 
3006OsiClpSolverInterface::setObjective(const double * array)
3007{
3008  // Say can't gurantee optimal basis etc
3009  lastAlgorithm_=999;
3010  modelPtr_->whatsChanged_ &= (0xffff&~64);
3011  int n = modelPtr_->numberColumns() ;
3012  if (fakeMinInSimplex_) {
3013    std::transform(array,array+n,
3014                   modelPtr_->objective(),std::negate<double>()) ;
3015  } else {
3016    CoinMemcpyN(array,n,modelPtr_->objective());
3017  }
3018}
3019/* Set the lower bounds for all columns
3020    array [getNumCols()] is an array of values for the objective.
3021    This defaults to a series of set operations and is here for speed.
3022*/
3023void 
3024OsiClpSolverInterface::setColLower(const double * array)
3025{
3026  // Say can't gurantee optimal basis etc
3027  lastAlgorithm_=999;
3028  modelPtr_->whatsChanged_ &= (0x1ffff&128);
3029  CoinMemcpyN(array,modelPtr_->numberColumns(),
3030                    modelPtr_->columnLower());
3031}
3032/* Set the upper bounds for all columns
3033    array [getNumCols()] is an array of values for the objective.
3034    This defaults to a series of set operations and is here for speed.
3035*/
3036void 
3037OsiClpSolverInterface::setColUpper(const double * array)
3038{
3039  // Say can't gurantee optimal basis etc
3040  lastAlgorithm_=999;
3041  modelPtr_->whatsChanged_ &= (0x1ffff&256);
3042  CoinMemcpyN(array,modelPtr_->numberColumns(),
3043                    modelPtr_->columnUpper());
3044}
3045//-----------------------------------------------------------------------------
3046void OsiClpSolverInterface::setColSolution(const double * cs) 
3047{
3048  // Say can't gurantee optimal basis etc
3049  lastAlgorithm_=999;
3050  CoinDisjointCopyN(cs,modelPtr_->numberColumns(),
3051                    modelPtr_->primalColumnSolution());
3052  if (modelPtr_->solveType()==2) {
3053    // directly into code as well
3054    CoinDisjointCopyN(cs,modelPtr_->numberColumns(),
3055                      modelPtr_->solutionRegion(1));
3056  }
3057  // compute row activity
3058  memset(modelPtr_->primalRowSolution(),0,modelPtr_->numberRows()*sizeof(double));
3059  modelPtr_->times(1.0,modelPtr_->primalColumnSolution(),modelPtr_->primalRowSolution());
3060}
3061//-----------------------------------------------------------------------------
3062void OsiClpSolverInterface::setRowPrice(const double * rs) 
3063{
3064  CoinDisjointCopyN(rs,modelPtr_->numberRows(),
3065                    modelPtr_->dualRowSolution());
3066  if (modelPtr_->solveType()==2) {
3067    // directly into code as well (? sign )
3068    CoinDisjointCopyN(rs,modelPtr_->numberRows(),
3069                      modelPtr_->djRegion(0));
3070  }
3071  // compute reduced costs
3072  memcpy(modelPtr_->dualColumnSolution(),modelPtr_->objective(),
3073         modelPtr_->numberColumns()*sizeof(double));
3074  modelPtr_->transposeTimes(-1.0,modelPtr_->dualRowSolution(),modelPtr_->dualColumnSolution());
3075}
3076
3077//#############################################################################
3078// Problem modifying methods (matrix)
3079//#############################################################################
3080void 
3081OsiClpSolverInterface::addCol(const CoinPackedVectorBase& vec,
3082                              const double collb, const double colub,   
3083                              const double obj)
3084{
3085  int numberColumns = modelPtr_->numberColumns();
3086  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3087  modelPtr_->resize(modelPtr_->numberRows(),numberColumns+1);
3088  linearObjective_ = modelPtr_->objective();
3089  basis_.resize(modelPtr_->numberRows(),numberColumns+1);
3090  setColBounds(numberColumns,collb,colub);
3091  setObjCoeff(numberColumns,obj);
3092  if (!modelPtr_->clpMatrix())
3093    modelPtr_->createEmptyMatrix();
3094  modelPtr_->matrix()->appendCol(vec);
3095  if (integerInformation_) {
3096    char * temp = new char[numberColumns+1];
3097    CoinMemcpyN(integerInformation_,numberColumns,temp);
3098    delete [] integerInformation_;
3099    integerInformation_ = temp;
3100    integerInformation_[numberColumns]=0;
3101  }
3102  freeCachedResults();
3103}
3104//-----------------------------------------------------------------------------
3105/* Add a column (primal variable) to the problem. */
3106void 
3107OsiClpSolverInterface::addCol(int numberElements, const int * rows, const double * elements,
3108                           const double collb, const double colub,   
3109                           const double obj) 
3110{
3111  CoinPackedVector column(numberElements, rows, elements);
3112  addCol(column,collb,colub,obj);
3113}
3114// Add a named column (primal variable) to the problem.
3115void 
3116OsiClpSolverInterface::addCol(const CoinPackedVectorBase& vec,
3117                      const double collb, const double colub,   
3118                      const double obj, std::string name) 
3119{
3120  int ndx = getNumCols() ;
3121  addCol(vec,collb,colub,obj) ;
3122  setColName(ndx,name) ;
3123}
3124//-----------------------------------------------------------------------------
3125void 
3126OsiClpSolverInterface::addCols(const int numcols,
3127                               const CoinPackedVectorBase * const * cols,
3128                               const double* collb, const double* colub,   
3129                               const double* obj)
3130{
3131  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3132  int numberColumns = modelPtr_->numberColumns();
3133  modelPtr_->resize(modelPtr_->numberRows(),numberColumns+numcols);
3134  linearObjective_ = modelPtr_->objective();
3135  basis_.resize(modelPtr_->numberRows(),numberColumns+numcols);
3136  double * lower = modelPtr_->columnLower()+numberColumns;
3137  double * upper = modelPtr_->columnUpper()+numberColumns;
3138  double * objective = modelPtr_->objective()+numberColumns;
3139  int iCol;
3140  if (collb) {
3141    for (iCol = 0; iCol < numcols; iCol++) {
3142      lower[iCol]= forceIntoRange(collb[iCol], -OsiClpInfinity, OsiClpInfinity);
3143      if (lower[iCol]<-1.0e27)
3144        lower[iCol]=-COIN_DBL_MAX;
3145    }
3146  } else {
3147    CoinFillN ( lower, numcols,0.0);
3148  }
3149  if (colub) {
3150    for (iCol = 0; iCol < numcols; iCol++) {
3151      upper[iCol]= forceIntoRange(colub[iCol], -OsiClpInfinity, OsiClpInfinity);
3152      if (upper[iCol]>1.0e27)
3153        upper[iCol]=COIN_DBL_MAX;
3154    }
3155  } else {
3156    CoinFillN ( upper, numcols,COIN_DBL_MAX);
3157  }
3158  if (obj) {
3159    for (iCol = 0; iCol < numcols; iCol++) {
3160      objective[iCol] = obj[iCol];
3161    }
3162  } else {
3163    CoinFillN ( objective, numcols,0.0);
3164  }
3165  if (!modelPtr_->clpMatrix())
3166    modelPtr_->createEmptyMatrix();
3167  modelPtr_->matrix()->appendCols(numcols,cols);
3168  if (integerInformation_) {
3169    char * temp = new char[numberColumns+numcols];
3170    CoinMemcpyN(integerInformation_,numberColumns,temp);
3171    delete [] integerInformation_;
3172    integerInformation_ = temp;
3173    for (iCol = 0; iCol < numcols; iCol++) 
3174      integerInformation_[numberColumns+iCol]=0;
3175  }
3176  freeCachedResults();
3177}
3178void 
3179OsiClpSolverInterface::addCols(const int numcols,
3180                               const int * columnStarts, const int * rows, const double * elements,
3181                               const double* collb, const double* colub,   
3182                               const double* obj)
3183{
3184  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3185  int numberColumns = modelPtr_->numberColumns();
3186  modelPtr_->resize(modelPtr_->numberRows(),numberColumns+numcols);
3187  linearObjective_ = modelPtr_->objective();
3188  basis_.resize(modelPtr_->numberRows(),numberColumns+numcols);
3189  double * lower = modelPtr_->columnLower()+numberColumns;
3190  double * upper = modelPtr_->columnUpper()+numberColumns;
3191  double * objective = modelPtr_->objective()+numberColumns;
3192  int iCol;
3193  if (collb) {
3194    for (iCol = 0; iCol < numcols; iCol++) {
3195      lower[iCol]= forceIntoRange(collb[iCol], -OsiClpInfinity, OsiClpInfinity);
3196      if (lower[iCol]<-1.0e27)
3197        lower[iCol]=-COIN_DBL_MAX;
3198    }
3199  } else {
3200    CoinFillN ( lower, numcols,0.0);
3201  }
3202  if (colub) {
3203    for (iCol = 0; iCol < numcols; iCol++) {
3204      upper[iCol]= forceIntoRange(colub[iCol], -OsiClpInfinity, OsiClpInfinity);
3205      if (upper[iCol]>1.0e27)
3206        upper[iCol]=COIN_DBL_MAX;
3207    }
3208  } else {
3209    CoinFillN ( upper, numcols,COIN_DBL_MAX);
3210  }
3211  if (obj) {
3212    for (iCol = 0; iCol < numcols; iCol++) {
3213      objective[iCol] = obj[iCol];
3214    }
3215  } else {
3216    CoinFillN ( objective, numcols,0.0);
3217  }
3218  if (!modelPtr_->clpMatrix())
3219    modelPtr_->createEmptyMatrix();
3220  modelPtr_->matrix()->appendCols(numcols,columnStarts,rows,elements);
3221  if (integerInformation_) {
3222    char * temp = new char[numberColumns+numcols];
3223    CoinMemcpyN(integerInformation_,numberColumns,temp);
3224    delete [] integerInformation_;
3225    integerInformation_ = temp;
3226    for (iCol = 0; iCol < numcols; iCol++) 
3227      integerInformation_[numberColumns+iCol]=0;
3228  }
3229  freeCachedResults();
3230}
3231//-----------------------------------------------------------------------------
3232void 
3233OsiClpSolverInterface::deleteCols(const int num, const int * columnIndices)
3234{
3235  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|8|64|128|256));
3236  findIntegers(false);
3237  deleteBranchingInfo(num,columnIndices);
3238  modelPtr_->deleteColumns(num,columnIndices);
3239  int nameDiscipline;
3240  getIntParam(OsiNameDiscipline,nameDiscipline) ;
3241  if (num&&nameDiscipline) {
3242    // Very clumsy (and inefficient) - need to sort and then go backwards in ? chunks
3243    int * indices = CoinCopyOfArray(columnIndices,num);
3244    std::sort(indices,indices+num);
3245    int num2=num;
3246    while(num2) {
3247      int next = indices[num2-1];
3248      int firstDelete = num2-1;
3249      int i;
3250      for (i=num2-2;i>=0;i--) {
3251        if (indices[i]+1==next) {
3252          next --;
3253          firstDelete=i;
3254        } else {
3255          break;
3256        }
3257      }
3258      OsiSolverInterface::deleteColNames(indices[firstDelete],num2-firstDelete);
3259      num2 = firstDelete;
3260      assert (num2>=0);
3261    }
3262    delete [] indices;
3263  }
3264  // synchronize integers (again)
3265  if (integerInformation_) {
3266    int numberColumns = modelPtr_->numberColumns();
3267    for (int i=0;i<numberColumns;i++) {
3268      if (modelPtr_->isInteger(i))
3269        integerInformation_[i]=1;
3270      else
3271        integerInformation_[i]=0;
3272    }
3273  }
3274  basis_.deleteColumns(num,columnIndices);
3275  linearObjective_ = modelPtr_->objective();
3276  freeCachedResults();
3277}
3278//-----------------------------------------------------------------------------
3279void 
3280OsiClpSolverInterface::addRow(const CoinPackedVectorBase& vec,
3281                              const double rowlb, const double rowub)
3282{
3283  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3284  freeCachedResults0();
3285  int numberRows = modelPtr_->numberRows();
3286  modelPtr_->resize(numberRows+1,modelPtr_->numberColumns());
3287  basis_.resize(numberRows+1,modelPtr_->numberColumns());
3288  setRowBounds(numberRows,rowlb,rowub);
3289  if (!modelPtr_->clpMatrix())
3290    modelPtr_->createEmptyMatrix();
3291  modelPtr_->matrix()->appendRow(vec);
3292  freeCachedResults1();
3293}
3294//-----------------------------------------------------------------------------
3295void OsiClpSolverInterface::addRow(const CoinPackedVectorBase& vec,
3296                                const double rowlb, const double rowub,
3297                                std::string name)
3298{
3299  int ndx = getNumRows() ;
3300  addRow(vec,rowlb,rowub) ;
3301  setRowName(ndx,name) ;
3302}
3303//-----------------------------------------------------------------------------
3304void 
3305OsiClpSolverInterface::addRow(const CoinPackedVectorBase& vec,
3306                              const char rowsen, const double rowrhs,   
3307                              const double rowrng)
3308{
3309  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3310  freeCachedResults0();
3311  int numberRows = modelPtr_->numberRows();
3312  modelPtr_->resize(numberRows+1,modelPtr_->numberColumns());
3313  basis_.resize(numberRows+1,modelPtr_->numberColumns());
3314  double rowlb = 0, rowub = 0;
3315  convertSenseToBound(rowsen, rowrhs, rowrng, rowlb, rowub);
3316  setRowBounds(numberRows,rowlb,rowub);
3317  if (!modelPtr_->clpMatrix())
3318    modelPtr_->createEmptyMatrix();
3319  modelPtr_->matrix()->appendRow(vec);
3320  freeCachedResults1();
3321}
3322//-----------------------------------------------------------------------------
3323void 
3324OsiClpSolverInterface::addRow(int numberElements, const int * columns, const double * elements,
3325                           const double rowlb, const double rowub) 
3326{
3327  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3328  freeCachedResults0();
3329  int numberRows = modelPtr_->numberRows();
3330  modelPtr_->resize(numberRows+1,modelPtr_->numberColumns());
3331  basis_.resize(numberRows+1,modelPtr_->numberColumns());
3332  setRowBounds(numberRows,rowlb,rowub);
3333  if (!modelPtr_->clpMatrix())
3334    modelPtr_->createEmptyMatrix();
3335  modelPtr_->matrix()->appendRow(numberElements, columns, elements);
3336  CoinBigIndex starts[2];
3337  starts[0]=0;
3338  starts[1]=numberElements;
3339  redoScaleFactors( 1,starts, columns, elements);
3340  freeCachedResults1();
3341}
3342//-----------------------------------------------------------------------------
3343void 
3344OsiClpSolverInterface::addRows(const int numrows,
3345                               const CoinPackedVectorBase * const * rows,
3346                               const double* rowlb, const double* rowub)
3347{
3348  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3349  freeCachedResults0();
3350  int numberRows = modelPtr_->numberRows();
3351  modelPtr_->resize(numberRows+numrows,modelPtr_->numberColumns());
3352  basis_.resize(numberRows+numrows,modelPtr_->numberColumns());
3353  double * lower = modelPtr_->rowLower()+numberRows;
3354  double * upper = modelPtr_->rowUpper()+numberRows;
3355  int iRow;
3356  for (iRow = 0; iRow < numrows; iRow++) {
3357    if (rowlb) 
3358      lower[iRow]= forceIntoRange(rowlb[iRow], -OsiClpInfinity, OsiClpInfinity);
3359    else 
3360      lower[iRow]=-OsiClpInfinity;
3361    if (rowub) 
3362      upper[iRow]= forceIntoRange(rowub[iRow], -OsiClpInfinity, OsiClpInfinity);
3363    else 
3364      upper[iRow]=OsiClpInfinity;
3365    if (lower[iRow]<-1.0e27)
3366      lower[iRow]=-COIN_DBL_MAX;
3367    if (upper[iRow]>1.0e27)
3368      upper[iRow]=COIN_DBL_MAX;
3369  }
3370  if (!modelPtr_->clpMatrix())
3371    modelPtr_->createEmptyMatrix();
3372  modelPtr_->matrix()->appendRows(numrows,rows);
3373  freeCachedResults1();
3374}
3375//-----------------------------------------------------------------------------
3376void 
3377OsiClpSolverInterface::addRows(const int numrows,
3378                               const CoinPackedVectorBase * const * rows,
3379                               const char* rowsen, const double* rowrhs,   
3380                               const double* rowrng)
3381{
3382  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3383  freeCachedResults0();
3384  int numberRows = modelPtr_->numberRows();
3385  modelPtr_->resize(numberRows+numrows,modelPtr_->numberColumns());
3386  basis_.resize(numberRows+numrows,modelPtr_->numberColumns());
3387  double * lower = modelPtr_->rowLower()+numberRows;
3388  double * upper = modelPtr_->rowUpper()+numberRows;
3389  int iRow;
3390  for (iRow = 0; iRow < numrows; iRow++) {
3391    double rowlb = 0, rowub = 0;
3392    convertSenseToBound(rowsen[iRow], rowrhs[iRow], rowrng[iRow], 
3393                        rowlb, rowub);
3394    lower[iRow]= forceIntoRange(rowlb, -OsiClpInfinity, OsiClpInfinity);
3395    upper[iRow]= forceIntoRange(rowub, -OsiClpInfinity, OsiClpInfinity);
3396    if (lower[iRow]<-1.0e27)
3397      lower[iRow]=-COIN_DBL_MAX;
3398    if (upper[iRow]>1.0e27)
3399      upper[iRow]=COIN_DBL_MAX;
3400  }
3401  if (!modelPtr_->clpMatrix())
3402    modelPtr_->createEmptyMatrix();
3403  modelPtr_->matrix()->appendRows(numrows,rows);
3404  freeCachedResults1();
3405}
3406void 
3407OsiClpSolverInterface::addRows(const int numrows,
3408                               const int * rowStarts, const int * columns, const double * element,
3409                               const double* rowlb, const double* rowub)
3410{
3411  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3412  freeCachedResults0();
3413  int numberRows = modelPtr_->numberRows();
3414  modelPtr_->resize(numberRows+numrows,modelPtr_->numberColumns());
3415  basis_.resize(numberRows+numrows,modelPtr_->numberColumns());
3416  double * lower = modelPtr_->rowLower()+numberRows;
3417  double * upper = modelPtr_->rowUpper()+numberRows;
3418  int iRow;
3419  for (iRow = 0; iRow < numrows; iRow++) {
3420    if (rowlb) 
3421      lower[iRow]= forceIntoRange(rowlb[iRow], -OsiClpInfinity, OsiClpInfinity);
3422    else 
3423      lower[iRow]=-OsiClpInfinity;
3424    if (rowub) 
3425      upper[iRow]= forceIntoRange(rowub[iRow], -OsiClpInfinity, OsiClpInfinity);
3426    else 
3427      upper[iRow]=OsiClpInfinity;
3428    if (lower[iRow]<-1.0e27)
3429      lower[iRow]=-COIN_DBL_MAX;
3430    if (upper[iRow]>1.0e27)
3431      upper[iRow]=COIN_DBL_MAX;
3432  }
3433  if (!modelPtr_->clpMatrix())
3434    modelPtr_->createEmptyMatrix();
3435  modelPtr_->matrix()->appendRows(numrows,rowStarts,columns,element);
3436  redoScaleFactors( numrows,rowStarts, columns, element);
3437  freeCachedResults1();
3438}
3439//-----------------------------------------------------------------------------
3440void 
3441OsiClpSolverInterface::deleteRows(const int num, const int * rowIndices)
3442{
3443  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
3444  // will still be optimal if all rows basic
3445  bool allBasic=true;
3446  int numBasis = basis_.getNumArtificial();
3447  for (int i=0;i<num;i++) {
3448    int iRow = rowIndices[i];
3449    if (iRow<numBasis) {
3450      if (basis_.getArtifStatus(iRow)!=CoinWarmStartBasis::basic) {
3451        allBasic=false;
3452        break;
3453      }
3454    }
3455  }
3456  int saveAlgorithm = allBasic ? lastAlgorithm_ : 999;
3457  modelPtr_->deleteRows(num,rowIndices);
3458  int nameDiscipline;
3459  getIntParam(OsiNameDiscipline,nameDiscipline) ;
3460  if (num&&nameDiscipline) {
3461    // Very clumsy (and inefficient) - need to sort and then go backwards in ? chunks
3462    int * indices = CoinCopyOfArray(rowIndices,num);
3463    std::sort(indices,indices+num);
3464    int num2=num;
3465    while(num2) {
3466      int next = indices[num2-1];
3467      int firstDelete = num2-1;
3468      int i;
3469      for (i=num2-2;i>=0;i--) {
3470        if (indices[i]+1==next) {
3471          next --;
3472          firstDelete=i;
3473        } else {
3474          break;
3475        }
3476      }
3477      OsiSolverInterface::deleteRowNames(indices[firstDelete],num2-firstDelete);
3478      num2 = firstDelete;
3479      assert (num2>=0);
3480    }
3481    delete [] indices;
3482  }
3483  basis_.deleteRows(num,rowIndices);
3484  CoinPackedMatrix * saveRowCopy = matrixByRow_;
3485  matrixByRow_=NULL;
3486  freeCachedResults();
3487  modelPtr_->setNewRowCopy(NULL);
3488  delete modelPtr_->scaledMatrix_;
3489  modelPtr_->scaledMatrix_=NULL;
3490  if (saveRowCopy) {
3491    matrixByRow_=saveRowCopy;
3492    matrixByRow_->deleteRows(num,rowIndices);
3493    if (matrixByRow_->getNumElements()!=modelPtr_->clpMatrix()->getNumElements()) {
3494      delete matrixByRow_; // odd type matrix
3495      matrixByRow_=NULL;
3496    }
3497  }
3498  lastAlgorithm_ = saveAlgorithm;
3499  if ((specialOptions_&131072)!=0) 
3500    lastNumberRows_=modelPtr_->numberRows();
3501}
3502
3503//#############################################################################
3504// Methods to input a problem
3505//#############################################################################
3506
3507void
3508OsiClpSolverInterface::loadProblem(const CoinPackedMatrix& matrix,
3509                                   const double* collb, const double* colub,   
3510                                   const double* obj,
3511                                   const double* rowlb, const double* rowub)
3512{
3513  modelPtr_->whatsChanged_ = 0;
3514  // Get rid of integer information (modelPtr will get rid of its copy)
3515  delete [] integerInformation_;
3516  integerInformation_=NULL;
3517  modelPtr_->loadProblem(matrix, collb, colub, obj, rowlb, rowub);
3518  linearObjective_ = modelPtr_->objective();
3519  freeCachedResults();
3520  basis_=CoinWarmStartBasis();
3521  if (ws_) {
3522     delete ws_;
3523     ws_ = 0;
3524  }
3525}
3526
3527//-----------------------------------------------------------------------------
3528
3529/*
3530  Expose the method that takes ClpMatrixBase. User request. Can't hurt, given
3531  the number of non-OSI methods already here.
3532*/
3533void OsiClpSolverInterface::loadProblem (const ClpMatrixBase& matrix,
3534                                         const double* collb,
3535                                         const double* colub,   
3536                                         const double* obj,
3537                                         const double* rowlb,
3538                                         const double* rowub)
3539{
3540  modelPtr_->whatsChanged_ = 0;
3541  // Get rid of integer information (modelPtr will get rid of its copy)
3542  delete [] integerInformation_;
3543  integerInformation_=NULL;
3544  modelPtr_->loadProblem(matrix,collb,colub,obj,rowlb,rowub);
3545  linearObjective_ = modelPtr_->objective();
3546  freeCachedResults();
3547  basis_=CoinWarmStartBasis();
3548  if (ws_) {
3549     delete ws_;
3550     ws_ = 0;
3551  }
3552}
3553
3554//-----------------------------------------------------------------------------
3555
3556void
3557OsiClpSolverInterface::assignProblem(CoinPackedMatrix*& matrix,
3558                                     double*& collb, double*& colub,
3559                                     double*& obj,
3560                                     double*& rowlb, double*& rowub)
3561{
3562  modelPtr_->whatsChanged_ = 0;
3563  // Get rid of integer information (modelPtr will get rid of its copy)
3564  loadProblem(*matrix, collb, colub, obj, rowlb, rowub);
3565  delete matrix;   matrix = NULL;
3566  delete[] collb;  collb = NULL;
3567  delete[] colub;  colub = NULL;
3568  delete[] obj;    obj = NULL;
3569  delete[] rowlb;  rowlb = NULL;
3570  delete[] rowub;  rowub = NULL;
3571}
3572
3573//-----------------------------------------------------------------------------
3574
3575void
3576OsiClpSolverInterface::loadProblem(const CoinPackedMatrix& matrix,
3577                                   const double* collb, const double* colub,
3578                                   const double* obj,
3579                                   const char* rowsen, const double* rowrhs,   
3580                                   const double* rowrng)
3581{
3582  modelPtr_->whatsChanged_ = 0;
3583  // Get rid of integer information (modelPtr will get rid of its copy)
3584  // assert( rowsen != NULL );
3585  // assert( rowrhs != NULL );
3586  // If any of Rhs NULLs then create arrays
3587  int numrows = matrix.getNumRows();
3588  const char * rowsenUse = rowsen;
3589  if (!rowsen) {
3590    char * rowsen = new char [numrows];
3591    for (int i=0;i<numrows;i++)
3592      rowsen[i]='G';
3593    rowsenUse = rowsen;
3594  } 
3595  const double * rowrhsUse = rowrhs;
3596  if (!rowrhs) {
3597    double * rowrhs = new double [numrows];
3598    for (int i=0;i<numrows;i++)
3599      rowrhs[i]=0.0;
3600    rowrhsUse = rowrhs;
3601  }
3602  const double * rowrngUse = rowrng;
3603  if (!rowrng) {
3604    double * rowrng = new double [numrows];
3605    for (int i=0;i<numrows;i++)
3606      rowrng[i]=0.0;
3607    rowrngUse = rowrng;
3608  }
3609  double * rowlb = new double[numrows];
3610  double * rowub = new double[numrows];
3611  for (int i = numrows-1; i >= 0; --i) {   
3612    convertSenseToBound(rowsenUse[i],rowrhsUse[i],rowrngUse[i],rowlb[i],rowub[i]);
3613  }
3614  if (rowsen!=rowsenUse)
3615    delete [] rowsenUse;
3616  if (rowrhs!=rowrhsUse)
3617    delete [] rowrhsUse;
3618  if (rowrng!=rowrngUse)
3619    delete [] rowrngUse;
3620  loadProblem(matrix, collb, colub, obj, rowlb, rowub);
3621  delete [] rowlb;
3622  delete [] rowub;
3623}
3624
3625//-----------------------------------------------------------------------------
3626
3627void
3628OsiClpSolverInterface::assignProblem(CoinPackedMatrix*& matrix,
3629                                     double*& collb, double*& colub,
3630                                     double*& obj,
3631                                     char*& rowsen, double*& rowrhs,
3632                                     double*& rowrng)
3633{
3634  modelPtr_->whatsChanged_ = 0;
3635  // Get rid of integer information (modelPtr will get rid of its copy)
3636  loadProblem(*matrix, collb, colub, obj, rowsen, rowrhs, rowrng);
3637  delete matrix;   matrix = NULL;
3638  delete[] collb;  collb = NULL;
3639  delete[] colub;  colub = NULL;
3640  delete[] obj;    obj = NULL;
3641  delete[] rowsen; rowsen = NULL;
3642  delete[] rowrhs; rowrhs = NULL;
3643  delete[] rowrng; rowrng = NULL;
3644}
3645
3646//-----------------------------------------------------------------------------
3647
3648void
3649OsiClpSolverInterface::loadProblem(const int numcols, const int numrows,
3650                                   const CoinBigIndex * start, const int* index,
3651                                   const double* value,
3652                                   const double* collb, const double* colub,
3653                                   const double* obj,
3654                                   const double* rowlb, const double* rowub)
3655{
3656  modelPtr_->whatsChanged_ = 0;
3657  // Get rid of integer information (modelPtr will get rid of its copy)
3658  delete [] integerInformation_;
3659  integerInformation_=NULL;
3660  modelPtr_->loadProblem(numcols, numrows, start,  index,
3661            value, collb, colub, obj,
3662            rowlb,  rowub);
3663  linearObjective_ = modelPtr_->objective();
3664  freeCachedResults();
3665  basis_=CoinWarmStartBasis();
3666  if (ws_) {
3667     delete ws_;
3668     ws_ = 0;
3669  }
3670}
3671//-----------------------------------------------------------------------------
3672
3673void
3674OsiClpSolverInterface::loadProblem(const int numcols, const int numrows,
3675                                   const CoinBigIndex * start, const int* index,
3676                                   const double* value,
3677                                   const double* collb, const double* colub,
3678                                   const double* obj,
3679                                   const char* rowsen, const double* rowrhs,   
3680                                   const double* rowrng)
3681{
3682  modelPtr_->whatsChanged_ = 0;
3683  // Get rid of integer information (modelPtr will get rid of its copy)
3684  // If any of Rhs NULLs then create arrays
3685  const char * rowsenUse = rowsen;
3686  if (!rowsen) {
3687    char * rowsen = new char [numrows];
3688    for (int i=0;i<numrows;i++)
3689      rowsen[i]='G';
3690    rowsenUse = rowsen;
3691  } 
3692  const double * rowrhsUse = rowrhs;
3693  if (!rowrhs) {
3694    double * rowrhs = new double [numrows];
3695    for (int i=0;i<numrows;i++)
3696      rowrhs[i]=0.0;
3697    rowrhsUse = rowrhs;
3698  }
3699  const double * rowrngUse = rowrng;
3700  if (!rowrng) {
3701    double * rowrng = new double [numrows];
3702    for (int i=0;i<numrows;i++)
3703      rowrng[i]=0.0;
3704    rowrngUse = rowrng;
3705  }
3706  double * rowlb = new double[numrows];
3707  double * rowub = new double[numrows];
3708  for (int i = numrows-1; i >= 0; --i) {   
3709    convertSenseToBound(rowsenUse[i],rowrhsUse[i],rowrngUse[i],rowlb[i],rowub[i]);
3710  }
3711  if (rowsen!=rowsenUse)
3712    delete [] rowsenUse;
3713  if (rowrhs!=rowrhsUse)
3714    delete [] rowrhsUse;
3715  if (rowrng!=rowrngUse)
3716    delete [] rowrngUse;
3717  loadProblem(numcols, numrows, start,  index, value, collb, colub, obj,
3718              rowlb,  rowub);
3719  delete[] rowlb;
3720  delete[] rowub;
3721}
3722// This loads a model from a coinModel object - returns number of errors
3723int 
3724OsiClpSolverInterface::loadFromCoinModel (  CoinModel & modelObject, bool keepSolution)
3725{
3726  modelPtr_->whatsChanged_ = 0;
3727  int numberErrors = 0;
3728  // Set arrays for normal use
3729  double * rowLower = modelObject.rowLowerArray();
3730  double * rowUpper = modelObject.rowUpperArray();
3731  double * columnLower = modelObject.columnLowerArray();
3732  double * columnUpper = modelObject.columnUpperArray();
3733  double * objective = modelObject.objectiveArray();
3734  int * integerType = modelObject.integerTypeArray();
3735  double * associated = modelObject.associatedArray();
3736  // If strings then do copies
3737  if (modelObject.stringsExist()) {
3738    numberErrors = modelObject.createArrays(rowLower, rowUpper, columnLower, columnUpper,
3739                                            objective, integerType,associated);
3740  }
3741  CoinPackedMatrix matrix;
3742  modelObject.createPackedMatrix(matrix,associated);
3743  int numberRows = modelObject.numberRows();
3744  int numberColumns = modelObject.numberColumns();
3745  CoinWarmStart * ws = getWarmStart();
3746  bool restoreBasis = keepSolution && numberRows&&numberRows==getNumRows()&&
3747    numberColumns==getNumCols();
3748  loadProblem(matrix, 
3749              columnLower, columnUpper, objective, rowLower, rowUpper);
3750  if (restoreBasis)
3751    setWarmStart(ws);
3752  delete ws;
3753  // Do names if wanted
3754  int numberItems;
3755  numberItems = modelObject.rowNames()->numberItems();
3756  if (numberItems) {
3757    const char *const * rowNames=modelObject.rowNames()->names();
3758    modelPtr_->copyRowNames(rowNames,0,numberItems);
3759  }
3760  numberItems = modelObject.columnNames()->numberItems();
3761  if (numberItems) {
3762    const char *const * columnNames=modelObject.columnNames()->names();
3763    modelPtr_->copyColumnNames(columnNames,0,numberItems);
3764  }
3765  // Do integers if wanted
3766  assert(integerType);
3767  for (int iColumn=0;iColumn<numberColumns;iColumn++) {
3768    if (integerType[iColumn])
3769      setInteger(iColumn);
3770  }
3771  if (rowLower!=modelObject.rowLowerArray()||
3772      columnLower!=modelObject.columnLowerArray()) {
3773    delete [] rowLower;
3774    delete [] rowUpper;
3775    delete [] columnLower;
3776    delete [] columnUpper;
3777    delete [] objective;
3778    delete [] integerType;
3779    delete [] associated;
3780    //if (numberErrors)
3781    //  handler_->message(CLP_BAD_STRING_VALUES,messages_)
3782    //    <<numberErrors
3783    //    <<CoinMessageEol;
3784  }
3785  modelPtr_->optimizationDirection_ = modelObject.optimizationDirection(); 
3786  return numberErrors;
3787}
3788
3789//-----------------------------------------------------------------------------
3790// Write mps files
3791//-----------------------------------------------------------------------------
3792
3793void OsiClpSolverInterface::writeMps(const char * filename,
3794                                     const char * extension,
3795                                     double objSense) const
3796{
3797  std::string f(filename);
3798  std::string e(extension);
3799  std::string fullname;
3800  if (e!="") {
3801    fullname = f + "." + e;
3802  } else {
3803    // no extension so no trailing period
3804    fullname = f;
3805  }
3806  // get names
3807  const char * const * const rowNames = modelPtr_->rowNamesAsChar();
3808  const char * const * const columnNames = modelPtr_->columnNamesAsChar();
3809  // Fall back on Osi version - possibly with names
3810  OsiSolverInterface::writeMpsNative(fullname.c_str(), 
3811                                     const_cast<const char **>(rowNames),
3812                                     const_cast<const char **>(columnNames),0,2,objSense,
3813                                     numberSOS_,setInfo_);
3814  if (rowNames) {
3815    modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_+1);
3816    modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
3817  }
3818}
3819
3820int 
3821OsiClpSolverInterface::writeMpsNative(const char *filename, 
3822                  const char ** rowNames, const char ** columnNames,
3823                  int formatType,int numberAcross,double objSense) const 
3824{
3825  return OsiSolverInterface::writeMpsNative(filename, rowNames, columnNames,
3826                               formatType, numberAcross,objSense,
3827                                            numberSOS_,setInfo_);
3828}
3829
3830//#############################################################################
3831// CLP specific public interfaces
3832//#############################################################################
3833
3834ClpSimplex * OsiClpSolverInterface::getModelPtr() const
3835{
3836  int saveAlgorithm = lastAlgorithm_;
3837  //freeCachedResults();
3838  lastAlgorithm_ = saveAlgorithm;
3839  //bool inCbcOrOther = (modelPtr_->specialOptions()&0x03000000)!=0;
3840  return modelPtr_;
3841}
3842
3843//-------------------------------------------------------------------
3844
3845//#############################################################################
3846// Constructors, destructors clone and assignment
3847//#############################################################################
3848//-------------------------------------------------------------------
3849// Default Constructor
3850//-------------------------------------------------------------------
3851OsiClpSolverInterface::OsiClpSolverInterface ()
3852:
3853OsiSolverInterface(),
3854rowsense_(NULL),
3855rhs_(NULL),
3856rowrange_(NULL),
3857ws_(NULL),
3858rowActivity_(NULL),
3859columnActivity_(NULL),
3860numberSOS_(0),
3861setInfo_(NULL),
3862smallModel_(NULL),
3863factorization_(NULL),
3864smallestElementInCut_(1.0e-15),
3865smallestChangeInCut_(1.0e-10),
3866largestAway_(-1.0),
3867spareArrays_(NULL),
3868matrixByRow_(NULL),
3869matrixByRowAtContinuous_(NULL), 
3870integerInformation_(NULL),
3871whichRange_(NULL),
3872fakeMinInSimplex_(false),
3873linearObjective_(NULL),
3874cleanupScaling_(0),
3875specialOptions_(0x80000000),
3876baseModel_(NULL),
3877lastNumberRows_(0),
3878continuousModel_(NULL),
3879fakeObjective_(NULL)
3880{
3881  //printf("in default %x\n",this);
3882  modelPtr_=NULL;
3883  notOwned_=false;
3884  disasterHandler_ = new OsiClpDisasterHandler();
3885  reset();
3886}
3887
3888//-------------------------------------------------------------------
3889// Clone
3890//-------------------------------------------------------------------
3891OsiSolverInterface * OsiClpSolverInterface::clone(bool CopyData) const
3892{
3893  //printf("in clone %x\n",this);
3894  OsiClpSolverInterface * newSolver;
3895  if (CopyData) {
3896    newSolver = new OsiClpSolverInterface(*this);
3897  } else {
3898    newSolver = new OsiClpSolverInterface();
3899  }
3900#if 0
3901  const double * obj = newSolver->getObjCoefficients();
3902  const double * oldObj = getObjCoefficients();
3903  if(newSolver->getNumCols()>3787)
3904    printf("%x - obj %x (from %x) val %g\n",newSolver,obj,oldObj,obj[3787]);
3905#endif
3906  return newSolver;
3907}
3908
3909
3910//-------------------------------------------------------------------
3911// Copy constructor
3912//-------------------------------------------------------------------
3913OsiClpSolverInterface::OsiClpSolverInterface (
3914                  const OsiClpSolverInterface & rhs)
3915: OsiSolverInterface(rhs),
3916rowsense_(NULL),
3917rhs_(NULL),
3918rowrange_(NULL),
3919ws_(NULL),
3920rowActivity_(NULL),
3921columnActivity_(NULL),
3922  stuff_(rhs.stuff_),
3923numberSOS_(rhs.numberSOS_),
3924setInfo_(NULL),
3925smallModel_(NULL),
3926factorization_(NULL),
3927smallestElementInCut_(rhs.smallestElementInCut_),
3928smallestChangeInCut_(rhs.smallestChangeInCut_),
3929largestAway_(-1.0),
3930spareArrays_(NULL),
3931basis_(),
3932itlimOrig_(9999999),
3933lastAlgorithm_(0),
3934notOwned_(false),
3935matrixByRow_(NULL),
3936matrixByRowAtContinuous_(NULL), 
3937integerInformation_(NULL),
3938whichRange_(NULL),
3939fakeMinInSimplex_(rhs.fakeMinInSimplex_)
3940{
3941  //printf("in copy %x - > %x\n",&rhs,this);
3942  if ( rhs.modelPtr_  ) 
3943    modelPtr_ = new ClpSimplex(*rhs.modelPtr_);
3944  else
3945    modelPtr_ = new ClpSimplex();
3946  if ( rhs.baseModel_  ) 
3947    baseModel_ = new ClpSimplex(*rhs.baseModel_);
3948  else
3949    baseModel_ = NULL;
3950  if ( rhs.continuousModel_  ) 
3951    continuousModel_ = new ClpSimplex(*rhs.continuousModel_);
3952  else
3953    continuousModel_ = NULL;
3954  if (rhs.matrixByRowAtContinuous_)
3955    matrixByRowAtContinuous_ = new CoinPackedMatrix(*rhs.matrixByRowAtContinuous_);
3956  if ( rhs.disasterHandler_  ) 
3957    disasterHandler_ = dynamic_cast<OsiClpDisasterHandler *>(rhs.disasterHandler_->clone());
3958  else
3959    disasterHandler_ = NULL;
3960  if (rhs.fakeObjective_)
3961    fakeObjective_ = new ClpLinearObjective(*rhs.fakeObjective_);
3962  else
3963    fakeObjective_ = NULL;
3964  linearObjective_ = modelPtr_->objective();
3965  if ( rhs.ws_ ) 
3966    ws_ = new CoinWarmStartBasis(*rhs.ws_);
3967  basis_ = rhs.basis_;
3968  if (rhs.integerInformation_) {
3969    int numberColumns = modelPtr_->numberColumns();
3970    integerInformation_ = new char[numberColumns];
3971    CoinMemcpyN(rhs.integerInformation_,        numberColumns,integerInformation_);
3972  }
3973  saveData_ = rhs.saveData_;
3974  solveOptions_  = rhs.solveOptions_;
3975  cleanupScaling_ = rhs.cleanupScaling_;
3976  specialOptions_ = rhs.specialOptions_;
3977  lastNumberRows_ = rhs.lastNumberRows_;
3978  rowScale_ = rhs.rowScale_;
3979  columnScale_ = rhs.columnScale_;
3980  fillParamMaps();
3981  messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
3982  if (numberSOS_) {
3983    setInfo_ = new CoinSet[numberSOS_];
3984    for (int i=0;i<numberSOS_;i++)
3985      setInfo_[i]=rhs.setInfo_[i];
3986  }
3987}
3988
3989// Borrow constructor - only delete one copy
3990OsiClpSolverInterface::OsiClpSolverInterface (ClpSimplex * rhs,
3991                                              bool reallyOwn)
3992:
3993OsiSolverInterface(),
3994rowsense_(NULL),
3995rhs_(NULL),
3996rowrange_(NULL),
3997ws_(NULL),
3998rowActivity_(NULL),
3999columnActivity_(NULL),
4000numberSOS_(0),
4001setInfo_(NULL),
4002smallModel_(NULL),
4003factorization_(NULL),
4004smallestElementInCut_(1.0e-15),
4005smallestChangeInCut_(1.0e-10),
4006largestAway_(-1.0),
4007spareArrays_(NULL),
4008basis_(),
4009itlimOrig_(9999999),
4010lastAlgorithm_(0),
4011notOwned_(false),
4012matrixByRow_(NULL),
4013matrixByRowAtContinuous_(NULL), 
4014integerInformation_(NULL),
4015whichRange_(NULL),
4016fakeMinInSimplex_(false),
4017cleanupScaling_(0),
4018specialOptions_(0x80000000),
4019baseModel_(NULL),
4020lastNumberRows_(0),
4021continuousModel_(NULL),
4022fakeObjective_(NULL)
4023{
4024  disasterHandler_ = new OsiClpDisasterHandler();
4025  //printf("in borrow %x - > %x\n",&rhs,this);
4026  modelPtr_ = rhs;
4027  basis_.resize(modelPtr_->numberRows(),modelPtr_->numberColumns());
4028  linearObjective_ = modelPtr_->objective();
4029  if (rhs) {
4030    notOwned_=!reallyOwn;
4031
4032    if (rhs->integerInformation()) {
4033      int numberColumns = modelPtr_->numberColumns();
4034      integerInformation_ = new char[numberColumns];
4035      CoinMemcpyN(rhs->integerInformation(),    numberColumns,integerInformation_);
4036    }
4037  }
4038  fillParamMaps();
4039}
4040   
4041// Releases so won't error
4042void 
4043OsiClpSolverInterface::releaseClp()
4044{
4045  modelPtr_=NULL;
4046  notOwned_=false;
4047}
4048   
4049//-------------------------------------------------------------------
4050// Destructor
4051//-------------------------------------------------------------------
4052OsiClpSolverInterface::~OsiClpSolverInterface ()
4053{
4054  //printf("in destructor %x\n",this);
4055  freeCachedResults();
4056  if (!notOwned_)
4057    delete modelPtr_;
4058  delete baseModel_;
4059  delete continuousModel_;
4060  delete disasterHandler_;
4061  delete fakeObjective_;
4062  delete ws_;
4063  delete [] rowActivity_;
4064  delete [] columnActivity_;
4065  delete [] setInfo_;
4066#ifdef KEEP_SMALL
4067  if (smallModel_) {
4068    delete [] spareArrays_;
4069    spareArrays_ = NULL;
4070    delete smallModel_;
4071    smallModel_=NULL;
4072  }
4073#endif
4074  assert(smallModel_==NULL);
4075  assert(factorization_==NULL);
4076  assert(spareArrays_==NULL);
4077  delete [] integerInformation_;
4078  delete matrixByRowAtContinuous_;
4079  delete matrixByRow_;
4080}
4081
4082//-------------------------------------------------------------------
4083// Assignment operator
4084//-------------------------------------------------------------------
4085OsiClpSolverInterface &
4086OsiClpSolverInterface::operator=(const OsiClpSolverInterface& rhs)
4087{
4088  if (this != &rhs) {   
4089    //printf("in = %x - > %x\n",&rhs,this);
4090    OsiSolverInterface::operator=(rhs);
4091    freeCachedResults();
4092    if (!notOwned_)
4093      delete modelPtr_;
4094    delete ws_;
4095    if ( rhs.modelPtr_  ) 
4096      modelPtr_ = new ClpSimplex(*rhs.modelPtr_);
4097    delete baseModel_;
4098    if ( rhs.baseModel_  ) 
4099      baseModel_ = new ClpSimplex(*rhs.baseModel_);
4100    else
4101      baseModel_ = NULL;
4102    delete continuousModel_;
4103    if ( rhs.continuousModel_  ) 
4104      continuousModel_ = new ClpSimplex(*rhs.continuousModel_);
4105    else
4106      continuousModel_ = NULL;
4107    delete matrixByRowAtContinuous_;
4108    delete matrixByRow_;
4109    matrixByRow_=NULL;
4110    if (rhs.matrixByRowAtContinuous_)
4111      matrixByRowAtContinuous_ = new CoinPackedMatrix(*rhs.matrixByRowAtContinuous_);
4112    else
4113      matrixByRowAtContinuous_=NULL;
4114    delete disasterHandler_;
4115    if ( rhs.disasterHandler_  ) 
4116      disasterHandler_ = dynamic_cast<OsiClpDisasterHandler *>(rhs.disasterHandler_->clone());
4117    else
4118      disasterHandler_ = NULL;
4119    delete fakeObjective_;
4120    if (rhs.fakeObjective_)
4121      fakeObjective_ = new ClpLinearObjective(*rhs.fakeObjective_);
4122    else
4123      fakeObjective_ = NULL;
4124    notOwned_=false;
4125    linearObjective_ = modelPtr_->objective();
4126    saveData_ = rhs.saveData_;
4127    solveOptions_  = rhs.solveOptions_;
4128    cleanupScaling_ = rhs.cleanupScaling_;
4129    specialOptions_ = rhs.specialOptions_;
4130    lastNumberRows_ = rhs.lastNumberRows_;
4131    rowScale_ = rhs.rowScale_;
4132    columnScale_ = rhs.columnScale_;
4133    basis_ = rhs.basis_;
4134    stuff_ = rhs.stuff_;
4135    if (rhs.integerInformation_) {
4136      int numberColumns = modelPtr_->numberColumns();
4137      integerInformation_ = new char[numberColumns];
4138      CoinMemcpyN(rhs.integerInformation_,      numberColumns,integerInformation_);
4139    }
4140    if ( rhs.ws_ ) 
4141      ws_ = new CoinWarmStartBasis(*rhs.ws_);
4142    else
4143      ws_=NULL;
4144    delete [] rowActivity_;
4145    delete [] columnActivity_;
4146    rowActivity_=NULL;
4147    columnActivity_=NULL;
4148    delete [] setInfo_;
4149    numberSOS_ = rhs.numberSOS_;
4150    setInfo_=NULL;
4151    if (numberSOS_) {
4152      setInfo_ = new CoinSet[numberSOS_];
4153      for (int i=0;i<numberSOS_;i++)
4154        setInfo_[i]=rhs.setInfo_[i];
4155    }
4156    assert(smallModel_==NULL);
4157    assert(factorization_==NULL);
4158    smallestElementInCut_ = rhs.smallestElementInCut_;
4159    smallestChangeInCut_ = rhs.smallestChangeInCut_;
4160    largestAway_ = -1.0;
4161    assert(spareArrays_==NULL);
4162    basis_ = rhs.basis_;
4163    fillParamMaps();
4164    messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4165  }
4166  return *this;
4167}
4168
4169//#############################################################################
4170// Applying cuts
4171//#############################################################################
4172
4173void OsiClpSolverInterface::applyRowCut( const OsiRowCut & rowCut )
4174{
4175  applyRowCuts(1, &rowCut);
4176}
4177/* Apply a collection of row cuts which are all effective.
4178   applyCuts seems to do one at a time which seems inefficient.
4179*/
4180void 
4181OsiClpSolverInterface::applyRowCuts(int numberCuts, const OsiRowCut * cuts)
4182{
4183  if (numberCuts) {
4184    // Say can't gurantee optimal basis etc
4185    lastAlgorithm_=999;
4186
4187    // Thanks to js
4188    const OsiRowCut * * cutsp = new const OsiRowCut * [numberCuts];
4189    for (int i=0;i<numberCuts;i++) 
4190      cutsp[i] = &cuts[i];
4191   
4192    applyRowCuts(numberCuts, cutsp);
4193   
4194    delete [] cutsp;
4195  }
4196}
4197/* Apply a collection of row cuts which are all effective.
4198   applyCuts seems to do one at a time which seems inefficient.
4199*/
4200void 
4201OsiClpSolverInterface::applyRowCuts(int numberCuts, const OsiRowCut ** cuts)
4202{
4203  int i;
4204  if (!numberCuts)
4205    return;
4206  modelPtr_->whatsChanged_ &= (0xffff&~(1|2|4|16|32));
4207  CoinPackedMatrix * saveRowCopy = matrixByRow_;
4208  matrixByRow_=NULL;
4209#if 0 // was #ifndef NDEBUG
4210  int nameDiscipline;
4211  getIntParam(OsiNameDiscipline,nameDiscipline) ;
4212  assert (!nameDiscipline);
4213#endif
4214  freeCachedResults0();
4215  // Say can't gurantee optimal basis etc
4216  lastAlgorithm_=999;
4217  int numberRows = modelPtr_->numberRows();
4218  modelPtr_->resize(numberRows+numberCuts,modelPtr_->numberColumns());
4219  basis_.resize(numberRows+numberCuts,modelPtr_->numberColumns());
4220  // redo as relaxed - use modelPtr_-> addRows with starts etc
4221  int size = 0;
4222  for (i=0;i<numberCuts;i++) 
4223    size += cuts[i]->row().getNumElements();
4224  CoinBigIndex * starts = new CoinBigIndex [numberCuts+1];
4225  int * indices = new int[size];
4226  double * elements = new double[size];
4227  double * lower = modelPtr_->rowLower()+numberRows;
4228  double * upper = modelPtr_->rowUpper()+numberRows;
4229  const double * columnLower = modelPtr_->columnLower();
4230  const double * columnUpper = modelPtr_->columnUpper();
4231  size=0;
4232  for (i=0;i<numberCuts;i++) {
4233    double rowLb = cuts[i]->lb();
4234    double rowUb = cuts[i]->ub();
4235    int n=cuts[i]->row().getNumElements();
4236    const int * index = cuts[i]->row().getIndices();
4237    const double * elem = cuts[i]->row().getElements();
4238    starts[i]=size;
4239    for (int j=0;j<n;j++) {
4240      double value = elem[j];
4241      int column = index[j];
4242      if (fabs(value)>=smallestChangeInCut_) {
4243        // always take
4244        indices[size]=column;
4245        elements[size++]=value;
4246      } else if (fabs(value)>=smallestElementInCut_) {
4247        double lowerValue = columnLower[column];
4248        double upperValue = columnUpper[column];
4249        double difference = upperValue-lowerValue;
4250        if (difference<1.0e20&&difference*fabs(value)<smallestChangeInCut_&&
4251            (rowLb<-1.0e20||rowUb>1.0e20)) {
4252          // Take out and adjust to relax
4253          //printf("small el %g adjusted\n",value);
4254          if (rowLb>-1.0e20) {
4255            // just lower bound on row
4256            if (value>0.0) {
4257              // pretend at upper
4258              rowLb -= value*upperValue;
4259            } else {
4260              // pretend at lower
4261              rowLb -= value*lowerValue;
4262            }
4263          } else {
4264            // just upper bound on row
4265            if (value>0.0) {
4266              // pretend at lower
4267              rowUb -= value*lowerValue;
4268            } else {
4269              // pretend at upper
4270              rowUb -= value*upperValue;
4271            }
4272          }
4273        } else {
4274          // take (unwillingly)
4275          indices[size]=column;
4276          elements[size++]=value;
4277        }
4278      } else {
4279        //printf("small el %g ignored\n",value);
4280      }
4281    }
4282    lower[i]= forceIntoRange(rowLb, -OsiClpInfinity, OsiClpInfinity);
4283    upper[i]= forceIntoRange(rowUb, -OsiClpInfinity, OsiClpInfinity);
4284    if (lower[i]<-1.0e27)
4285      lower[i]=-COIN_DBL_MAX;
4286    if (upper[i]>1.0e27)
4287      upper[i]=COIN_DBL_MAX;
4288  }
4289  starts[numberCuts]=size;
4290 if (!modelPtr_->clpMatrix())
4291    modelPtr_->createEmptyMatrix();
4292  //modelPtr_->matrix()->appendRows(numberCuts,rows);
4293  modelPtr_->clpMatrix()->appendMatrix(numberCuts,0,starts,indices,elements);
4294  modelPtr_->setNewRowCopy(NULL);
4295  modelPtr_->setClpScaledMatrix(NULL);
4296  freeCachedResults1();
4297  redoScaleFactors( numberCuts,starts, indices, elements);
4298  if (saveRowCopy) {
4299#if 1
4300    matrixByRow_=saveRowCopy;
4301    matrixByRow_->appendRows(numberCuts,starts,indices,elements,0);
4302    if (matrixByRow_->getNumElements()!=modelPtr_->clpMatrix()->getNumElements()) {
4303      delete matrixByRow_; // odd type matrix
4304      matrixByRow_=NULL;
4305    }
4306#else
4307    delete saveRowCopy;
4308#endif
4309  }
4310  delete [] starts;
4311  delete [] indices;
4312  delete [] elements;
4313
4314}
4315//#############################################################################
4316// Apply Cuts
4317//#############################################################################
4318
4319OsiSolverInterface::ApplyCutsReturnCode
4320OsiClpSolverInterface::applyCuts( const OsiCuts & cs, double effectivenessLb ) 
4321{
4322  OsiSolverInterface::ApplyCutsReturnCode retVal;
4323  int i;
4324
4325  // Loop once for each column cut
4326  for ( i=0; i<cs.sizeColCuts(); i ++ ) {
4327    if ( cs.colCut(i).effectiveness() < effectivenessLb ) {
4328      retVal.incrementIneffective();
4329      continue;
4330    }
4331    if ( !cs.colCut(i).consistent() ) {
4332      retVal.incrementInternallyInconsistent();
4333      continue;
4334    }
4335    if ( !cs.colCut(i).consistent(*this) ) {
4336      retVal.incrementExternallyInconsistent();
4337      continue;
4338    }
4339    if ( cs.colCut(i).infeasible(*this) ) {
4340      retVal.incrementInfeasible();
4341      continue;
4342    }
4343    applyColCut( cs.colCut(i) );
4344    retVal.incrementApplied();
4345  }
4346
4347  // Loop once for each row cut
4348  const OsiRowCut ** addCuts = new const OsiRowCut* [cs.sizeRowCuts()];
4349  int nAdd=0;
4350  for ( i=0; i<cs.sizeRowCuts(); i ++ ) {
4351    if ( cs.rowCut(i).effectiveness() < effectivenessLb ) {
4352      retVal.incrementIneffective();
4353      continue;
4354    }
4355    if ( !cs.rowCut(i).consistent() ) {
4356      retVal.incrementInternallyInconsistent();
4357      continue;
4358    }
4359    if ( !cs.rowCut(i).consistent(*this) ) {
4360      retVal.incrementExternallyInconsistent();
4361      continue;
4362    }
4363    if ( cs.rowCut(i).infeasible(*this) ) {
4364      retVal.incrementInfeasible();
4365      continue;
4366    }
4367    addCuts[nAdd++] = cs.rowCutPtr(i);
4368    retVal.incrementApplied();
4369  }
4370  // now apply
4371  applyRowCuts(nAdd,addCuts);
4372  delete [] addCuts;
4373 
4374  return retVal;
4375}
4376// Extend scale factors
4377void 
4378OsiClpSolverInterface::redoScaleFactors(int numberAdd,const CoinBigIndex * starts,
4379                                        const int * indices, const double * elements)
4380{
4381  if ((specialOptions_&131072)!=0) {
4382    int numberRows = modelPtr_->numberRows()-numberAdd;
4383    assert (lastNumberRows_==numberRows); // ???
4384    int iRow;
4385    int newNumberRows = numberRows + numberAdd;
4386    rowScale_.extend(static_cast<int>(2*newNumberRows*sizeof(double)));
4387    double * rowScale = rowScale_.array();
4388    double * oldInverseScale = rowScale + lastNumberRows_;
4389    double * inverseRowScale = rowScale + newNumberRows;
4390    for (iRow=lastNumberRows_-1;iRow>=0;iRow--)
4391      inverseRowScale[iRow] = oldInverseScale[iRow] ;
4392    //int numberColumns = baseModel_->numberColumns();
4393    const double * columnScale = columnScale_.array();
4394    //const double * inverseColumnScale = columnScale + numberColumns;
4395    // Geometric mean on row scales
4396    // adjust arrays
4397    rowScale += lastNumberRows_;
4398    inverseRowScale += lastNumberRows_;
4399    for (iRow=0;iRow<numberAdd;iRow++) {
4400      CoinBigIndex j;
4401      double largest=1.0e-20;
4402      double smallest=1.0e50;
4403      for (j=starts[iRow];j<starts[iRow+1];j++) {
4404        int iColumn = indices[j];
4405        double value = fabs(elements[j]);
4406        // Don't bother with tiny elements
4407        if (value>1.0e-20) {
4408          value *= columnScale[iColumn];
4409          largest = CoinMax(largest,value);
4410          smallest = CoinMin(smallest,value);
4411        }
4412      }
4413      double scale=sqrt(smallest*largest);
4414      scale=CoinMax(1.0e-10,CoinMin(1.0e10,scale));
4415      inverseRowScale[iRow]=scale;
4416      rowScale[iRow]=1.0/scale;
4417    }
4418    lastNumberRows_=newNumberRows;
4419  }
4420}
4421// Delete all scale factor stuff and reset option
4422void OsiClpSolverInterface::deleteScaleFactors()
4423{
4424  delete baseModel_;
4425  baseModel_=NULL;
4426  lastNumberRows_=0;
4427  specialOptions_ &= ~131072;
4428}
4429//-----------------------------------------------------------------------------
4430
4431void OsiClpSolverInterface::applyColCut( const OsiColCut & cc )
4432{
4433  modelPtr_->whatsChanged_ &= (0x1ffff&~(128|256));
4434  // Say can't gurantee optimal basis etc
4435  lastAlgorithm_=999;
4436  double * lower = modelPtr_->columnLower();
4437  double * upper = modelPtr_->columnUpper();
4438  const CoinPackedVector & lbs = cc.lbs();
4439  const CoinPackedVector & ubs = cc.ubs();
4440  int i;
4441
4442  for ( i=0; i<lbs.getNumElements(); i++ ) {
4443    int iCol = lbs.getIndices()[i];
4444    double value = lbs.getElements()[i];
4445    if ( value > lower[iCol] )
4446      lower[iCol]= value;
4447  }
4448  for ( i=0; i<ubs.getNumElements(); i++ ) {
4449    int iCol = ubs.getIndices()[i];
4450    double value = ubs.getElements()[i];
4451    if ( value < upper[iCol] )
4452      upper[iCol]= value;
4453  }
4454}
4455//#############################################################################
4456// Private methods
4457//#############################################################################
4458
4459
4460//-------------------------------------------------------------------
4461
4462void OsiClpSolverInterface::freeCachedResults() const
4463{ 
4464  // Say can't gurantee optimal basis etc
4465  lastAlgorithm_=999;
4466  delete [] rowsense_;
4467  delete [] rhs_;
4468  delete [] rowrange_;
4469  delete matrixByRow_;
4470  if (modelPtr_&&modelPtr_->scaledMatrix_) {
4471    delete modelPtr_->scaledMatrix_;
4472    modelPtr_->scaledMatrix_=NULL;
4473  }
4474  //delete ws_;
4475  rowsense_=NULL;
4476  rhs_=NULL;
4477  rowrange_=NULL;
4478  matrixByRow_=NULL;
4479  //ws_ = NULL;
4480  if (modelPtr_&&modelPtr_->clpMatrix()) {
4481    modelPtr_->clpMatrix()->refresh(modelPtr_); // make sure all clean
4482#ifndef NDEBUG
4483    ClpPackedMatrix * clpMatrix = dynamic_cast<ClpPackedMatrix *> (modelPtr_->clpMatrix());
4484    if (clpMatrix) {
4485      if (clpMatrix->getNumCols())
4486        assert (clpMatrix->getNumRows()==modelPtr_->getNumRows());
4487      if (clpMatrix->getNumRows())
4488        assert (clpMatrix->getNumCols()==modelPtr_->getNumCols());
4489    }
4490#endif
4491  }
4492}
4493
4494//-------------------------------------------------------------------
4495
4496void OsiClpSolverInterface::freeCachedResults0() const
4497{ 
4498  delete [] rowsense_;
4499  delete [] rhs_;
4500  delete [] rowrange_;
4501  rowsense_=NULL;
4502  rhs_=NULL;
4503  rowrange_=NULL;
4504}
4505
4506//-------------------------------------------------------------------
4507
4508void OsiClpSolverInterface::freeCachedResults1() const
4509{ 
4510  // Say can't gurantee optimal basis etc
4511  lastAlgorithm_=999;
4512  delete matrixByRow_;
4513  matrixByRow_=NULL;
4514  //ws_ = NULL;
4515  if (modelPtr_&&modelPtr_->clpMatrix()) {
4516    delete modelPtr_->scaledMatrix_;
4517    modelPtr_->scaledMatrix_=NULL;
4518    modelPtr_->clpMatrix()->refresh(modelPtr_); // make sure all clean
4519#ifndef NDEBUG
4520    ClpPackedMatrix * clpMatrix = dynamic_cast<ClpPackedMatrix *> (modelPtr_->clpMatrix());
4521    if (clpMatrix) {
4522      assert (clpMatrix->getNumRows()==modelPtr_->getNumRows());
4523      assert (clpMatrix->getNumCols()==modelPtr_->getNumCols());
4524    }
4525#endif
4526  }
4527}
4528
4529//------------------------------------------------------------------
4530void OsiClpSolverInterface::extractSenseRhsRange() const
4531{
4532  if (rowsense_ == NULL) {
4533    // all three must be NULL
4534    assert ((rhs_ == NULL) && (rowrange_ == NULL));
4535   
4536    int nr=modelPtr_->numberRows();
4537    if ( nr!=0 ) {
4538      rowsense_ = new char[nr];
4539      rhs_ = new double[nr];
4540      rowrange_ = new double[nr];
4541      std::fill(rowrange_,rowrange_+nr,0.0);
4542     
4543      const double * lb = modelPtr_->rowLower();
4544      const double * ub = modelPtr_->rowUpper();
4545     
4546      int i;
4547      for ( i=0; i<nr; i++ ) {
4548        convertBoundToSense(lb[i], ub[i], rowsense_[i], rhs_[i], rowrange_[i]);
4549      }
4550    }
4551  }
4552}
4553// Set language
4554void 
4555OsiClpSolverInterface::newLanguage(CoinMessages::Language language)
4556{
4557  modelPtr_->newLanguage(language);
4558  OsiSolverInterface::newLanguage(language);
4559}
4560//#############################################################################
4561
4562void
4563OsiClpSolverInterface::fillParamMaps()
4564{
4565  assert (static_cast<int> (OsiMaxNumIteration)==        static_cast<int>(ClpMaxNumIteration));
4566  assert (static_cast<int> (OsiMaxNumIterationHotStart)==static_cast<int>(ClpMaxNumIterationHotStart));
4567  //assert (static_cast<int> (OsiLastIntParam)==           static_cast<int>(ClpLastIntParam));
4568 
4569  assert (static_cast<int> (OsiDualObjectiveLimit)==  static_cast<int>(ClpDualObjectiveLimit));
4570  assert (static_cast<int> (OsiPrimalObjectiveLimit)==static_cast<int>(ClpPrimalObjectiveLimit));
4571  assert (static_cast<int> (OsiDualTolerance)==       static_cast<int>(ClpDualTolerance));
4572  assert (static_cast<int> (OsiPrimalTolerance)==     static_cast<int>(ClpPrimalTolerance));
4573  assert (static_cast<int> (OsiObjOffset)==           static_cast<int>(ClpObjOffset));
4574  //assert (static_cast<int> (OsiLastDblParam)==        static_cast<int>(ClpLastDblParam));
4575 
4576  assert (static_cast<int> (OsiProbName)==    static_cast<int> (ClpProbName));
4577  //strParamMap_[OsiLastStrParam] = ClpLastStrParam;
4578}
4579// Sets up basis
4580void 
4581OsiClpSolverInterface::setBasis ( const CoinWarmStartBasis & basis)
4582{
4583  setBasis(basis,modelPtr_);
4584  setWarmStart(&basis); 
4585}
4586// Warm start
4587CoinWarmStartBasis
4588OsiClpSolverInterface::getBasis(ClpSimplex * model) const
4589{
4590  int iRow,iColumn;
4591  int numberRows = model->numberRows();
4592  int numberColumns = model->numberColumns();
4593  CoinWarmStartBasis basis;
4594  basis.setSize(numberColumns,numberRows);
4595  if (model->statusExists()) {
4596    // Flip slacks
4597    int lookupA[]={0,1,3,2,0,2};
4598    for (iRow=0;iRow<numberRows;iRow++) {
4599      int iStatus = model->getRowStatus(iRow);
4600      iStatus = lookupA[iStatus];
4601      basis.setArtifStatus(iRow,static_cast<CoinWarmStartBasis::Status> (iStatus));
4602    }
4603    int lookupS[]={0,1,2,3,0,3};
4604    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4605      int iStatus = model->getColumnStatus(iColumn);
4606      iStatus = lookupS[iStatus];
4607      basis.setStructStatus(iColumn,static_cast<CoinWarmStartBasis::Status> (iStatus));
4608    }
4609  }
4610  //basis.print();
4611  return basis;
4612}
4613// Warm start from statusArray
4614CoinWarmStartBasis * 
4615OsiClpSolverInterface::getBasis(const unsigned char * statusArray) const 
4616{
4617  int iRow,iColumn;
4618  int numberRows = modelPtr_->numberRows();
4619  int numberColumns = modelPtr_->numberColumns();
4620  CoinWarmStartBasis * basis = new CoinWarmStartBasis();
4621  basis->setSize(numberColumns,numberRows);
4622  // Flip slacks
4623  int lookupA[]={0,1,3,2,0,2};
4624  for (iRow=0;iRow<numberRows;iRow++) {
4625    int iStatus = statusArray[numberColumns+iRow]&7;
4626    iStatus = lookupA[iStatus];
4627    basis->setArtifStatus(iRow,static_cast<CoinWarmStartBasis::Status> (iStatus));
4628  }
4629  int lookupS[]={0,1,2,3,0,3};
4630  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4631    int iStatus = statusArray[iColumn]&7;
4632    iStatus = lookupS[iStatus];
4633    basis->setStructStatus(iColumn,static_cast<CoinWarmStartBasis::Status> (iStatus));
4634  }
4635  //basis->print();
4636  return basis;
4637}
4638// Sets up basis
4639void 
4640OsiClpSolverInterface::setBasis ( const CoinWarmStartBasis & basis,
4641                                  ClpSimplex * model)
4642{
4643  // Say can't gurantee optimal basis etc
4644  lastAlgorithm_=999;
4645  // transform basis to status arrays
4646  int iRow,iColumn;
4647  int numberRows = model->numberRows();
4648  int numberColumns = model->numberColumns();
4649  if (!model->statusExists()) {
4650    /*
4651      get status arrays
4652      ClpBasis would seem to have overheads and we will need
4653      extra bits anyway.
4654    */
4655    model->createStatus();
4656  }
4657  if (basis.getNumArtificial()!=numberRows||
4658      basis.getNumStructural()!=numberColumns) {
4659    CoinWarmStartBasis basis2 = basis;
4660    // resize
4661    basis2.resize(numberRows,numberColumns);
4662    // move status
4663    model->createStatus();
4664    // For rows lower and upper are flipped
4665    for (iRow=0;iRow<numberRows;iRow++) {
4666      int stat = basis2.getArtifStatus(iRow);
4667      if (stat>1)
4668        stat = 5 - stat; // so 2->3 and 3->2
4669      model->setRowStatus(iRow, static_cast<ClpSimplex::Status> (stat));
4670    }
4671    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4672      model->setColumnStatus(iColumn,
4673                             static_cast<ClpSimplex::Status> (basis2.getStructStatus(iColumn)));
4674    }
4675  } else {
4676    // move status
4677    model->createStatus();
4678    // For rows lower and upper are flipped
4679    for (iRow=0;iRow<numberRows;iRow++) {
4680      int stat = basis.getArtifStatus(iRow);
4681      if (stat>1)
4682        stat = 5 - stat; // so 2->3 and 3->2
4683      model->setRowStatus(iRow, static_cast<ClpSimplex::Status> (stat));
4684    }
4685    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4686      model->setColumnStatus(iColumn,
4687                             static_cast<ClpSimplex::Status> (basis.getStructStatus(iColumn)));
4688    }
4689  }
4690}
4691// Warm start difference from basis_ to statusArray
4692CoinWarmStartDiff * 
4693OsiClpSolverInterface::getBasisDiff(const unsigned char * statusArray) const
4694{
4695  int iRow,iColumn;
4696  int numberRows = modelPtr_->numberRows();
4697  int numberColumns = modelPtr_->numberColumns();
4698  CoinWarmStartBasis basis;
4699  basis.setSize(numberColumns,numberRows);
4700  assert (modelPtr_->statusExists());
4701  int lookupS[]={0,1,2,3,0,3};
4702  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4703    int iStatus = statusArray[iColumn]&7;
4704    iStatus = lookupS[iStatus];
4705    basis.setStructStatus(iColumn,static_cast<CoinWarmStartBasis::Status> (iStatus));
4706  }
4707  statusArray += numberColumns;
4708  // Flip slacks
4709  int lookupA[]={0,1,3,2,0,2};
4710  for (iRow=0;iRow<numberRows;iRow++) {
4711    int iStatus = statusArray[iRow]&7;
4712    iStatus = lookupA[iStatus];
4713    basis.setArtifStatus(iRow,static_cast<CoinWarmStartBasis::Status> (iStatus));
4714  }
4715  // Now basis is what we want while basis_ is old
4716  CoinWarmStartDiff * difference = basis.generateDiff(&basis_);
4717  return difference;
4718}
4719/*
4720  Read an mps file from the given filename - returns number of errors
4721  (see CoinMpsIO class)
4722*/
4723int 
4724OsiClpSolverInterface::readMps(const char *filename,
4725                               const char *extension ) 
4726{
4727  // Get rid of integer stuff
4728  delete [] integerInformation_;
4729  integerInformation_=NULL;
4730  freeCachedResults();
4731 
4732  CoinMpsIO m;
4733  m.setInfinity(getInfinity());
4734  m.passInMessageHandler(modelPtr_->messageHandler());
4735  *m.messagesPointer()=modelPtr_->coinMessages();
4736
4737  delete [] setInfo_;
4738  setInfo_=NULL;
4739  numberSOS_=0;
4740  CoinSet ** sets=NULL;
4741  // Temporarily reduce log level to get CoinMpsIO to shut up.
4742  int saveLogLevel = modelPtr_->messageHandler()->logLevel() ;
4743  modelPtr_->messageHandler()->setLogLevel(0) ;
4744  int numberErrors = m.readMps(filename,extension,numberSOS_,sets);
4745  modelPtr_->messageHandler()->setLogLevel(saveLogLevel) ;
4746  if (numberSOS_) {
4747    setInfo_ = new CoinSet[numberSOS_];
4748    for (int i=0;i<numberSOS_;i++) {
4749      setInfo_[i]=*sets[i];
4750      delete sets[i];
4751    }
4752    delete [] sets;
4753  }
4754  handler_->message(COIN_SOLVER_MPS,messages_)
4755    <<m.getProblemName()<< numberErrors <<CoinMessageEol;
4756  if (!numberErrors) {
4757
4758    // set objective function offest
4759    setDblParam(OsiObjOffset,m.objectiveOffset());
4760
4761    // set problem name
4762    setStrParam(OsiProbName,m.getProblemName());
4763   
4764    // no errors
4765    loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
4766                m.getObjCoefficients(),m.getRowSense(),m.getRightHandSide(),
4767                m.getRowRange());
4768    const char * integer = m.integerColumns();
4769    int nCols=m.getNumCols();
4770    int nRows=m.getNumRows();
4771    if (integer) {
4772      int i,n=0;
4773      int * index = new int [nCols];
4774      for (i=0;i<nCols;i++) {
4775        if (integer[i]) {
4776          index[n++]=i;
4777        }
4778      }
4779      setInteger(index,n);
4780      delete [] index;
4781      if (n) 
4782        modelPtr_->copyInIntegerInformation(integer);
4783    }
4784
4785    // set objective name
4786    setObjName(m.getObjectiveName());
4787
4788    // Always keep names
4789    int nameDiscipline;
4790    getIntParam(OsiNameDiscipline,nameDiscipline) ;
4791    int iRow;
4792    std::vector<std::string> rowNames = std::vector<std::string> ();
4793    std::vector<std::string> columnNames = std::vector<std::string> ();
4794    rowNames.reserve(nRows);
4795    for (iRow=0;iRow<nRows;iRow++) {
4796      const char * name = m.rowName(iRow);
4797      rowNames.push_back(name);
4798      if (nameDiscipline) 
4799        OsiSolverInterface::setRowName(iRow,name) ;
4800    }
4801   
4802    int iColumn;
4803    columnNames.reserve(nCols);
4804    for (iColumn=0;iColumn<nCols;iColumn++) {
4805      const char * name = m.columnName(iColumn);
4806      columnNames.push_back(name);
4807      if (nameDiscipline) 
4808        OsiSolverInterface::setColName(iColumn,name) ;
4809    }
4810    modelPtr_->copyNames(rowNames,columnNames);
4811  }
4812  return numberErrors;
4813}
4814int 
4815OsiClpSolverInterface::readMps(const char *filename, const char*extension,
4816                            int & numberSets, CoinSet ** & sets)
4817{
4818  int numberErrors = readMps(filename,extension);
4819  numberSets= numberSOS_;
4820  sets = &setInfo_;
4821  return numberErrors;
4822}
4823/* Read an mps file from the given filename returns
4824   number of errors (see OsiMpsReader class) */
4825int 
4826OsiClpSolverInterface::readMps(const char *filename,bool keepNames,bool allowErrors)
4827{
4828  // Get rid of integer stuff
4829  delete [] integerInformation_;
4830  integerInformation_=NULL;
4831  freeCachedResults();
4832 
4833  CoinMpsIO m;
4834  m.setInfinity(getInfinity());
4835  m.passInMessageHandler(modelPtr_->messageHandler());
4836  *m.messagesPointer()=modelPtr_->coinMessages();
4837  m.setSmallElementValue(CoinMax(modelPtr_->getSmallElementValue(), 
4838                                 m.getSmallElementValue()));
4839
4840  delete [] setInfo_;
4841  setInfo_=NULL;
4842  numberSOS_=0;
4843  CoinSet ** sets=NULL;
4844  int numberErrors = m.readMps(filename,"",numberSOS_,sets);
4845  if (numberSOS_) {
4846    setInfo_ = new CoinSet[numberSOS_];
4847    for (int i=0;i<numberSOS_;i++) {
4848      setInfo_[i]=*sets[i];
4849      delete sets[i];
4850    }
4851    delete [] sets;
4852  }
4853  handler_->message(COIN_SOLVER_MPS,messages_)
4854    <<m.getProblemName()<< numberErrors <<CoinMessageEol;
4855  if (!numberErrors||((numberErrors>0&&numberErrors<100000)&&allowErrors)) {
4856
4857    // set objective function offest
4858    setDblParam(OsiObjOffset,m.objectiveOffset());
4859
4860    // set problem name
4861    setStrParam(OsiProbName,m.getProblemName());
4862
4863    // set objective name
4864    setObjName(m.getObjectiveName());
4865
4866    // no errors
4867    loadProblem(*m.getMatrixByCol(),m.getColLower(),m.getColUpper(),
4868                m.getObjCoefficients(),m.getRowSense(),m.getRightHandSide(),
4869                m.getRowRange());
4870    int nCols=m.getNumCols();
4871    // get quadratic part
4872    if (m.reader()->whichSection (  ) == COIN_QUAD_SECTION ) {
4873      int * start=NULL;
4874      int * column = NULL;
4875      double * element = NULL;
4876      int status=m.readQuadraticMps(NULL,start,column,element,2);
4877      if (!status) 
4878        modelPtr_->loadQuadraticObjective(nCols,start,column,element);
4879      delete [] start;
4880      delete [] column;
4881      delete [] element;
4882    }
4883    const char * integer = m.integerColumns();
4884    int nRows=m.getNumRows();
4885    if (integer) {
4886      int i,n=0;
4887      int * index = new int [nCols];
4888      for (i=0;i<nCols;i++) {
4889        if (integer[i]) {
4890          index[n++]=i;
4891        }
4892      }
4893      setInteger(index,n);
4894      delete [] index;
4895      if (n) 
4896        modelPtr_->copyInIntegerInformation(integer);
4897    }
4898    if (keepNames) {
4899      // keep names
4900      int nameDiscipline;
4901      getIntParam(OsiNameDiscipline,nameDiscipline) ;
4902      int iRow;
4903      std::vector<std::string> rowNames = std::vector<std::string> ();
4904      std::vector<std::string> columnNames = std::vector<std::string> ();
4905      rowNames.reserve(nRows);
4906      for (iRow=0;iRow<nRows;iRow++) {
4907        const char * name = m.rowName(iRow);
4908        rowNames.push_back(name);
4909        if (nameDiscipline) 
4910          OsiSolverInterface::setRowName(iRow,name) ;
4911      }
4912     
4913      int iColumn;
4914      columnNames.reserve(nCols);
4915      for (iColumn=0;iColumn<nCols;iColumn++) {
4916        const char * name = m.columnName(iColumn);
4917        columnNames.push_back(name);
4918        if (nameDiscipline) 
4919          OsiSolverInterface::setColName(iColumn,name) ;
4920      }
4921      modelPtr_->copyNames(rowNames,columnNames);
4922    }
4923  }
4924  return numberErrors;
4925}
4926// Read file in LP format (with names)
4927int 
4928OsiClpSolverInterface::readLp(const char *filename, const double epsilon )
4929{
4930  CoinLpIO m;
4931  m.passInMessageHandler(modelPtr_->messageHandler());
4932  *m.messagesPointer()=modelPtr_->coinMessages();
4933  m.readLp(filename, epsilon);
4934  freeCachedResults();
4935
4936  // set objective function offest
4937  setDblParam(OsiObjOffset, 0);
4938
4939  // set problem name
4940  setStrParam(OsiProbName, m.getProblemName());
4941
4942  // set objective name
4943  setObjName(m.getObjName());
4944
4945  // no errors
4946  loadProblem(*m.getMatrixByRow(), m.getColLower(), m.getColUpper(),
4947              m.getObjCoefficients(), m.getRowLower(), m.getRowUpper());
4948
4949  const char *integer = m.integerColumns();
4950  int nCols = m.getNumCols();
4951  int nRows = m.getNumRows();
4952  if (integer) {
4953    int i, n = 0;
4954    int *index = new int [nCols];
4955    for (i=0; i<nCols; i++) {
4956      if (integer[i]) {
4957        index[n++] = i;
4958      }
4959    }
4960    setInteger(index,n);
4961    delete [] index;
4962  }
4963  // Always keep names
4964  int nameDiscipline;
4965  getIntParam(OsiNameDiscipline,nameDiscipline) ;
4966  int iRow;
4967  std::vector<std::string> rowNames = std::vector<std::string> ();
4968  std::vector<std::string> columnNames = std::vector<std::string> ();
4969  rowNames.reserve(nRows);
4970  for (iRow=0;iRow<nRows;iRow++) {
4971    const char * name = m.rowName(iRow);
4972    rowNames.push_back(name);
4973    if (nameDiscipline) 
4974      OsiSolverInterface::setRowName(iRow,name) ;
4975  }
4976 
4977  int iColumn;
4978  columnNames.reserve(nCols);
4979  for (iColumn=0;iColumn<nCols;iColumn++) {
4980    const char * name = m.columnName(iColumn);
4981    columnNames.push_back(name);
4982    if (nameDiscipline) 
4983      OsiSolverInterface::setColName(iColumn,name) ;
4984  }
4985  modelPtr_->copyNames(rowNames,columnNames);
4986  return(0);
4987}
4988/* Write the problem into an Lp file of the given filename.
4989   If objSense is non zero then -1.0 forces the code to write a
4990   maximization objective and +1.0 to write a minimization one.
4991   If 0.0 then solver can do what it wants.
4992   This version calls writeLpNative with names */
4993void 
4994OsiClpSolverInterface::writeLp(const char *filename,
4995                               const char *extension ,
4996                               double epsilon ,
4997                               int numberAcross ,
4998                               int decimals ,
4999                               double objSense ,
5000                               bool changeNameOnRange) const
5001{
5002  std::string f(filename);
5003  std::string e(extension);
5004  std::string fullname;
5005  if (e!="") {
5006    fullname = f + "." + e;
5007  } else {
5008    // no extension so no trailing period
5009    fullname = f;
5010  }
5011  // get names
5012  const char * const * const rowNames = modelPtr_->rowNamesAsChar();
5013  const char * const * const columnNames = modelPtr_->columnNamesAsChar();
5014  // Fall back on Osi version - possibly with names
5015  OsiSolverInterface::writeLpNative(fullname.c_str(), 
5016                                    rowNames,columnNames, epsilon, numberAcross,
5017                                    decimals, objSense,changeNameOnRange);
5018  if (rowNames) {
5019    modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_+1);
5020    modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
5021  }
5022}
5023void 
5024OsiClpSolverInterface::writeLp(FILE * fp,
5025                               double epsilon ,
5026                               int numberAcross ,
5027                               int decimals ,
5028                               double objSense ,
5029                               bool changeNameOnRange) const
5030{
5031  // get names
5032  const char * const * const rowNames = modelPtr_->rowNamesAsChar();
5033  const char * const * const columnNames = modelPtr_->columnNamesAsChar();
5034  // Fall back on Osi version - possibly with names
5035  OsiSolverInterface::writeLpNative(fp,
5036                                    rowNames,columnNames, epsilon, numberAcross,
5037                                    decimals, objSense,changeNameOnRange);
5038  if (rowNames) {
5039    modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_+1);
5040    modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
5041  }
5042}
5043/*
5044  I (JJF) am getting incredibly annoyed because I can't just replace a matrix.
5045  The default behavior of this is do nothing so only use where that would not matter
5046  e.g. strengthening a matrix for MIP
5047*/
5048void 
5049OsiClpSolverInterface::replaceMatrixOptional(const CoinPackedMatrix & matrix)
5050{
5051  modelPtr_->whatsChanged_ &= (0xffff&~(2|4|8));
5052  replaceMatrix(matrix);
5053}
5054// And if it does matter (not used at present)
5055void 
5056OsiClpSolverInterface::replaceMatrix(const CoinPackedMatrix & matrix)
5057{
5058  modelPtr_->whatsChanged_ &= (0xffff&~(2|4|8));
5059  delete modelPtr_->matrix_;
5060  delete modelPtr_->rowCopy_;
5061  modelPtr_->rowCopy_=NULL;
5062  if (matrix.isColOrdered()) {
5063    modelPtr_->matrix_=new ClpPackedMatrix(matrix);
5064  } else {
5065    CoinPackedMatrix matrix2;
5066    matrix2.setExtraGap(0.0);
5067    matrix2.setExtraMajor(0.0);
5068    matrix2.reverseOrderedCopyOf(matrix);
5069    modelPtr_->matrix_=new ClpPackedMatrix(matrix2);
5070  }   
5071  modelPtr_->matrix_->setDimensions(modelPtr_->numberRows_,modelPtr_->numberColumns_);
5072  freeCachedResults();
5073}
5074// Get pointer to array[getNumCols()] of primal solution vector
5075const double * 
5076OsiClpSolverInterface::getColSolution() const 
5077{ 
5078  if (modelPtr_->solveType()!=2) {
5079    return modelPtr_->primalColumnSolution();
5080  } else {
5081    // simplex interface
5082    return modelPtr_->solutionRegion(1);
5083  }
5084}
5085 
5086// Get pointer to array[getNumRows()] of dual prices
5087const double * 
5088OsiClpSolverInterface::getRowPrice() const
5089{ 
5090  if (modelPtr_->solveType()!=2) {
5091    return modelPtr_->dualRowSolution();
5092  } else {
5093    // simplex interface
5094    //return modelPtr_->djRegion(0);
5095    return modelPtr_->dualRowSolution();
5096  }
5097}
5098 
5099// Get a pointer to array[getNumCols()] of reduced costs
5100const double * 
5101OsiClpSolverInterface::getReducedCost() const 
5102{ 
5103  if (modelPtr_->solveType()!=2) {
5104    return modelPtr_->dualColumnSolution();
5105  } else {
5106    // simplex interface
5107    return modelPtr_->djRegion(1);
5108  }
5109}
5110
5111/* Get pointer to array[getNumRows()] of row activity levels (constraint
5112   matrix times the solution vector */
5113const double * 
5114OsiClpSolverInterface::getRowActivity() const 
5115{ 
5116  if (modelPtr_->solveType()!=2) {
5117    return modelPtr_->primalRowSolution();
5118  } else {
5119    // simplex interface
5120    return modelPtr_->solutionRegion(0);
5121  }
5122}
5123double 
5124OsiClpSolverInterface::getObjValue() const 
5125{
5126  if (modelPtr_->numberIterations()||modelPtr_->upperIn_!=-COIN_DBL_MAX) {
5127    // This does not pass unitTest when getObjValue is called before solve.
5128    //printf("obj a %g %g\n",modelPtr_->objectiveValue(),
5129    //     OsiSolverInterface::getObjValue());
5130    if (fakeMinInSimplex_)
5131      return -modelPtr_->objectiveValue() ;
5132    else
5133      return modelPtr_->objectiveValue();
5134  } else {
5135    return OsiSolverInterface::getObjValue();
5136  }
5137}
5138
5139/* Set an objective function coefficient */
5140void 
5141OsiClpSolverInterface::setObjCoeff( int elementIndex, double elementValue )
5142{
5143  modelPtr_->whatsChanged_ &= 0xffff;
5144  // Say can't gurantee optimal basis etc
5145  lastAlgorithm_=999;
5146#ifndef NDEBUG
5147  int n = modelPtr_->numberColumns();
5148  if (elementIndex<0||elementIndex>=n) {
5149    indexError(elementIndex,"setObjCoeff");
5150  }
5151#endif
5152  modelPtr_->setObjectiveCoefficient(elementIndex,
5153    ((fakeMinInSimplex_)?-elementValue:elementValue));
5154}
5155
5156/* Set a single column lower bound<br>
5157   Use -DBL_MAX for -infinity. */
5158void 
5159OsiClpSolverInterface::setColLower( int index, double elementValue )
5160{
5161  modelPtr_->whatsChanged_ &= 0x1ffff;
5162#ifndef NDEBUG
5163  int n = modelPtr_->numberColumns();
5164  if (index<0||index>=n) {
5165    indexError(index,"setColLower");
5166  }
5167#endif
5168  double currentValue = modelPtr_->columnActivity_[index];
5169  bool changed=(currentValue<elementValue-modelPtr_->primalTolerance()||
5170                index>=basis_.getNumStructural()||
5171                basis_.getStructStatus(index)==CoinWarmStartBasis::atLowerBound);
5172  // Say can't gurantee optimal basis etc
5173  if (changed)
5174    lastAlgorithm_=999;
5175  if (!modelPtr_->lower_)
5176    modelPtr_->whatsChanged_ &= ~0xffff; // switch off
5177  modelPtr_->setColumnLower(index,elementValue);
5178}
5179     
5180/* Set a single column upper bound<br>
5181   Use DBL_MAX for infinity. */
5182void 
5183OsiClpSolverInterface::setColUpper( int index, double elementValue )
5184{
5185  modelPtr_->whatsChanged_ &= 0x1ffff;
5186#ifndef NDEBUG
5187  int n = modelPtr_->numberColumns();
5188  if (index<0||index>=n) {
5189    indexError(index,"setColUpper");
5190  }
5191#endif
5192  double currentValue = modelPtr_->columnActivity_[index];
5193  bool changed=(currentValue>elementValue+modelPtr_->primalTolerance()||
5194                index>=basis_.getNumStructural()||
5195                basis_.getStructStatus(index)==CoinWarmStartBasis::atUpperBound);
5196  // Say can't gurantee optimal basis etc
5197  if (changed)
5198    lastAlgorithm_=999;
5199  if (!modelPtr_->upper_)
5200    modelPtr_->whatsChanged_ &= ~0xffff; // switch off
5201  modelPtr_->setColumnUpper(index,elementValue);
5202}
5203
5204/* Set a single column lower and upper bound */
5205void 
5206OsiClpSolverInterface::setColBounds( int elementIndex,
5207                                     double lower, double upper )
5208{
5209  modelPtr_->whatsChanged_ &= 0x1ffff;
5210  // Say can't gurantee optimal basis etc
5211  lastAlgorithm_=999;
5212#ifndef NDEBUG
5213  int n = modelPtr_->numberColumns();
5214  if (elementIndex<0||elementIndex>=n) {
5215    indexError(elementIndex,"setColBounds");
5216  }
5217#endif
5218  if (!modelPtr_->lower_)
5219    modelPtr_->whatsChanged_ &= ~0xffff; // switch off
5220  modelPtr_->setColumnBounds(elementIndex,lower,upper);
5221}
5222void OsiClpSolverInterface::setColSetBounds(const int* indexFirst,
5223                                            const int* indexLast,
5224                                            const double* boundList)
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  const int * indexFirst2=indexFirst;
5232  while (indexFirst2 != indexLast) {
5233    const int iColumn=*indexFirst2++;
5234    if (iColumn<0||iColumn>=n) {
5235      indexError(iColumn,"setColSetBounds");
5236    }
5237  }
5238#endif
5239  modelPtr_->setColSetBounds(indexFirst,indexLast,boundList);
5240}
5241//------------------------------------------------------------------
5242/* Set a single row lower bound<br>
5243   Use -DBL_MAX for -infinity. */
5244void 
5245OsiClpSolverInterface::setRowLower( int elementIndex, double elementValue ) {
5246  // Say can't gurantee optimal basis etc
5247  lastAlgorithm_=999;
5248  modelPtr_->whatsChanged_ &= 0xffff;
5249#ifndef NDEBUG
5250  int n = modelPtr_->numberRows();
5251  if (elementIndex<0||elementIndex>=n) {
5252    indexError(elementIndex,"setRowLower");
5253  }
5254#endif
5255  modelPtr_->setRowLower(elementIndex , elementValue);
5256  if (rowsense_!=NULL) {
5257    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5258    convertBoundToSense(modelPtr_->rowLower_[elementIndex], 
5259                        modelPtr_->rowUpper_[elementIndex],
5260                        rowsense_[elementIndex], rhs_[elementIndex], rowrange_[elementIndex]);
5261  }
5262}
5263     
5264/* Set a single row upper bound<br>
5265   Use DBL_MAX for infinity. */
5266void 
5267OsiClpSolverInterface::setRowUpper( int elementIndex, double elementValue ) {
5268  modelPtr_->whatsChanged_ &= 0xffff;
5269  // Say can't guarantee optimal basis etc
5270  lastAlgorithm_=999;
5271#ifndef NDEBUG
5272  int n = modelPtr_->numberRows();
5273  if (elementIndex<0||elementIndex>=n) {
5274    indexError(elementIndex,"setRowUpper");
5275  }
5276#endif
5277  modelPtr_->setRowUpper(elementIndex , elementValue);
5278  if (rowsense_!=NULL) {
5279    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5280    convertBoundToSense(modelPtr_->rowLower_[elementIndex], 
5281                        modelPtr_->rowUpper_[elementIndex],
5282                        rowsense_[elementIndex], rhs_[elementIndex], rowrange_[elementIndex]);
5283  }
5284}
5285   
5286/* Set a single row lower and upper bound */
5287void 
5288OsiClpSolverInterface::setRowBounds( int elementIndex,
5289              double lower, double upper ) {
5290  modelPtr_->whatsChanged_ &= 0xffff;
5291  // Say can't gurantee optimal basis etc
5292  lastAlgorithm_=999;
5293#ifndef NDEBUG
5294  int n = modelPtr_->numberRows();
5295  if (elementIndex<0||elementIndex>=n) {
5296    indexError(elementIndex,"setRowBounds");
5297  }
5298#endif
5299  modelPtr_->setRowBounds(elementIndex,lower,upper);
5300  if (rowsense_!=NULL) {
5301    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5302    convertBoundToSense(modelPtr_->rowLower_[elementIndex], 
5303                        modelPtr_->rowUpper_[elementIndex],
5304                        rowsense_[elementIndex], rhs_[elementIndex], rowrange_[elementIndex]);
5305  }
5306}
5307//-----------------------------------------------------------------------------
5308void
5309OsiClpSolverInterface::setRowType(int i, char sense, double rightHandSide,
5310                                  double range)
5311{
5312  modelPtr_->whatsChanged_ &= 0xffff;
5313  // Say can't gurantee optimal basis etc
5314  lastAlgorithm_=999;
5315#ifndef NDEBUG
5316  int n = modelPtr_->numberRows();
5317  if (i<0||i>=n) {
5318    indexError(i,"setRowType");
5319  }
5320#endif
5321  double lower = 0, upper = 0;
5322  convertSenseToBound(sense, rightHandSide, range, lower, upper);
5323  setRowBounds(i, lower, upper);
5324  // If user is using sense then set
5325  if (rowsense_) {
5326    rowsense_[i] = sense;
5327    rhs_[i] = rightHandSide;
5328    rowrange_[i] = range;
5329  }
5330}
5331// Set name of row
5332void 
5333//OsiClpSolverInterface::setRowName(int rowIndex, std::string & name)
5334OsiClpSolverInterface::setRowName(int rowIndex, std::string name) 
5335{
5336  if (rowIndex>=0&&rowIndex<modelPtr_->numberRows()) {
5337    int nameDiscipline;
5338    getIntParam(OsiNameDiscipline,nameDiscipline) ;
5339    if (nameDiscipline) {
5340      modelPtr_->setRowName(rowIndex,name);
5341      OsiSolverInterface::setRowName(rowIndex,name) ;
5342    }
5343  }
5344}
5345// Return name of row if one exists or Rnnnnnnn
5346// we ignore maxLen
5347std::string
5348OsiClpSolverInterface::getRowName(int rowIndex, unsigned int /*maxLen*/) const
5349{ 
5350  if (rowIndex == getNumRows())
5351    return getObjName();
5352  int useNames;
5353  getIntParam (OsiNameDiscipline,useNames);
5354  if (useNames)
5355    return modelPtr_->getRowName(rowIndex);
5356  else
5357    return dfltRowColName('r',rowIndex);
5358}
5359   
5360// Set name of col
5361void 
5362//OsiClpSolverInterface::setColName(int colIndex, std::string & name)
5363OsiClpSolverInterface::setColName(int colIndex, std::string name) 
5364{
5365  if (colIndex>=0&&colIndex<modelPtr_->numberColumns()) {
5366    int nameDiscipline;
5367    getIntParam(OsiNameDiscipline,nameDiscipline) ;
5368    if (nameDiscipline) {
5369      modelPtr_->setColumnName(colIndex,name);
5370      OsiSolverInterface::setColName(colIndex,name) ;
5371    }
5372  }
5373}
5374// Return name of col if one exists or Rnnnnnnn
5375std::string
5376OsiClpSolverInterface::getColName(int colIndex, unsigned int /*maxLen*/) const
5377{
5378  int useNames;
5379  getIntParam (OsiNameDiscipline,useNames);
5380  if (useNames)
5381    return modelPtr_->getColumnName(colIndex);
5382  else
5383    return dfltRowColName('c',colIndex);
5384}
5385   
5386   
5387//-----------------------------------------------------------------------------
5388void OsiClpSolverInterface::setRowSetBounds(const int* indexFirst,
5389                                            const int* indexLast,
5390                                            const double* boundList)
5391{
5392  modelPtr_->whatsChanged_ &= 0xffff;
5393  // Say can't gurantee optimal basis etc
5394  lastAlgorithm_=999;
5395#ifndef NDEBUG
5396  int n = modelPtr_->numberRows();
5397  const int * indexFirst2=indexFirst;
5398  while (indexFirst2 != indexLast) {
5399    const int iColumn=*indexFirst2++;
5400    if (iColumn<0||iColumn>=n) {
5401      indexError(iColumn,"setColumnSetBounds");
5402    }
5403  }
5404#endif
5405  modelPtr_->setRowSetBounds(indexFirst,indexLast,boundList);
5406  if (rowsense_ != NULL) {
5407    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5408    double * lower = modelPtr_->rowLower();
5409    double * upper = modelPtr_->rowUpper();
5410    while (indexFirst != indexLast) {
5411      const int iRow=*indexFirst++;
5412      convertBoundToSense(lower[iRow], upper[iRow],
5413                          rowsense_[iRow], rhs_[iRow], rowrange_[iRow]);
5414    }
5415  }
5416}
5417//-----------------------------------------------------------------------------
5418void
5419OsiClpSolverInterface::setRowSetTypes(const int* indexFirst,
5420                                      const int* indexLast,
5421                                      const char* senseList,
5422                                      const double* rhsList,
5423                                      const double* rangeList)
5424{
5425  modelPtr_->whatsChanged_ &= 0xffff;
5426  // Say can't gurantee optimal basis etc
5427  lastAlgorithm_=999;
5428#ifndef NDEBUG
5429  int n = modelPtr_->numberRows();
5430#endif
5431  const int len = static_cast<int>(indexLast - indexFirst);
5432  while (indexFirst != indexLast) {
5433    const int iRow= *indexFirst++;
5434#ifndef NDEBUG
5435    if (iRow<0||iRow>=n) {
5436      indexError(iRow,"isContinuous");
5437    }
5438#endif
5439    double lowerValue = 0;
5440    double upperValue = 0;
5441    if (rangeList){
5442      convertSenseToBound(*senseList++, *rhsList++, *rangeList++,
5443                          lowerValue, upperValue);
5444    } else {
5445      convertSenseToBound(*senseList++, *rhsList++, 0,
5446                          lowerValue, upperValue);
5447    }
5448    modelPtr_->setRowBounds(iRow,lowerValue,upperValue);
5449  }
5450  if (rowsense_ != NULL) {
5451    assert ((rhs_ != NULL) && (rowrange_ != NULL));
5452    indexFirst -= len;
5453    senseList -= len;
5454    rhsList -= len;
5455    if (rangeList)
5456       rangeList -= len;
5457    while (indexFirst != indexLast) {
5458      const int iRow=*indexFirst++;
5459      rowsense_[iRow] = *senseList++;
5460      rhs_[iRow] = *rhsList++;
5461      if (rangeList)
5462         rowrange_[iRow] = *rangeList++;
5463    }
5464  }
5465}
5466
5467/*
5468  Clp's copy-in/copy-out design paradigm is a challenge for the simplex modes.
5469  Normal operation goes like this:
5470    * startup() loads clp's work arrays, performing scaling for numerical
5471      stability and compensating for max.
5472    * clp solves the problem
5473    * finish() unloads the work arrays into answer arrays, undoing scaling
5474      and max compensation.
5475  There are two solutions: undo scaling and max on demand, or make them
5476  into noops. The various getBInv* methods undo scaling on demand (but
5477  see special option 512) and do not need to worry about max. Other get
5478  solution methods are not coded to do this, so the second approach is
5479  used. For simplex modes, turn off scaling (necessary for both primal and
5480  dual solutions) and temporarily convert max to min (necessary for dual
5481  solution). This makes the unscaling in getBInv* superfluous, but don't
5482  remove it. Arguably the better solution here would be to go through and
5483  add unscaling and max compensation to the get solution methods. Look for
5484  fakeMinInSimplex to see the places this propagates to.
5485
5486  TODO: setRowPrice never has worked properly, and I didn't try to fix it in
5487        this go-round.
5488
5489  As of 100907, change applied to [enable|disable]Factorization (mode 1).
5490  Limitation of [enable|disable]SimplexInterface (mode 2) noted in
5491  documentation.  -- lh, 100907 --
5492*/
5493/*
5494  Enables normal operation of subsequent functions.  This method is supposed
5495  to ensure that all typical things (like reduced costs, etc.) are updated
5496  when individual pivots are executed and can be queried by other methods
5497*/
5498void 
5499OsiClpSolverInterface::enableSimplexInterface(bool doingPrimal)
5500{
5501  modelPtr_->whatsChanged_ &= 0xffff;
5502  if (modelPtr_->solveType()==2)
5503    return;
5504  assert (modelPtr_->solveType()==1);
5505  int saveIts = modelPtr_->numberIterations_;
5506  modelPtr_->setSolveType(2);
5507  if (doingPrimal)
5508    modelPtr_->setAlgorithm(1);
5509  else
5510    modelPtr_->setAlgorithm(-1);
5511  // Do initialization
5512  saveData_ = modelPtr_->saveData();
5513  saveData_.scalingFlag_=modelPtr_->scalingFlag();
5514  modelPtr_->scaling(0);
5515  specialOptions_ = 0x80000000;
5516  // set infeasibility cost up
5517  modelPtr_->setInfeasibilityCost(1.0e12);
5518  ClpDualRowDantzig dantzig;
5519  modelPtr_->setDualRowPivotAlgorithm(dantzig);
5520  ClpPrimalColumnDantzig dantzigP;
5521  dantzigP.saveWeights(modelPtr_,0); // set modelPtr
5522  modelPtr_->setPrimalColumnPivotAlgorithm(dantzigP);
5523  int saveOptions = modelPtr_->specialOptions_;
5524  modelPtr_->specialOptions_ &= ~262144;
5525  delete modelPtr_->scaledMatrix_;
5526  modelPtr_->scaledMatrix_=NULL;
5527  // make sure using standard factorization
5528  modelPtr_->factorization()->forceOtherFactorization(4);
5529#ifdef NDEBUG
5530  modelPtr_->startup(0);
5531#else
5532  int returnCode=modelPtr_->startup(0);
5533  assert (!returnCode||returnCode==2);
5534#endif
5535  modelPtr_->specialOptions_=saveOptions;
5536  modelPtr_->numberIterations_=saveIts;
5537}
5538
5539//Undo whatever setting changes the above method had to make
5540void 
5541OsiClpSolverInterface::disableSimplexInterface()
5542{
5543  modelPtr_->whatsChanged_ &= 0xffff;
5544  assert (modelPtr_->solveType()==2);
5545  // declare optimality anyway  (for message handler)
5546  modelPtr_->setProblemStatus(0);
5547  modelPtr_->setSolveType(1);
5548  // message will not appear anyway
5549  int saveMessageLevel=modelPtr_->messageHandler()->logLevel();
5550  modelPtr_->messageHandler()->setLogLevel(0);
5551  modelPtr_->finish();
5552  modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
5553  modelPtr_->restoreData(saveData_);
5554  modelPtr_->scaling(saveData_.scalingFlag_);
5555  ClpDualRowSteepest steepest;
5556  modelPtr_->setDualRowPivotAlgorithm(steepest);
5557  ClpPrimalColumnSteepest steepestP;
5558  modelPtr_->setPrimalColumnPivotAlgorithm(steepestP);
5559  basis_ = getBasis(modelPtr_);
5560  modelPtr_->setSolveType(1);
5561}
5562
5563/*
5564  Force scaling off. If the client thinks we're maximising, arrange it so
5565  that clp sees minimisation while the client still sees maximisation. In
5566  keeping with the spirit of the getBInv methods, special option 512 will
5567  leave all work to the client.
5568*/
5569void 
5570OsiClpSolverInterface::enableFactorization() const
5571{
5572  saveData_.specialOptions_=specialOptions_;
5573  // Try to preserve work regions, reuse factorization
5574  if ((specialOptions_&(1+8))!=1+8)
5575    setSpecialOptionsMutable((1+8)|specialOptions_);
5576  // Are we allowed to make the output sensible to mere mortals?
5577  if ((specialOptions_&512)==0) {
5578    // Force scaling to off
5579    saveData_.scalingFlag_ = modelPtr_->scalingFlag() ;
5580    modelPtr_->scaling(0) ;
5581    // Temporarily force to min but keep a copy of original objective.
5582    if (getObjSense() < 0.0) {
5583      fakeMinInSimplex_ = true ;
5584      modelPtr_->setOptimizationDirection(1.0) ;
5585      double *c = modelPtr_->objective() ;
5586      int n = getNumCols() ;
5587      linearObjective_ = new double[n] ;
5588      CoinMemcpyN(c,n,linearObjective_) ;
5589      std::transform(c,c+n,c,std::negate<double>()) ;
5590    }
5591  }
5592  int saveStatus = modelPtr_->problemStatus_;
5593#ifdef NDEBUG
5594  modelPtr_->startup(0);
5595#else
5596  int returnCode=modelPtr_->startup(0);
5597  assert (!returnCode||returnCode==2);
5598#endif
5599  modelPtr_->problemStatus_=saveStatus;
5600}
5601
5602/*
5603  Undo enableFactorization. Retrieve the special options and scaling and
5604  remove the temporary objective used to fake minimisation in clp.
5605*/
5606void 
5607OsiClpSolverInterface::disableFactorization() const
5608{
5609  specialOptions_=saveData_.specialOptions_;
5610  // declare optimality anyway  (for message handler)
5611  modelPtr_->setProblemStatus(0);
5612  // message will not appear anyway
5613  int saveMessageLevel=modelPtr_->messageHandler()->logLevel();
5614  modelPtr_->messageHandler()->setLogLevel(0);
5615  modelPtr_->finish();
5616  modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
5617  // Client asked for transforms on the way in, so back out.
5618  if ((specialOptions_&512)==0) {
5619    modelPtr_->scaling(saveData_.scalingFlag_) ;
5620    if (fakeMinInSimplex_ == true) {
5621      fakeMinInSimplex_ = false ;
5622      modelPtr_->setOptimizationDirection(-1.0) ;
5623      double *c = modelPtr_->objective() ;
5624      int n = getNumCols() ;
5625      std::transform(c,c+n,c,std::negate<double>()) ;
5626      delete[] linearObjective_ ;
5627    }
5628  }
5629}
5630
5631
5632
5633/* The following two methods may be replaced by the
5634   methods of OsiSolverInterface using OsiWarmStartBasis if:
5635   1. OsiWarmStartBasis resize operation is implemented
5636   more efficiently and
5637   2. It is ensured that effects on the solver are the same
5638   
5639   Returns a basis status of the structural/artificial variables
5640*/
5641void 
5642OsiClpSolverInterface::getBasisStatus(int* cstat, int* rstat) const
5643{
5644  int iRow,iColumn;
5645  int numberRows = modelPtr_->numberRows();
5646  int numberColumns = modelPtr_->numberColumns();
5647  const double * pi = modelPtr_->dualRowSolution();
5648  const double * dj = modelPtr_->dualColumnSolution();
5649  double multiplier = modelPtr_->optimizationDirection();
5650  // Flip slacks
5651  int lookupA[]={0,1,3,2,0,3};
5652  for (iRow=0;iRow<numberRows;iRow++) {
5653    int iStatus = modelPtr_->getRowStatus(iRow);
5654    if (iStatus==5) {
5655      // Fixed - look at reduced cost
5656      if (pi[iRow]*multiplier>1.0e-7)
5657        iStatus = 3;