source: stable/1.14/Clp/src/OsiClp/OsiClpSolverInterface.cpp @ 1886

Last change on this file since 1886 was 1886, checked in by forrest, 7 years ago

add proximity search to possibilities

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