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

Last change on this file since 1677 was 1677, checked in by forrest, 9 years ago

gub

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