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

Last change on this file since 1602 was 1602, checked in by lou, 9 years ago

Make simplex mode 1 (tableau access) work properly with maximisation. Note
that mode 2 (single pivot) is still minimisation only, though it shouldn't be
too hard to extend the fix.

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