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

Last change on this file since 1729 was 1729, checked in by forrest, 8 years ago

messages for Cbc fathoming and allow perturbation and some safety?

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