source: branches/devel/Cbc/src/CoinSolve.cpp @ 452

Last change on this file since 452 was 452, checked in by forrest, 13 years ago

modify CoinSolve? so will solve all 44 with -miplib

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 171.5 KB
Line 
1// Copyright (C) 2004, International Business Machines
2// Corporation and others.  All Rights Reserved.
3   
4#include "CbcConfig.h"
5#include "CoinPragma.hpp"
6
7#include <cassert>
8#include <cstdio>
9#include <cmath>
10#include <cfloat>
11#include <cstring>
12#include <iostream>
13
14
15#include "CoinPragma.hpp"
16#include "CoinHelperFunctions.hpp"
17// Same version as CBC
18#define CBCVERSION "1.02.00"
19
20#include "CoinMpsIO.hpp"
21
22#include "ClpFactorization.hpp"
23#include "CoinTime.hpp"
24#include "ClpSimplex.hpp"
25#include "ClpSimplexOther.hpp"
26#include "ClpSolve.hpp"
27#include "ClpPackedMatrix.hpp"
28#include "ClpPlusMinusOneMatrix.hpp"
29#include "ClpNetworkMatrix.hpp"
30#include "ClpDualRowSteepest.hpp"
31#include "ClpDualRowDantzig.hpp"
32#include "ClpLinearObjective.hpp"
33#include "ClpPrimalColumnSteepest.hpp"
34#include "ClpPrimalColumnDantzig.hpp"
35#include "ClpPresolve.hpp"
36#include "CbcOrClpParam.hpp"
37#include "OsiRowCutDebugger.hpp"
38#include "OsiChooseVariable.hpp"
39#ifdef DMALLOC
40#include "dmalloc.h"
41#endif
42#ifdef WSSMP_BARRIER
43#define FOREIGN_BARRIER
44#endif
45#ifdef UFL_BARRIER
46#define FOREIGN_BARRIER
47#endif
48#ifdef TAUCS_BARRIER
49#define FOREIGN_BARRIER
50#endif
51#include "CoinWarmStartBasis.hpp"
52
53#include "OsiSolverInterface.hpp"
54#include "OsiCuts.hpp"
55#include "OsiRowCut.hpp"
56#include "OsiColCut.hpp"
57
58#include "CglPreProcess.hpp"
59#include "CglCutGenerator.hpp"
60#include "CglGomory.hpp"
61#include "CglProbing.hpp"
62#include "CglKnapsackCover.hpp"
63#include "CglRedSplit.hpp"
64#include "CglClique.hpp"
65#include "CglFlowCover.hpp"
66#include "CglMixedIntegerRounding2.hpp"
67#include "CglTwomir.hpp"
68#include "CglDuplicateRow.hpp"
69#include "CglStored.hpp"
70#include "CglLandP.hpp"
71
72#include "CbcModel.hpp"
73#include "CbcHeuristic.hpp"
74#include "CbcHeuristicLocal.hpp"
75#include "CbcHeuristicGreedy.hpp"
76#include "CbcHeuristicFPump.hpp"
77#include "CbcTreeLocal.hpp"
78#include "CbcCompareActual.hpp"
79#include "CbcBranchActual.hpp"
80#include  "CbcOrClpParam.hpp"
81#include  "CbcCutGenerator.hpp"
82#include  "CbcStrategy.hpp"
83
84#include "OsiClpSolverInterface.hpp"
85#ifdef COIN_HAS_ASL
86#include "Cbc_ampl.h"
87static bool usingAmpl=false;
88#endif
89static double totalTime=0.0;
90static void statistics(ClpSimplex * originalModel, ClpSimplex * model);
91static bool maskMatches(const int * starts, char ** masks,
92                        std::string & check);
93static void generateCode(CbcModel * model, const char * fileName,int type,int preProcess);
94#ifdef NDEBUG
95#undef NDEBUG
96#endif
97//#############################################################################
98// To use USERCBC uncomment the following define and add in your fake main program here
99//#define USER_HAS_FAKE_MAIN
100//  Start any fake main program
101#ifdef USER_HAS_FAKE_MAIN
102#endif
103//  End any fake main program
104//#############################################################################
105
106// Allow for interrupts
107// But is this threadsafe ? (so switched off by option)
108
109#include "CoinSignal.hpp"
110static CbcModel * currentBranchModel = NULL;
111
112extern "C" {
113   static void signal_handler(int whichSignal)
114   {
115      if (currentBranchModel!=NULL) 
116         currentBranchModel->setMaximumNodes(0); // stop at next node
117      return;
118   }
119}
120
121int mainTest (int argc, const char *argv[],int algorithm,
122              ClpSimplex empty, bool doPresolve,int switchOff,bool doVector);
123void CbcClpUnitTest (const CbcModel & saveModel);
124int CbcOrClpRead_mode=1;
125FILE * CbcOrClpReadCommand=stdin;
126static bool noPrinting=false;
127static int * analyze(OsiClpSolverInterface * solverMod, int & numberChanged, double & increment,
128                     bool changeInt)
129{
130  OsiSolverInterface * solver = solverMod->clone();
131  if (0) {
132    // just get increment
133    CbcModel model(*solver);
134    model.analyzeObjective();
135    double increment2=model.getCutoffIncrement();
136    printf("initial cutoff increment %g\n",increment2);
137  }
138  const double *objective = solver->getObjCoefficients() ;
139  const double *lower = solver->getColLower() ;
140  const double *upper = solver->getColUpper() ;
141  int numberColumns = solver->getNumCols() ;
142  int numberRows = solver->getNumRows();
143  double direction = solver->getObjSense();
144  int iRow,iColumn;
145
146  // Row copy
147  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
148  const double * elementByRow = matrixByRow.getElements();
149  const int * column = matrixByRow.getIndices();
150  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
151  const int * rowLength = matrixByRow.getVectorLengths();
152
153  // Column copy
154  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
155  const double * element = matrixByCol.getElements();
156  const int * row = matrixByCol.getIndices();
157  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
158  const int * columnLength = matrixByCol.getVectorLengths();
159
160  const double * rowLower = solver->getRowLower();
161  const double * rowUpper = solver->getRowUpper();
162
163  char * ignore = new char [numberRows];
164  int * changed = new int[numberColumns];
165  int * which = new int[numberRows];
166  double * changeRhs = new double[numberRows];
167  memset(changeRhs,0,numberRows*sizeof(double));
168  memset(ignore,0,numberRows);
169  numberChanged=0;
170  int numberInteger=0;
171  for (iColumn=0;iColumn<numberColumns;iColumn++) {
172    if (upper[iColumn] > lower[iColumn]+1.0e-8&&solver->isInteger(iColumn)) 
173      numberInteger++;
174  }
175  bool finished=false;
176  while (!finished) {
177    int saveNumberChanged = numberChanged;
178    for (iRow=0;iRow<numberRows;iRow++) {
179      int numberContinuous=0;
180      double value1=0.0,value2=0.0;
181      bool allIntegerCoeff=true;
182      double sumFixed=0.0;
183      int jColumn1=-1,jColumn2=-1;
184      for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
185        int jColumn = column[j];
186        double value = elementByRow[j];
187        if (upper[jColumn] > lower[jColumn]+1.0e-8) {
188          if (!solver->isInteger(jColumn)) {
189            if (numberContinuous==0) {
190              jColumn1=jColumn;
191              value1=value;
192            } else {
193              jColumn2=jColumn;
194              value2=value;
195            }
196            numberContinuous++;
197          } else {
198            if (fabs(value-floor(value+0.5))>1.0e-12)
199              allIntegerCoeff=false;
200          }
201        } else {
202          sumFixed += lower[jColumn]*value;
203        }
204      }
205      double low = rowLower[iRow];
206      if (low>-1.0e20) {
207        low -= sumFixed;
208        if (fabs(low-floor(low+0.5))>1.0e-12)
209          allIntegerCoeff=false;
210      }
211      double up = rowUpper[iRow];
212      if (up<1.0e20) {
213        up -= sumFixed;
214        if (fabs(up-floor(up+0.5))>1.0e-12)
215          allIntegerCoeff=false;
216      }
217      if (!allIntegerCoeff)
218        continue; // can't do
219      if (numberContinuous==1) {
220        // see if really integer
221        // This does not allow for complicated cases
222        if (low==up) {
223          if (fabs(value1)>1.0e-3) {
224            value1 = 1.0/value1;
225            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
226              // integer
227              changed[numberChanged++]=jColumn1;
228              solver->setInteger(jColumn1);
229              if (upper[jColumn1]>1.0e20)
230                solver->setColUpper(jColumn1,1.0e20);
231              if (lower[jColumn1]<-1.0e20)
232                solver->setColLower(jColumn1,-1.0e20);
233            }
234          }
235        } else {
236          if (fabs(value1)>1.0e-3) {
237            value1 = 1.0/value1;
238            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
239              // This constraint will not stop it being integer
240              ignore[iRow]=1;
241            }
242          }
243        }
244      } else if (numberContinuous==2) {
245        if (low==up) {
246          /* need general theory - for now just look at 2 cases -
247             1 - +- 1 one in column and just costs i.e. matching objective
248             2 - +- 1 two in column but feeds into G/L row which will try and minimize
249          */
250          if (fabs(value1)==1.0&&value1*value2==-1.0&&!lower[jColumn1]
251              &&!lower[jColumn2]) {
252            int n=0;
253            int i;
254            double objChange=direction*(objective[jColumn1]+objective[jColumn2]);
255            double bound = CoinMin(upper[jColumn1],upper[jColumn2]);
256            bound = CoinMin(bound,1.0e20);
257            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
258              int jRow = row[i];
259              double value = element[i];
260              if (jRow!=iRow) {
261                which[n++]=jRow;
262                changeRhs[jRow]=value;
263              }
264            }
265            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
266              int jRow = row[i];
267              double value = element[i];
268              if (jRow!=iRow) {
269                if (!changeRhs[jRow]) {
270                  which[n++]=jRow;
271                  changeRhs[jRow]=value;
272                } else {
273                  changeRhs[jRow]+=value;
274                }
275              }
276            }
277            if (objChange>=0.0) {
278              // see if all rows OK
279              bool good=true;
280              for (i=0;i<n;i++) {
281                int jRow = which[i];
282                double value = changeRhs[jRow];
283                if (value) {
284                  value *= bound;
285                  if (rowLength[jRow]==1) {
286                    if (value>0.0) {
287                      double rhs = rowLower[jRow];
288                      if (rhs>0.0) {
289                        double ratio =rhs/value;
290                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
291                          good=false;
292                      }
293                    } else {
294                      double rhs = rowUpper[jRow];
295                      if (rhs<0.0) {
296                        double ratio =rhs/value;
297                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
298                          good=false;
299                      }
300                    }
301                  } else if (rowLength[jRow]==2) {
302                    if (value>0.0) {
303                      if (rowLower[jRow]>-1.0e20)
304                        good=false;
305                    } else {
306                      if (rowUpper[jRow]<1.0e20)
307                        good=false;
308                    }
309                  } else {
310                    good=false;
311                  }
312                }
313              }
314              if (good) {
315                // both can be integer
316                changed[numberChanged++]=jColumn1;
317                solver->setInteger(jColumn1);
318                if (upper[jColumn1]>1.0e20)
319                  solver->setColUpper(jColumn1,1.0e20);
320                if (lower[jColumn1]<-1.0e20)
321                  solver->setColLower(jColumn1,-1.0e20);
322                changed[numberChanged++]=jColumn2;
323                solver->setInteger(jColumn2);
324                if (upper[jColumn2]>1.0e20)
325                  solver->setColUpper(jColumn2,1.0e20);
326                if (lower[jColumn2]<-1.0e20)
327                  solver->setColLower(jColumn2,-1.0e20);
328              }
329            }
330            // clear
331            for (i=0;i<n;i++) {
332              changeRhs[which[i]]=0.0;
333            }
334          }
335        }
336      }
337    }
338    for (iColumn=0;iColumn<numberColumns;iColumn++) {
339      if (upper[iColumn] > lower[iColumn]+1.0e-8&&!solver->isInteger(iColumn)) {
340        double value;
341        value = upper[iColumn];
342        if (value<1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
343          continue;
344        value = lower[iColumn];
345        if (value>-1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
346          continue;
347        bool integer=true;
348        for (CoinBigIndex j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
349          int iRow = row[j];
350          if (!ignore[iRow]) {
351            integer=false;
352            break;
353          }
354        }
355        if (integer) {
356          // integer
357          changed[numberChanged++]=iColumn;
358          solver->setInteger(iColumn);
359          if (upper[iColumn]>1.0e20)
360            solver->setColUpper(iColumn,1.0e20);
361          if (lower[iColumn]<-1.0e20)
362            solver->setColLower(iColumn,-1.0e20);
363        }
364      }
365    }
366    finished = numberChanged==saveNumberChanged;
367  }
368  delete [] which;
369  delete [] changeRhs;
370  delete [] ignore;
371  if (numberInteger&&!noPrinting)
372    printf("%d integer variables",numberInteger);
373  if (changeInt) {
374    if (!noPrinting) {
375      if (numberChanged)
376        printf(" and %d variables made integer\n",numberChanged);
377      else
378        printf("\n");
379    }
380    delete [] ignore;
381    //increment=0.0;
382    if (!numberChanged) {
383      delete [] changed;
384      delete solver;
385      return NULL;
386    } else {
387      for (iColumn=0;iColumn<numberColumns;iColumn++) {
388        if (solver->isInteger(iColumn))
389          solverMod->setInteger(iColumn);
390      }
391      delete solver;
392      return changed;
393    }
394  } else {
395    if (!noPrinting) {
396      if (numberChanged)
397        printf(" and %d variables could be made integer\n",numberChanged);
398      else
399        printf("\n");
400    }
401    // just get increment
402    CbcModel model(*solver);
403    if (noPrinting)
404      model.setLogLevel(0);
405    model.analyzeObjective();
406    double increment2=model.getCutoffIncrement();
407    if (increment2>increment) {
408      if (!noPrinting)
409        printf("cutoff increment increased from %g to %g\n",increment,increment2);
410      increment=increment2;
411    }
412    delete solver;
413    numberChanged=0;
414    delete [] changed;
415    return NULL;
416  }
417}
418static int outDupRow(OsiSolverInterface * solver) 
419{
420  CglDuplicateRow dupCuts(solver);
421  CglTreeInfo info;
422  info.level = 0;
423  info.pass = 0;
424  int numberRows = solver->getNumRows();
425  info.formulation_rows = numberRows;
426  info.inTree = false;
427  info.strengthenRow= NULL;
428  info.pass = 0;
429  OsiCuts cs;
430  dupCuts.generateCuts(*solver,cs,info);
431  const int * duplicate = dupCuts.duplicate();
432  // Get rid of duplicate rows
433  int * which = new int[numberRows]; 
434  int numberDrop=0;
435  for (int iRow=0;iRow<numberRows;iRow++) {
436    if (duplicate[iRow]==-2||duplicate[iRow]>=0) 
437      which[numberDrop++]=iRow;
438  }
439  if (numberDrop) {
440    solver->deleteRows(numberDrop,which);
441  }
442  delete [] which;
443  // see if we have any column cuts
444  int numberColumnCuts = cs.sizeColCuts() ;
445  const double * columnLower = solver->getColLower();
446  const double * columnUpper = solver->getColUpper();
447  for (int k = 0;k<numberColumnCuts;k++) {
448    OsiColCut * thisCut = cs.colCutPtr(k) ;
449    const CoinPackedVector & lbs = thisCut->lbs() ;
450    const CoinPackedVector & ubs = thisCut->ubs() ;
451    int j ;
452    int n ;
453    const int * which ;
454    const double * values ;
455    n = lbs.getNumElements() ;
456    which = lbs.getIndices() ;
457    values = lbs.getElements() ;
458    for (j = 0;j<n;j++) {
459      int iColumn = which[j] ;
460      if (values[j]>columnLower[iColumn]) 
461        solver->setColLower(iColumn,values[j]) ;
462    }
463    n = ubs.getNumElements() ;
464    which = ubs.getIndices() ;
465    values = ubs.getElements() ;
466    for (j = 0;j<n;j++) {
467      int iColumn = which[j] ;
468      if (values[j]<columnUpper[iColumn]) 
469        solver->setColUpper(iColumn,values[j]) ;
470    }
471  }
472  return numberDrop;
473}
474void checkSOS(CbcModel * babModel, const OsiSolverInterface * solver)
475{
476#ifdef COIN_DEVELOP
477  //const double *objective = solver->getObjCoefficients() ;
478  const double *columnLower = solver->getColLower() ;
479  const double * columnUpper = solver->getColUpper() ;
480  const double * solution = solver->getColSolution();
481  //int numberColumns = solver->getNumCols() ;
482  //int numberRows = solver->getNumRows();
483  //double direction = solver->getObjSense();
484  //int iRow,iColumn;
485
486  // Row copy
487  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
488  //const double * elementByRow = matrixByRow.getElements();
489  //const int * column = matrixByRow.getIndices();
490  //const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
491  const int * rowLength = matrixByRow.getVectorLengths();
492
493  // Column copy
494  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
495  const double * element = matrixByCol.getElements();
496  const int * row = matrixByCol.getIndices();
497  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
498  const int * columnLength = matrixByCol.getVectorLengths();
499
500  const double * rowLower = solver->getRowLower();
501  const double * rowUpper = solver->getRowUpper();
502  OsiObject ** objects = babModel->objects();
503  int numberObjects = babModel->numberObjects();
504  for (int iObj = 0;iObj<numberObjects;iObj++) {
505    CbcSOS * objSOS =
506      dynamic_cast <CbcSOS *>(objects[iObj]) ;
507    if (objSOS) {
508      int n=objSOS->numberMembers();
509      const int * which = objSOS->members();
510      const double * weight = objSOS->weights();
511      int type = objSOS->sosType();
512      // convexity row?
513      int iColumn;
514      iColumn=which[0];
515      int j;
516      int convex=-1;
517      for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
518        int iRow = row[j];
519        double value = element[j];
520        if (rowLower[iRow]==1.0&&rowUpper[iRow]==1.0&&
521            value==1.0) {
522          // possible
523          if (rowLength[iRow]==n) {
524            if (convex==-1)
525              convex=iRow;
526            else
527              convex=-2;
528          }
529        }
530      }
531      printf ("set %d of type %d has %d members - possible convexity row %d\n",
532              iObj,type,n,convex);
533      for (int i=0;i<n;i++) {
534        iColumn = which[i];
535        int convex2=-1;
536        for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
537          int iRow = row[j];
538          if (iRow==convex) {
539            double value = element[j];
540            if (value==1.0) {
541              convex2=iRow;
542            }
543          }
544        }
545        if (convex2<0&&convex>=0) {
546          printf("odd convexity row\n");
547          convex=-2;
548        }
549        printf("col %d has weight %g and value %g, bounds %g %g\n",
550               iColumn,weight[i],solution[iColumn],columnLower[iColumn],
551               columnUpper[iColumn]);
552      }
553    }
554  }
555#endif
556}
557int main (int argc, const char *argv[])
558{
559  /* Note
560     This is meant as a stand-alone executable to do as much of coin as possible.
561     It should only have one solver known to it.
562  */
563  {
564    double time1 = CoinCpuTime(),time2;
565    bool goodModel=false;
566    CoinSighandler_t saveSignal=SIG_DFL;
567    // register signal handler
568    saveSignal = signal(SIGINT,signal_handler);
569    // Set up all non-standard stuff
570    OsiClpSolverInterface solver1;
571    CbcModel model(solver1);
572    CbcModel * babModel = NULL;
573    model.setNumberBeforeTrust(21);
574    int cutPass=-1234567;
575    int tunePreProcess=5;
576    OsiSolverInterface * solver = model.solver();
577    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
578    ClpSimplex * lpSolver = clpSolver->getModelPtr();
579    clpSolver->messageHandler()->setLogLevel(0) ;
580    model.messageHandler()->setLogLevel(1);
581    // For priorities etc
582    int * priorities=NULL;
583    int * branchDirection=NULL;
584    double * pseudoDown=NULL;
585    double * pseudoUp=NULL;
586    double * solutionIn = NULL;
587    int * prioritiesIn = NULL;
588    int numberSOS = 0;
589    int * sosStart = NULL;
590    int * sosIndices = NULL;
591    char * sosType = NULL;
592    double * sosReference = NULL;
593    int * sosPriority=NULL;
594#ifdef COIN_HAS_ASL
595    ampl_info info;
596    memset(&info,0,sizeof(info));
597    if (argc>2&&!strcmp(argv[2],"-AMPL")) {
598      usingAmpl=true;
599      int returnCode = readAmpl(&info,argc,const_cast<char **>(argv));
600      if (returnCode)
601        return returnCode;
602      CbcOrClpRead_mode=2; // so will start with parameters
603      // see if log in list
604      noPrinting=true;
605      for (int i=1;i<info.numberArguments;i++) {
606        if (!strcmp(info.arguments[i],"log")) {
607          if (i<info.numberArguments-1&&atoi(info.arguments[i+1])>0)
608            noPrinting=false;
609          break;
610        }
611      }
612      if (noPrinting) {
613        model.messageHandler()->setLogLevel(0);
614        setCbcOrClpPrinting(false);
615      }
616      if (!noPrinting)
617        printf("%d rows, %d columns and %d elements\n",
618               info.numberRows,info.numberColumns,info.numberElements);
619      solver->loadProblem(info.numberColumns,info.numberRows,info.starts,
620                          info.rows,info.elements,
621                          info.columnLower,info.columnUpper,info.objective,
622                          info.rowLower,info.rowUpper);
623      // If we had a solution use it
624      if (info.primalSolution) {
625        solver->setColSolution(info.primalSolution);
626      }
627      // status
628      if (info.rowStatus) {
629        unsigned char * statusArray = lpSolver->statusArray();
630        int i;
631        for (i=0;i<info.numberColumns;i++)
632          statusArray[i]=(char)info.columnStatus[i];
633        statusArray+=info.numberColumns;
634        for (i=0;i<info.numberRows;i++)
635          statusArray[i]=(char)info.rowStatus[i];
636        CoinWarmStartBasis * basis = lpSolver->getBasis();
637        solver->setWarmStart(basis);
638        delete basis;
639      }
640      freeArrays1(&info);
641      // modify objective if necessary
642      solver->setObjSense(info.direction);
643      solver->setDblParam(OsiObjOffset,info.offset);
644      // Set integer variables
645      for (int i=info.numberColumns-info.numberBinary-info.numberIntegers;
646           i<info.numberColumns;i++)
647        solver->setInteger(i);
648      goodModel=true;
649      // change argc etc
650      argc = info.numberArguments;
651      argv = const_cast<const char **>(info.arguments);
652    }
653#endif   
654    // default action on import
655    int allowImportErrors=0;
656    int keepImportNames=1;
657    int doIdiot=-1;
658    int outputFormat=2;
659    int slpValue=-1;
660    int cppValue=-1;
661    int printOptions=0;
662    int printMode=0;
663    int presolveOptions=0;
664    int substitution=3;
665    int dualize=0;
666    int doCrash=0;
667    int doVector=0;
668    int doSprint=-1;
669    int doScaling=1;
670    // set reasonable defaults
671    int preSolve=5;
672    int preProcess=1;
673    bool useStrategy=false;
674    bool preSolveFile=false;
675   
676    double djFix=1.0e100;
677    double gapRatio=1.0e100;
678    double tightenFactor=0.0;
679    lpSolver->setPerturbation(50);
680    lpSolver->messageHandler()->setPrefix(false);
681    const char dirsep =  CoinFindDirSeparator();
682    std::string directory = (dirsep == '/' ? "./" : ".\\");
683    std::string defaultDirectory = directory;
684    std::string importFile ="";
685    std::string exportFile ="default.mps";
686    std::string importBasisFile ="";
687    std::string importPriorityFile ="";
688    std::string debugFile="";
689    std::string printMask="";
690    double * debugValues = NULL;
691    int numberDebugValues = -1;
692    int basisHasValues=0;
693    std::string exportBasisFile ="default.bas";
694    std::string saveFile ="default.prob";
695    std::string restoreFile ="default.prob";
696    std::string solutionFile ="stdout";
697    std::string solutionSaveFile ="solution.file";
698#define CBCMAXPARAMETERS 200
699    CbcOrClpParam parameters[CBCMAXPARAMETERS];
700    int numberParameters ;
701    establishParams(numberParameters,parameters) ;
702    parameters[whichParam(BASISIN,numberParameters,parameters)].setStringValue(importBasisFile);
703    parameters[whichParam(PRIORITYIN,numberParameters,parameters)].setStringValue(importPriorityFile);
704    parameters[whichParam(BASISOUT,numberParameters,parameters)].setStringValue(exportBasisFile);
705    parameters[whichParam(DEBUG,numberParameters,parameters)].setStringValue(debugFile);
706    parameters[whichParam(PRINTMASK,numberParameters,parameters)].setStringValue(printMask);
707    parameters[whichParam(DIRECTORY,numberParameters,parameters)].setStringValue(directory);
708    parameters[whichParam(DUALBOUND,numberParameters,parameters)].setDoubleValue(lpSolver->dualBound());
709    parameters[whichParam(DUALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->dualTolerance());
710    parameters[whichParam(EXPORT,numberParameters,parameters)].setStringValue(exportFile);
711    parameters[whichParam(IDIOT,numberParameters,parameters)].setIntValue(doIdiot);
712    parameters[whichParam(IMPORT,numberParameters,parameters)].setStringValue(importFile);
713    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].setDoubleValue(1.0e-8);
714    int slog = whichParam(SOLVERLOGLEVEL,numberParameters,parameters);
715    int log = whichParam(LOGLEVEL,numberParameters,parameters);
716    parameters[slog].setIntValue(0);
717    parameters[log].setIntValue(1);
718    parameters[whichParam(MAXFACTOR,numberParameters,parameters)].setIntValue(lpSolver->factorizationFrequency());
719    parameters[whichParam(MAXITERATION,numberParameters,parameters)].setIntValue(lpSolver->maximumIterations());
720    parameters[whichParam(OUTPUTFORMAT,numberParameters,parameters)].setIntValue(outputFormat);
721    parameters[whichParam(PRESOLVEPASS,numberParameters,parameters)].setIntValue(preSolve);
722    parameters[whichParam(PERTVALUE,numberParameters,parameters)].setIntValue(lpSolver->perturbation());
723    parameters[whichParam(PRIMALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->primalTolerance());
724    parameters[whichParam(PRIMALWEIGHT,numberParameters,parameters)].setDoubleValue(lpSolver->infeasibilityCost());
725    parameters[whichParam(RESTORE,numberParameters,parameters)].setStringValue(restoreFile);
726    parameters[whichParam(SAVE,numberParameters,parameters)].setStringValue(saveFile);
727    //parameters[whichParam(TIMELIMIT,numberParameters,parameters)].setDoubleValue(1.0e8);
728    parameters[whichParam(TIMELIMIT_BAB,numberParameters,parameters)].setDoubleValue(1.0e8);
729    parameters[whichParam(SOLUTION,numberParameters,parameters)].setStringValue(solutionFile);
730    parameters[whichParam(SAVESOL,numberParameters,parameters)].setStringValue(solutionSaveFile);
731    parameters[whichParam(SPRINT,numberParameters,parameters)].setIntValue(doSprint);
732    parameters[whichParam(SUBSTITUTION,numberParameters,parameters)].setIntValue(substitution);
733    parameters[whichParam(DUALIZE,numberParameters,parameters)].setIntValue(dualize);
734    model.setNumberBeforeTrust(5);
735    parameters[whichParam(NUMBERBEFORE,numberParameters,parameters)].setIntValue(5);
736    parameters[whichParam(MAXNODES,numberParameters,parameters)].setIntValue(model.getMaximumNodes());
737    model.setNumberStrong(5);
738    parameters[whichParam(STRONGBRANCHING,numberParameters,parameters)].setIntValue(model.numberStrong());
739    parameters[whichParam(INFEASIBILITYWEIGHT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcInfeasibilityWeight));
740    parameters[whichParam(INTEGERTOLERANCE,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcIntegerTolerance));
741    parameters[whichParam(INCREMENT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcCutoffIncrement));
742    // Set up likely cut generators and defaults
743    parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("on");
744    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(128|64|1);
745    parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].setIntValue(-1);
746    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].setIntValue(100);
747    parameters[whichParam(CUTSSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
748    parameters[whichParam(HEURISTICSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
749    parameters[whichParam(NODESTRATEGY,numberParameters,parameters)].setCurrentOption("fewest");
750    int nodeStrategy=0;
751    int doSOS=1;
752    int verbose=0;
753    CglGomory gomoryGen;
754    // try larger limit
755    gomoryGen.setLimitAtRoot(512);
756    gomoryGen.setLimit(50);
757    // set default action (0=off,1=on,2=root)
758    int gomoryAction=3;
759    parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
760
761    CglProbing probingGen;
762    probingGen.setUsingObjective(true);
763    probingGen.setMaxPass(3);
764    probingGen.setMaxPassRoot(3);
765    // Number of unsatisfied variables to look at
766    probingGen.setMaxProbe(10);
767    probingGen.setMaxProbeRoot(50);
768    // How far to follow the consequences
769    probingGen.setMaxLook(10);
770    probingGen.setMaxLookRoot(50);
771    probingGen.setMaxLookRoot(10);
772    // Only look at rows with fewer than this number of elements
773    probingGen.setMaxElements(200);
774    probingGen.setRowCuts(3);
775    // set default action (0=off,1=on,2=root)
776    int probingAction=1;
777    parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
778
779    CglKnapsackCover knapsackGen;
780    //knapsackGen.switchOnExpensive();
781    // set default action (0=off,1=on,2=root)
782    int knapsackAction=3;
783    parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
784
785    CglRedSplit redsplitGen;
786    //redsplitGen.setLimit(100);
787    // set default action (0=off,1=on,2=root)
788    // Off as seems to give some bad cuts
789    int redsplitAction=2;
790    parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption("root");
791
792    CglClique cliqueGen(false,true);
793    cliqueGen.setStarCliqueReport(false);
794    cliqueGen.setRowCliqueReport(false);
795    cliqueGen.setMinViolation(0.1);
796    // set default action (0=off,1=on,2=root)
797    int cliqueAction=3;
798    parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption("ifmove");
799
800    CglMixedIntegerRounding2 mixedGen;
801    // set default action (0=off,1=on,2=root)
802    int mixedAction=3;
803    parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
804
805    CglFlowCover flowGen;
806    // set default action (0=off,1=on,2=root)
807    int flowAction=3;
808    parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
809
810    CglTwomir twomirGen;
811    twomirGen.setMaxElements(250);
812    // set default action (0=off,1=on,2=root)
813    int twomirAction=2;
814    parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption("root");
815    CglLandP landpGen;
816    // set default action (0=off,1=on,2=root)
817    int landpAction=0;
818    parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption("off");
819    // Stored cuts
820    bool storedCuts = false;
821
822    bool useRounding=true;
823    parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption("on");
824    bool useFpump=true;
825    parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption("on");
826    bool useGreedy=true;
827    parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption("on");
828    bool useCombine=true;
829    parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
830    bool useLocalTree=false;
831    parameters[whichParam(COSTSTRATEGY,numberParameters,parameters)].setCurrentOption("off");
832    int useCosts=0;
833    // don't use input solution
834    int useSolution=0;
835   
836    // total number of commands read
837    int numberGoodCommands=0;
838    // Set false if user does anything advanced
839    bool defaultSettings=true;
840
841    // Hidden stuff for barrier
842    int choleskyType = 0;
843    int gamma=0;
844    int scaleBarrier=0;
845    int doKKT=0;
846    int crossover=2; // do crossover unless quadratic
847    // For names
848    int lengthName = 0;
849    std::vector<std::string> rowNames;
850    std::vector<std::string> columnNames;
851   
852    std::string field;
853    if (!noPrinting) {
854      std::cout<<"Coin Cbc and Clp Solver version "<<CBCVERSION
855               <<", build "<<__DATE__<<std::endl;
856      // Print command line
857      if (argc>1) {
858        printf("command line - ");
859        for (int i=0;i<argc;i++)
860          printf("%s ",argv[i]);
861        printf("\n");
862      }
863    }
864    while (1) {
865      // next command
866      field=CoinReadGetCommand(argc,argv);
867      // exit if null or similar
868      if (!field.length()) {
869        if (numberGoodCommands==1&&goodModel) {
870          // we just had file name - do branch and bound
871          field="branch";
872        } else if (!numberGoodCommands) {
873          // let's give the sucker a hint
874          std::cout
875            <<"CoinSolver takes input from arguments ( - switches to stdin)"
876            <<std::endl
877            <<"Enter ? for list of commands or help"<<std::endl;
878          field="-";
879        } else {
880          break;
881        }
882      }
883     
884      // see if ? at end
885      int numberQuery=0;
886      if (field!="?"&&field!="???") {
887        int length = field.length();
888        int i;
889        for (i=length-1;i>0;i--) {
890          if (field[i]=='?') 
891            numberQuery++;
892          else
893            break;
894        }
895        field=field.substr(0,length-numberQuery);
896      }
897      // find out if valid command
898      int iParam;
899      int numberMatches=0;
900      int firstMatch=-1;
901      for ( iParam=0; iParam<numberParameters; iParam++ ) {
902        int match = parameters[iParam].matches(field);
903        if (match==1) {
904          numberMatches = 1;
905          firstMatch=iParam;
906          break;
907        } else {
908          if (match&&firstMatch<0)
909            firstMatch=iParam;
910          numberMatches += match>>1;
911        }
912      }
913      if (iParam<numberParameters&&!numberQuery) {
914        // found
915        CbcOrClpParam found = parameters[iParam];
916        CbcOrClpParameterType type = found.type();
917        int valid;
918        numberGoodCommands++;
919        if (type==BAB&&goodModel) {
920          // check if any integers
921#ifdef COIN_HAS_ASL
922          if (info.numberSos&&doSOS&&usingAmpl) {
923            // SOS
924            numberSOS = info.numberSos;
925          }
926#endif
927          if (!lpSolver->integerInformation()&&!numberSOS&&
928              !clpSolver->numberSOS())
929            type=DUALSIMPLEX;
930        }
931        if (type==GENERALQUERY) {
932          bool evenHidden=false;
933          if ((verbose&8)!=0) {
934            // even hidden
935            evenHidden = true;
936            verbose &= ~8;
937          }
938#ifdef COIN_HAS_ASL
939          if (verbose<4&&usingAmpl)
940            verbose +=4;
941#endif
942          if (verbose<4) {
943            std::cout<<"In argument list keywords have leading - "
944              ", -stdin or just - switches to stdin"<<std::endl;
945            std::cout<<"One command per line (and no -)"<<std::endl;
946            std::cout<<"abcd? gives list of possibilities, if only one + explanation"<<std::endl;
947            std::cout<<"abcd?? adds explanation, if only one fuller help"<<std::endl;
948            std::cout<<"abcd without value (where expected) gives current value"<<std::endl;
949            std::cout<<"abcd value sets value"<<std::endl;
950            std::cout<<"Commands are:"<<std::endl;
951          } else {
952            std::cout<<"Cbc options are set within AMPL with commands like:"<<std::endl<<std::endl;
953            std::cout<<"         option cbc_options \"cuts=root log=2 feas=on slog=1\""<<std::endl<<std::endl;
954            std::cout<<"only maximize, dual, primal, help and quit are recognized without ="<<std::endl;
955          }
956          int maxAcross=5;
957          if ((verbose%4)!=0)
958            maxAcross=1;
959          int limits[]={1,51,101,151,201,251,301,351,401};
960          std::vector<std::string> types;
961          types.push_back("Double parameters:");
962          types.push_back("Branch and Cut double parameters:");
963          types.push_back("Integer parameters:");
964          types.push_back("Branch and Cut integer parameters:");
965          types.push_back("Keyword parameters:");
966          types.push_back("Branch and Cut keyword parameters:");
967          types.push_back("Actions or string parameters:");
968          types.push_back("Branch and Cut actions:");
969          int iType;
970          for (iType=0;iType<8;iType++) {
971            int across=0;
972            if ((verbose%4)!=0)
973              std::cout<<std::endl;
974            std::cout<<types[iType]<<std::endl;
975            if ((verbose&2)!=0)
976              std::cout<<std::endl;
977            for ( iParam=0; iParam<numberParameters; iParam++ ) {
978              int type = parameters[iParam].type();
979              if ((parameters[iParam].displayThis()||evenHidden)&&
980                  type>=limits[iType]
981                  &&type<limits[iType+1]) {
982                // but skip if not useful for ampl (and in ampl mode)
983                if (verbose>=4&&(parameters[iParam].whereUsed()&4)==0)
984                  continue;
985                if (!across) {
986                  if ((verbose&2)==0) 
987                    std::cout<<"  ";
988                  else
989                    std::cout<<"Command ";
990                }
991                std::cout<<parameters[iParam].matchName()<<"  ";
992                across++;
993                if (across==maxAcross) {
994                  across=0;
995                  if ((verbose%4)!=0) {
996                    // put out description as well
997                    if ((verbose&1)!=0) 
998                      std::cout<<parameters[iParam].shortHelp();
999                    std::cout<<std::endl;
1000                    if ((verbose&2)!=0) {
1001                      std::cout<<"---- description"<<std::endl;
1002                      parameters[iParam].printLongHelp();
1003                      std::cout<<"----"<<std::endl<<std::endl;
1004                    }
1005                  } else {
1006                    std::cout<<std::endl;
1007                  }
1008                }
1009              }
1010            }
1011            if (across)
1012              std::cout<<std::endl;
1013          }
1014        } else if (type==FULLGENERALQUERY) {
1015          std::cout<<"Full list of commands is:"<<std::endl;
1016          int maxAcross=5;
1017          int limits[]={1,51,101,151,201,251,301,351,401};
1018          std::vector<std::string> types;
1019          types.push_back("Double parameters:");
1020          types.push_back("Branch and Cut double parameters:");
1021          types.push_back("Integer parameters:");
1022          types.push_back("Branch and Cut integer parameters:");
1023          types.push_back("Keyword parameters:");
1024          types.push_back("Branch and Cut keyword parameters:");
1025          types.push_back("Actions or string parameters:");
1026          types.push_back("Branch and Cut actions:");
1027          int iType;
1028          for (iType=0;iType<8;iType++) {
1029            int across=0;
1030            std::cout<<types[iType]<<"  ";
1031            for ( iParam=0; iParam<numberParameters; iParam++ ) {
1032              int type = parameters[iParam].type();
1033              if (type>=limits[iType]
1034                  &&type<limits[iType+1]) {
1035                if (!across)
1036                  std::cout<<"  ";
1037                std::cout<<parameters[iParam].matchName()<<"  ";
1038                across++;
1039                if (across==maxAcross) {
1040                  std::cout<<std::endl;
1041                  across=0;
1042                }
1043              }
1044            }
1045            if (across)
1046              std::cout<<std::endl;
1047          }
1048        } else if (type<101) {
1049          // get next field as double
1050          double value = CoinReadGetDoubleField(argc,argv,&valid);
1051          if (!valid) {
1052            if (type<51) {
1053              parameters[iParam].setDoubleParameter(lpSolver,value);
1054            } else if (type<81) {
1055              parameters[iParam].setDoubleParameter(model,value);
1056            } else {
1057              parameters[iParam].setDoubleParameter(lpSolver,value);
1058              switch(type) {
1059              case DJFIX:
1060                djFix=value;
1061                preSolve=5;
1062                defaultSettings=false; // user knows what she is doing
1063                break;
1064              case GAPRATIO:
1065                gapRatio=value;
1066                break;
1067              case TIGHTENFACTOR:
1068                tightenFactor=value;
1069                defaultSettings=false; // user knows what she is doing
1070                break;
1071              default:
1072                abort();
1073              }
1074            }
1075          } else if (valid==1) {
1076            abort();
1077          } else {
1078            std::cout<<parameters[iParam].name()<<" has value "<<
1079              parameters[iParam].doubleValue()<<std::endl;
1080          }
1081        } else if (type<201) {
1082          // get next field as int
1083          int value = CoinReadGetIntField(argc,argv,&valid);
1084          if (!valid) {
1085            if (type<151) {
1086              if (parameters[iParam].type()==PRESOLVEPASS)
1087                preSolve = value;
1088              else if (parameters[iParam].type()==IDIOT)
1089                doIdiot = value;
1090              else if (parameters[iParam].type()==SPRINT)
1091                doSprint = value;
1092              else if (parameters[iParam].type()==OUTPUTFORMAT)
1093                outputFormat = value;
1094              else if (parameters[iParam].type()==SLPVALUE)
1095                slpValue = value;
1096              else if (parameters[iParam].type()==CPP)
1097                cppValue = value;
1098              else if (parameters[iParam].type()==PRESOLVEOPTIONS)
1099                presolveOptions = value;
1100              else if (parameters[iParam].type()==PRINTOPTIONS)
1101                printOptions = value;
1102              else if (parameters[iParam].type()==SUBSTITUTION)
1103                substitution = value;
1104              else if (parameters[iParam].type()==DUALIZE)
1105                dualize = value;
1106              else if (parameters[iParam].type()==CUTPASS)
1107                cutPass = value;
1108              else if (parameters[iParam].type()==PROCESSTUNE)
1109                tunePreProcess = value;
1110              else if (parameters[iParam].type()==VERBOSE)
1111                verbose = value;
1112              else if (parameters[iParam].type()==FPUMPITS)
1113                { useFpump = true;parameters[iParam].setIntValue(value);}
1114              parameters[iParam].setIntParameter(lpSolver,value);
1115            } else {
1116              parameters[iParam].setIntParameter(model,value);
1117            }
1118          } else if (valid==1) {
1119            abort();
1120          } else {
1121            std::cout<<parameters[iParam].name()<<" has value "<<
1122              parameters[iParam].intValue()<<std::endl;
1123          }
1124        } else if (type<301) {
1125          // one of several strings
1126          std::string value = CoinReadGetString(argc,argv);
1127          int action = parameters[iParam].parameterOption(value);
1128          if (action<0) {
1129            if (value!="EOL") {
1130              // no match
1131              parameters[iParam].printOptions();
1132            } else {
1133              // print current value
1134              std::cout<<parameters[iParam].name()<<" has value "<<
1135                parameters[iParam].currentOption()<<std::endl;
1136            }
1137          } else {
1138            parameters[iParam].setCurrentOption(action,!noPrinting);
1139            // for now hard wired
1140            switch (type) {
1141            case DIRECTION:
1142              if (action==0)
1143                lpSolver->setOptimizationDirection(1);
1144              else if (action==1)
1145                lpSolver->setOptimizationDirection(-1);
1146              else
1147                lpSolver->setOptimizationDirection(0);
1148              break;
1149            case DUALPIVOT:
1150              if (action==0) {
1151                ClpDualRowSteepest steep(3);
1152                lpSolver->setDualRowPivotAlgorithm(steep);
1153              } else if (action==1) {
1154                //ClpDualRowDantzig dantzig;
1155                ClpDualRowSteepest dantzig(5);
1156                lpSolver->setDualRowPivotAlgorithm(dantzig);
1157              } else if (action==2) {
1158                // partial steep
1159                ClpDualRowSteepest steep(2);
1160                lpSolver->setDualRowPivotAlgorithm(steep);
1161              } else {
1162                ClpDualRowSteepest steep;
1163                lpSolver->setDualRowPivotAlgorithm(steep);
1164              }
1165              break;
1166            case PRIMALPIVOT:
1167              if (action==0) {
1168                ClpPrimalColumnSteepest steep(3);
1169                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1170              } else if (action==1) {
1171                ClpPrimalColumnSteepest steep(0);
1172                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1173              } else if (action==2) {
1174                ClpPrimalColumnDantzig dantzig;
1175                lpSolver->setPrimalColumnPivotAlgorithm(dantzig);
1176              } else if (action==3) {
1177                ClpPrimalColumnSteepest steep(2);
1178                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1179              } else if (action==4) {
1180                ClpPrimalColumnSteepest steep(1);
1181                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1182              } else if (action==5) {
1183                ClpPrimalColumnSteepest steep(4);
1184                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1185              } else if (action==6) {
1186                ClpPrimalColumnSteepest steep(10);
1187                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1188              }
1189              break;
1190            case SCALING:
1191              lpSolver->scaling(action);
1192              solver->setHintParam(OsiDoScale,action!=0,OsiHintTry);
1193              doScaling = 1-action;
1194              break;
1195            case AUTOSCALE:
1196              lpSolver->setAutomaticScaling(action!=0);
1197              break;
1198            case SPARSEFACTOR:
1199              lpSolver->setSparseFactorization((1-action)!=0);
1200              break;
1201            case BIASLU:
1202              lpSolver->factorization()->setBiasLU(action);
1203              break;
1204            case PERTURBATION:
1205              if (action==0)
1206                lpSolver->setPerturbation(50);
1207              else
1208                lpSolver->setPerturbation(100);
1209              break;
1210            case ERRORSALLOWED:
1211              allowImportErrors = action;
1212              break;
1213            case INTPRINT:
1214              printMode=action;
1215              break;
1216              //case ALGORITHM:
1217              //algorithm  = action;
1218              //defaultSettings=false; // user knows what she is doing
1219              //abort();
1220              //break;
1221            case KEEPNAMES:
1222              keepImportNames = 1-action;
1223              break;
1224            case PRESOLVE:
1225              if (action==0)
1226                preSolve = 5;
1227              else if (action==1)
1228                preSolve=0;
1229              else if (action==2)
1230                preSolve=10;
1231              else
1232                preSolveFile=true;
1233              break;
1234            case PFI:
1235              lpSolver->factorization()->setForrestTomlin(action==0);
1236              break;
1237            case CRASH:
1238              doCrash=action;
1239              break;
1240            case VECTOR:
1241              doVector=action;
1242              break;
1243            case MESSAGES:
1244              lpSolver->messageHandler()->setPrefix(action!=0);
1245              break;
1246            case CHOLESKY:
1247              choleskyType = action;
1248              break;
1249            case GAMMA:
1250              gamma=action;
1251              break;
1252            case BARRIERSCALE:
1253              scaleBarrier=action;
1254              break;
1255            case KKT:
1256              doKKT=action;
1257              break;
1258            case CROSSOVER:
1259              crossover=action;
1260              break;
1261            case SOS:
1262              doSOS=action;
1263              break;
1264            case GOMORYCUTS:
1265              defaultSettings=false; // user knows what she is doing
1266              gomoryAction = action;
1267              break;
1268            case PROBINGCUTS:
1269              defaultSettings=false; // user knows what she is doing
1270              probingAction = action;
1271              break;
1272            case KNAPSACKCUTS:
1273              defaultSettings=false; // user knows what she is doing
1274              knapsackAction = action;
1275              break;
1276            case REDSPLITCUTS:
1277              defaultSettings=false; // user knows what she is doing
1278              redsplitAction = action;
1279              break;
1280            case CLIQUECUTS:
1281              defaultSettings=false; // user knows what she is doing
1282              cliqueAction = action;
1283              break;
1284            case FLOWCUTS:
1285              defaultSettings=false; // user knows what she is doing
1286              flowAction = action;
1287              break;
1288            case MIXEDCUTS:
1289              defaultSettings=false; // user knows what she is doing
1290              mixedAction = action;
1291              break;
1292            case TWOMIRCUTS:
1293              defaultSettings=false; // user knows what she is doing
1294              twomirAction = action;
1295              break;
1296            case LANDPCUTS:
1297              defaultSettings=false; // user knows what she is doing
1298              landpAction = action;
1299              break;
1300            case ROUNDING:
1301              defaultSettings=false; // user knows what she is doing
1302              useRounding = action;
1303              break;
1304            case FPUMP:
1305              defaultSettings=false; // user knows what she is doing
1306              useFpump=action;
1307              break;
1308            case CUTSSTRATEGY:
1309              gomoryAction = action;
1310              probingAction = action;
1311              knapsackAction = action;
1312              cliqueAction = action;
1313              flowAction = action;
1314              mixedAction = action;
1315              twomirAction = action;
1316              //landpAction = action;
1317              parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption(action);
1318              parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption(action);
1319              parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption(action);
1320              if (!action) {
1321                redsplitAction = action;
1322                parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption(action);
1323              }
1324              parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption(action);
1325              parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption(action);
1326              parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption(action);
1327              parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption(action);
1328              if (!action) {
1329                landpAction = action;
1330                parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption(action);
1331              }
1332              break;
1333            case HEURISTICSTRATEGY:
1334              useRounding = action;
1335              useGreedy = action;
1336              useCombine = action;
1337              //useLocalTree = action;
1338              useFpump=action;
1339              parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption(action);
1340              parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption(action);
1341              parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption(action);
1342              //parameters[whichParam(LOCALTREE,numberParameters,parameters)].setCurrentOption(action);
1343              parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption(action);
1344              break;
1345            case GREEDY:
1346              defaultSettings=false; // user knows what she is doing
1347              useGreedy = action;
1348              break;
1349            case COMBINE:
1350              defaultSettings=false; // user knows what she is doing
1351              useCombine = action;
1352              break;
1353            case LOCALTREE:
1354              defaultSettings=false; // user knows what she is doing
1355              useLocalTree = action;
1356              break;
1357            case COSTSTRATEGY:
1358              useCosts=action;
1359              break;
1360            case NODESTRATEGY:
1361              nodeStrategy=action;
1362              break;
1363            case PREPROCESS:
1364              preProcess = action;
1365              break;
1366            case USESOLUTION:
1367              useSolution = action;
1368              break;
1369            default:
1370              abort();
1371            }
1372          }
1373        } else {
1374          // action
1375          if (type==EXIT) {
1376#ifdef COIN_HAS_ASL
1377            if(usingAmpl) {
1378              if (info.numberIntegers||info.numberBinary) {
1379                // integer
1380              } else {
1381                // linear
1382              }
1383              writeAmpl(&info);
1384              freeArrays2(&info);
1385              freeArgs(&info);
1386            }
1387#endif
1388            break; // stop all
1389          }
1390          switch (type) {
1391          case DUALSIMPLEX:
1392          case PRIMALSIMPLEX:
1393          case SOLVECONTINUOUS:
1394          case BARRIER:
1395            if (goodModel) {
1396              double objScale = 
1397                parameters[whichParam(OBJSCALE2,numberParameters,parameters)].doubleValue();
1398              if (objScale!=1.0) {
1399                int iColumn;
1400                int numberColumns=lpSolver->numberColumns();
1401                double * dualColumnSolution = 
1402                  lpSolver->dualColumnSolution();
1403                ClpObjective * obj = lpSolver->objectiveAsObject();
1404                assert(dynamic_cast<ClpLinearObjective *> (obj));
1405                double offset;
1406                double * objective = obj->gradient(NULL,NULL,offset,true);
1407                for (iColumn=0;iColumn<numberColumns;iColumn++) {
1408                  dualColumnSolution[iColumn] *= objScale;
1409                  objective[iColumn] *= objScale;;
1410                }
1411                int iRow;
1412                int numberRows=lpSolver->numberRows();
1413                double * dualRowSolution = 
1414                  lpSolver->dualRowSolution();
1415                for (iRow=0;iRow<numberRows;iRow++) 
1416                  dualRowSolution[iRow] *= objScale;
1417                lpSolver->setObjectiveOffset(objScale*lpSolver->objectiveOffset());
1418              }
1419              ClpSolve::SolveType method;
1420              ClpSolve::PresolveType presolveType;
1421              ClpSimplex * model2 = lpSolver;
1422              if (dualize) {
1423                model2 = ((ClpSimplexOther *) model2)->dualOfModel();
1424                printf("Dual of model has %d rows and %d columns\n",
1425                       model2->numberRows(),model2->numberColumns());
1426                model2->setOptimizationDirection(1.0);
1427              }
1428              if (noPrinting)
1429                lpSolver->setLogLevel(0);
1430              ClpSolve solveOptions;
1431              solveOptions.setPresolveActions(presolveOptions);
1432              solveOptions.setSubstitution(substitution);
1433              if (preSolve!=5&&preSolve) {
1434                presolveType=ClpSolve::presolveNumber;
1435                if (preSolve<0) {
1436                  preSolve = - preSolve;
1437                  if (preSolve<=100) {
1438                    presolveType=ClpSolve::presolveNumber;
1439                    printf("Doing %d presolve passes - picking up non-costed slacks\n",
1440                           preSolve);
1441                    solveOptions.setDoSingletonColumn(true);
1442                  } else {
1443                    preSolve -=100;
1444                    presolveType=ClpSolve::presolveNumberCost;
1445                    printf("Doing %d presolve passes - picking up costed slacks\n",
1446                           preSolve);
1447                  }
1448                } 
1449              } else if (preSolve) {
1450                presolveType=ClpSolve::presolveOn;
1451              } else {
1452                presolveType=ClpSolve::presolveOff;
1453              }
1454              solveOptions.setPresolveType(presolveType,preSolve);
1455              if (type==DUALSIMPLEX||type==SOLVECONTINUOUS) {
1456                method=ClpSolve::useDual;
1457              } else if (type==PRIMALSIMPLEX) {
1458                method=ClpSolve::usePrimalorSprint;
1459              } else {
1460                method = ClpSolve::useBarrier;
1461                if (crossover==1) {
1462                  method=ClpSolve::useBarrierNoCross;
1463                } else if (crossover==2) {
1464                  ClpObjective * obj = lpSolver->objectiveAsObject();
1465                  if (obj->type()>1) {
1466                    method=ClpSolve::useBarrierNoCross;
1467                    presolveType=ClpSolve::presolveOff;
1468                    solveOptions.setPresolveType(presolveType,preSolve);
1469                  } 
1470                }
1471              }
1472              solveOptions.setSolveType(method);
1473              if(preSolveFile)
1474                presolveOptions |= 0x40000000;
1475              solveOptions.setSpecialOption(4,presolveOptions);
1476              solveOptions.setSpecialOption(5,printOptions);
1477              if (doVector) {
1478                ClpMatrixBase * matrix = lpSolver->clpMatrix();
1479                if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
1480                  ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
1481                  clpMatrix->makeSpecialColumnCopy();
1482                }
1483              }
1484              if (method==ClpSolve::useDual) {
1485                // dual
1486                if (doCrash)
1487                  solveOptions.setSpecialOption(0,1,doCrash); // crash
1488                else if (doIdiot)
1489                  solveOptions.setSpecialOption(0,2,doIdiot); // possible idiot
1490              } else if (method==ClpSolve::usePrimalorSprint) {
1491                // primal
1492                // if slp turn everything off
1493                if (slpValue>0) {
1494                  doCrash=false;
1495                  doSprint=0;
1496                  doIdiot=-1;
1497                  solveOptions.setSpecialOption(1,10,slpValue); // slp
1498                  method=ClpSolve::usePrimal;
1499                }
1500                if (doCrash) {
1501                  solveOptions.setSpecialOption(1,1,doCrash); // crash
1502                } else if (doSprint>0) {
1503                  // sprint overrides idiot
1504                  solveOptions.setSpecialOption(1,3,doSprint); // sprint
1505                } else if (doIdiot>0) {
1506                  solveOptions.setSpecialOption(1,2,doIdiot); // idiot
1507                } else if (slpValue<=0) {
1508                  if (doIdiot==0) {
1509                    if (doSprint==0)
1510                      solveOptions.setSpecialOption(1,4); // all slack
1511                    else
1512                      solveOptions.setSpecialOption(1,9); // all slack or sprint
1513                  } else {
1514                    if (doSprint==0)
1515                      solveOptions.setSpecialOption(1,8); // all slack or idiot
1516                    else
1517                      solveOptions.setSpecialOption(1,7); // initiative
1518                  }
1519                }
1520                if (basisHasValues==-1)
1521                  solveOptions.setSpecialOption(1,11); // switch off values
1522              } else if (method==ClpSolve::useBarrier||method==ClpSolve::useBarrierNoCross) {
1523                int barrierOptions = choleskyType;
1524                if (scaleBarrier)
1525                  barrierOptions |= 8;
1526                if (doKKT)
1527                  barrierOptions |= 16;
1528                if (gamma)
1529                  barrierOptions |= 32*gamma;
1530                if (crossover==3) 
1531                  barrierOptions |= 256; // try presolve in crossover
1532                solveOptions.setSpecialOption(4,barrierOptions);
1533              }
1534              model2->initialSolve(solveOptions);
1535              basisHasValues=1;
1536              if (dualize) {
1537                int returnCode=((ClpSimplexOther *) lpSolver)->restoreFromDual(model2);
1538                delete model2;
1539                if (returnCode)
1540                  lpSolver->primal(1);
1541                model2=lpSolver;
1542              }
1543#ifdef COIN_HAS_ASL
1544              if (usingAmpl) {
1545                double value = model2->getObjValue()*model2->getObjSense();
1546                char buf[300];
1547                int pos=0;
1548                int iStat = model2->status();
1549                if (iStat==0) {
1550                  pos += sprintf(buf+pos,"optimal," );
1551                } else if (iStat==1) {
1552                  // infeasible
1553                  pos += sprintf(buf+pos,"infeasible,");
1554                } else if (iStat==2) {
1555                  // unbounded
1556                  pos += sprintf(buf+pos,"unbounded,");
1557                } else if (iStat==3) {
1558                  pos += sprintf(buf+pos,"stopped on iterations or time,");
1559                } else if (iStat==4) {
1560                  iStat = 7;
1561                  pos += sprintf(buf+pos,"stopped on difficulties,");
1562                } else if (iStat==5) {
1563                  iStat = 3;
1564                  pos += sprintf(buf+pos,"stopped on ctrl-c,");
1565                } else {
1566                  pos += sprintf(buf+pos,"status unknown,");
1567                  iStat=6;
1568                }
1569                info.problemStatus=iStat;
1570                info.objValue = value;
1571                pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
1572                               value);
1573                sprintf(buf+pos,"\n%d iterations",
1574                        model2->getIterationCount());
1575                free(info.primalSolution);
1576                int numberColumns=model2->numberColumns();
1577                info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
1578                CoinCopyN(model2->primalColumnSolution(),numberColumns,info.primalSolution);
1579                int numberRows = model2->numberRows();
1580                free(info.dualSolution);
1581                info.dualSolution = (double *) malloc(numberRows*sizeof(double));
1582                CoinCopyN(model2->dualRowSolution(),numberRows,info.dualSolution);
1583                CoinWarmStartBasis * basis = model2->getBasis();
1584                free(info.rowStatus);
1585                info.rowStatus = (int *) malloc(numberRows*sizeof(int));
1586                free(info.columnStatus);
1587                info.columnStatus = (int *) malloc(numberColumns*sizeof(int));
1588                // Put basis in
1589                int i;
1590                // free,basic,ub,lb are 0,1,2,3
1591                for (i=0;i<numberRows;i++) {
1592                  CoinWarmStartBasis::Status status = basis->getArtifStatus(i);
1593                  info.rowStatus[i]=status;
1594                }
1595                for (i=0;i<numberColumns;i++) {
1596                  CoinWarmStartBasis::Status status = basis->getStructStatus(i);
1597                  info.columnStatus[i]=status;
1598                }
1599                // put buffer into info
1600                strcpy(info.buffer,buf);
1601                delete basis;
1602              }
1603#endif
1604            } else {
1605              std::cout<<"** Current model not valid"<<std::endl;
1606            }
1607            break;
1608          case STATISTICS:
1609            if (goodModel) {
1610              // If presolve on look at presolved
1611              bool deleteModel2=false;
1612              ClpSimplex * model2 = lpSolver;
1613              if (preSolve) {
1614                ClpPresolve pinfo;
1615                int presolveOptions2 = presolveOptions&~0x40000000;
1616                if ((presolveOptions2&0xffff)!=0)
1617                  pinfo.setPresolveActions(presolveOptions2);
1618                pinfo.setSubstitution(substitution);
1619                if ((printOptions&1)!=0)
1620                  pinfo.statistics();
1621                double presolveTolerance = 
1622                  parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
1623                model2 = 
1624                  pinfo.presolvedModel(*lpSolver,presolveTolerance,
1625                                       true,preSolve);
1626                if (model2) {
1627                  printf("Statistics for presolved model\n");
1628                  deleteModel2=true;
1629                } else {
1630                  printf("Presolved model looks infeasible - will use unpresolved\n");
1631                  model2 = lpSolver;
1632                }
1633              } else {
1634                printf("Statistics for unpresolved model\n");
1635                model2 =  lpSolver;
1636              }
1637              statistics(lpSolver,model2);
1638              if (deleteModel2)
1639                delete model2;
1640            } else {
1641              std::cout<<"** Current model not valid"<<std::endl;
1642            }
1643            break;
1644          case TIGHTEN:
1645            if (goodModel) {
1646              int numberInfeasibilities = lpSolver->tightenPrimalBounds();
1647              if (numberInfeasibilities)
1648                std::cout<<"** Analysis indicates model infeasible"<<std::endl;
1649            } else {
1650              std::cout<<"** Current model not valid"<<std::endl;
1651            }
1652            break;
1653          case PLUSMINUS:
1654            if (goodModel) {
1655              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
1656              ClpPackedMatrix* clpMatrix =
1657                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
1658              if (clpMatrix) {
1659                ClpPlusMinusOneMatrix * newMatrix = new ClpPlusMinusOneMatrix(*(clpMatrix->matrix()));
1660                if (newMatrix->getIndices()) {
1661                  lpSolver->replaceMatrix(newMatrix);
1662                  delete saveMatrix;
1663                  std::cout<<"Matrix converted to +- one matrix"<<std::endl;
1664                } else {
1665                  std::cout<<"Matrix can not be converted to +- 1 matrix"<<std::endl;
1666                }
1667              } else {
1668                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
1669              }
1670            } else {
1671              std::cout<<"** Current model not valid"<<std::endl;
1672            }
1673            break;
1674          case OUTDUPROWS:
1675            if (goodModel) {
1676              int numberRows = clpSolver->getNumRows();
1677              //int nOut = outDupRow(clpSolver);
1678              CglDuplicateRow dupcuts(clpSolver);
1679              storedCuts = dupcuts.outDuplicates(clpSolver)!=0;
1680              int nOut = numberRows-clpSolver->getNumRows();
1681              if (nOut&&!noPrinting)
1682                printf("%d rows eliminated\n",nOut);
1683            } else {
1684              std::cout<<"** Current model not valid"<<std::endl;
1685            }
1686            break;
1687          case NETWORK:
1688            if (goodModel) {
1689              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
1690              ClpPackedMatrix* clpMatrix =
1691                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
1692              if (clpMatrix) {
1693                ClpNetworkMatrix * newMatrix = new ClpNetworkMatrix(*(clpMatrix->matrix()));
1694                if (newMatrix->getIndices()) {
1695                  lpSolver->replaceMatrix(newMatrix);
1696                  delete saveMatrix;
1697                  std::cout<<"Matrix converted to network matrix"<<std::endl;
1698                } else {
1699                  std::cout<<"Matrix can not be converted to network matrix"<<std::endl;
1700                }
1701              } else {
1702                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
1703              }
1704            } else {
1705              std::cout<<"** Current model not valid"<<std::endl;
1706            }
1707            break;
1708          case MIPLIB:
1709            // User can set options - main difference is lack of model and CglPreProcess
1710            goodModel=true;
1711/*
1712  Run branch-and-cut. First set a few options -- node comparison, scaling. If
1713  the solver is Clp, consider running some presolve code (not yet converted
1714  this to generic OSI) with branch-and-cut. If presolve is disabled, or the
1715  solver is not Clp, simply run branch-and-cut. Print elapsed time at the end.
1716*/
1717          case BAB: // branchAndBound
1718          case STRENGTHEN:
1719            if (goodModel) {
1720              bool miplib = type==MIPLIB;
1721              int logLevel = parameters[slog].intValue();
1722              // Reduce printout
1723              if (logLevel<=1)
1724                model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
1725              // Don't switch off all output
1726              {
1727                OsiSolverInterface * solver = model.solver();
1728                OsiClpSolverInterface * si =
1729                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
1730                assert (si != NULL);
1731                si->setSpecialOptions(0x40000000);
1732              }
1733              if (!miplib) {
1734                model.initialSolve();
1735                OsiSolverInterface * solver = model.solver();
1736                OsiClpSolverInterface * si =
1737                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
1738                ClpSimplex * clpSolver = si->getModelPtr();
1739                if (clpSolver->tightenPrimalBounds()!=0) {
1740                  std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
1741                  exit(1);
1742                }
1743                if (clpSolver->dualBound()==1.0e10) {
1744                  // user did not set - so modify
1745                  // get largest scaled away from bound
1746                  double largest=1.0e-12;
1747                  int numberRows = clpSolver->numberRows();
1748                  const double * rowPrimal = clpSolver->primalRowSolution();
1749                  const double * rowLower = clpSolver->rowLower();
1750                  const double * rowUpper = clpSolver->rowUpper();
1751                  const double * rowScale = clpSolver->rowScale();
1752                  int iRow;
1753                  for (iRow=0;iRow<numberRows;iRow++) {
1754                    double value = rowPrimal[iRow];
1755                    double above = value-rowLower[iRow];
1756                    double below = rowUpper[iRow]-value;
1757                    if (rowScale) {
1758                      double multiplier = rowScale[iRow];
1759                      above *= multiplier;
1760                      below *= multiplier;
1761                    }
1762                    if (above<1.0e12)
1763                      largest = CoinMax(largest,above);
1764                    if (below<1.0e12)
1765                      largest = CoinMax(largest,below);
1766                  }
1767                 
1768                  int numberColumns = clpSolver->numberColumns();
1769                  const double * columnPrimal = clpSolver->primalColumnSolution();
1770                  const double * columnLower = clpSolver->columnLower();
1771                  const double * columnUpper = clpSolver->columnUpper();
1772                  const double * columnScale = clpSolver->columnScale();
1773                  int iColumn;
1774                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
1775                    double value = columnPrimal[iColumn];
1776                    double above = value-columnLower[iColumn];
1777                    double below = columnUpper[iColumn]-value;
1778                    if (columnScale) {
1779                      double multiplier = 1.0/columnScale[iColumn];
1780                      above *= multiplier;
1781                      below *= multiplier;
1782                    }
1783                    if (above<1.0e12)
1784                      largest = CoinMax(largest,above);
1785                    if (below<1.0e12)
1786                      largest = CoinMax(largest,below);
1787                  }
1788                  //if (!noPrinting)
1789                  //std::cout<<"Largest (scaled) away from bound "<<largest<<std::endl;
1790                  clpSolver->setDualBound(CoinMax(1.0001e8,CoinMin(1000.0*largest,1.00001e10)));
1791                }
1792                clpSolver->dual();  // clean up
1793              }
1794              // If user made settings then use them
1795              if (!defaultSettings) {
1796                OsiSolverInterface * solver = model.solver();
1797                if (!doScaling)
1798                  solver->setHintParam(OsiDoScale,false,OsiHintTry);
1799                OsiClpSolverInterface * si =
1800                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
1801                assert (si != NULL);
1802                // get clp itself
1803                ClpSimplex * modelC = si->getModelPtr();
1804                //if (modelC->tightenPrimalBounds()!=0) {
1805                //std::cout<<"Problem is infeasible!"<<std::endl;
1806                //break;
1807                //}
1808                // bounds based on continuous
1809                if (tightenFactor) {
1810                  if (modelC->tightenPrimalBounds(tightenFactor)!=0) {
1811                    std::cout<<"Problem is infeasible!"<<std::endl;
1812                    break;
1813                  }
1814                }
1815                if (djFix<1.0e20) {
1816                  // do some fixing
1817                  int numberColumns = modelC->numberColumns();
1818                  int i;
1819                  const char * type = modelC->integerInformation();
1820                  double * lower = modelC->columnLower();
1821                  double * upper = modelC->columnUpper();
1822                  double * solution = modelC->primalColumnSolution();
1823                  double * dj = modelC->dualColumnSolution();
1824                  int numberFixed=0;
1825                  for (i=0;i<numberColumns;i++) {
1826                    if (type[i]) {
1827                      double value = solution[i];
1828                      if (value<lower[i]+1.0e-5&&dj[i]>djFix) {
1829                        solution[i]=lower[i];
1830                        upper[i]=lower[i];
1831                        numberFixed++;
1832                      } else if (value>upper[i]-1.0e-5&&dj[i]<-djFix) {
1833                        solution[i]=upper[i];
1834                        lower[i]=upper[i];
1835                        numberFixed++;
1836                      }
1837                    }
1838                  }
1839                  printf("%d columns fixed\n",numberFixed);
1840                }
1841              }
1842              // See if we want preprocessing
1843              OsiSolverInterface * saveSolver=NULL;
1844              CglPreProcess process;
1845              delete babModel;
1846              babModel = new CbcModel(model);
1847              OsiSolverInterface * solver3 = clpSolver->clone();
1848              babModel->assignSolver(solver3);
1849              OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
1850              int numberChanged=0;
1851              if (clpSolver2->messageHandler()->logLevel())
1852                clpSolver2->messageHandler()->setLogLevel(1);
1853              if (logLevel>-1)
1854                clpSolver2->messageHandler()->setLogLevel(logLevel);
1855              lpSolver = clpSolver2->getModelPtr();
1856              if (lpSolver->factorizationFrequency()==200&&!miplib) {
1857                // User did not touch preset
1858                int numberRows = lpSolver->numberRows();
1859                const int cutoff1=10000;
1860                const int cutoff2=100000;
1861                const int base=75;
1862                const int freq0 = 50;
1863                const int freq1=200;
1864                const int freq2=400;
1865                const int maximum=1000;
1866                int frequency;
1867                if (numberRows<cutoff1)
1868                  frequency=base+numberRows/freq0;
1869                else if (numberRows<cutoff2)
1870                  frequency=base+cutoff1/freq0 + (numberRows-cutoff1)/freq1;
1871                else
1872                  frequency=base+cutoff1/freq0 + (cutoff2-cutoff1)/freq1 + (numberRows-cutoff2)/freq2;
1873                lpSolver->setFactorizationFrequency(CoinMin(maximum,frequency));
1874              }
1875              time2 = CoinCpuTime();
1876              totalTime += time2-time1;
1877              time1 = time2;
1878              double timeLeft = babModel->getMaximumSeconds();
1879              int numberOriginalColumns = babModel->solver()->getNumCols();
1880              if (preProcess==7) {
1881                // use strategy instead
1882                preProcess=0;
1883                useStrategy=true;
1884              }
1885              if (preProcess&&type==BAB) {
1886                // See if sos from mps file
1887                if (numberSOS==0&&clpSolver->numberSOS()&&doSOS) {
1888                  // SOS
1889                  numberSOS = clpSolver->numberSOS();
1890                  const CoinSet * setInfo = clpSolver->setInfo();
1891                  sosStart = new int [numberSOS+1];
1892                  sosType = new char [numberSOS];
1893                  int i;
1894                  int nTotal=0;
1895                  sosStart[0]=0;
1896                  for ( i=0;i<numberSOS;i++) {
1897                    int type = setInfo[i].setType();
1898                    int n=setInfo[i].numberEntries();
1899                    sosType[i]=type;
1900                    nTotal += n;
1901                    sosStart[i+1] = nTotal;
1902                  }
1903                  sosIndices = new int[nTotal];
1904                  sosReference = new double [nTotal];
1905                  for (i=0;i<numberSOS;i++) {
1906                    int n=setInfo[i].numberEntries();
1907                    const int * which = setInfo[i].which();
1908                    const double * weights = setInfo[i].weights();
1909                    int base = sosStart[i];
1910                    for (int j=0;j<n;j++) {
1911                      int k=which[j];
1912                      sosIndices[j+base]=k;
1913                      sosReference[j+base] = weights ? weights[j] : (double) j;
1914                    }
1915                  }
1916                }
1917                saveSolver=babModel->solver()->clone();
1918                /* Do not try and produce equality cliques and
1919                   do up to 10 passes */
1920                OsiSolverInterface * solver2;
1921                {
1922                  // Tell solver we are in Branch and Cut
1923                  saveSolver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo) ;
1924                  // Default set of cut generators
1925                  CglProbing generator1;
1926                  generator1.setUsingObjective(true);
1927                  generator1.setMaxPass(3);
1928                  generator1.setMaxProbeRoot(saveSolver->getNumCols());
1929                  generator1.setMaxElements(100);
1930                  generator1.setMaxLookRoot(50);
1931                  generator1.setRowCuts(3);
1932                  // Add in generators
1933                  process.addCutGenerator(&generator1);
1934                  int translate[]={9999,0,0,-1,2,3,-2};
1935                  process.messageHandler()->setLogLevel(babModel->logLevel());
1936#ifdef COIN_HAS_ASL
1937                  if (info.numberSos&&doSOS&&usingAmpl) {
1938                    // SOS
1939                    numberSOS = info.numberSos;
1940                    sosStart = info.sosStart;
1941                    sosIndices = info.sosIndices;
1942                  }
1943#endif
1944                  if (numberSOS&&doSOS) {
1945                    // SOS
1946                    int numberColumns = saveSolver->getNumCols();
1947                    char * prohibited = new char[numberColumns];
1948                    memset(prohibited,0,numberColumns);
1949                    int n=sosStart[numberSOS];
1950                    for (int i=0;i<n;i++) {
1951                      int iColumn = sosIndices[i];
1952                      prohibited[iColumn]=1;
1953                    }
1954                    process.passInProhibited(prohibited,numberColumns);
1955                    delete [] prohibited;
1956                  }
1957                  solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],10,
1958                                                         tunePreProcess);
1959                  // Tell solver we are not in Branch and Cut
1960                  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
1961                  if (solver2)
1962                    solver2->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
1963                }
1964#ifdef COIN_HAS_ASL
1965                if (!solver2&&usingAmpl) {
1966                  // infeasible
1967                  info.problemStatus=1;
1968                  info.objValue = 1.0e100;
1969                  sprintf(info.buffer,"infeasible/unbounded by pre-processing");
1970                  info.primalSolution=NULL;
1971                  info.dualSolution=NULL;
1972                  break;
1973                }
1974#endif
1975                if (!noPrinting) {
1976                  if (!solver2) {
1977                    printf("Pre-processing says infeasible or unbounded\n");
1978                    break;
1979                  } else {
1980                    printf("processed model has %d rows, %d columns and %d elements\n",
1981                           solver2->getNumRows(),solver2->getNumCols(),solver2->getNumElements());
1982                  }
1983                }
1984                //solver2->resolve();
1985                if (preProcess==2) {
1986                  OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (solver2);
1987                  ClpSimplex * lpSolver = clpSolver2->getModelPtr();
1988                  lpSolver->writeMps("presolved.mps",0,1,lpSolver->optimizationDirection());
1989                  printf("Preprocessed model (minimization) on presolved.mps\n");
1990                }
1991                // we have to keep solver2 so pass clone
1992                solver2 = solver2->clone();
1993                babModel->assignSolver(solver2);
1994                babModel->initialSolve();
1995                babModel->setMaximumSeconds(timeLeft-(CoinCpuTime()-time1));
1996              }
1997              // now tighten bounds
1998              if (!miplib) {
1999                OsiClpSolverInterface * si =
2000                  dynamic_cast<OsiClpSolverInterface *>(babModel->solver()) ;
2001                assert (si != NULL);
2002                // get clp itself
2003                ClpSimplex * modelC = si->getModelPtr();
2004                if (noPrinting)
2005                  modelC->setLogLevel(0);
2006                if (modelC->tightenPrimalBounds()!=0) {
2007                  std::cout<<"Problem is infeasible!"<<std::endl;
2008                  break;
2009                }
2010                modelC->dual();
2011              }
2012              if (debugValues) {
2013                // for debug
2014                std::string problemName ;
2015                babModel->solver()->getStrParam(OsiProbName,problemName) ;
2016                babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
2017                twomirGen.probname_=strdup(problemName.c_str());
2018                // checking seems odd
2019                //redsplitGen.set_given_optsol(babModel->solver()->getRowCutDebuggerAlways()->optimalSolution(),
2020                //                         babModel->getNumCols());
2021              }
2022              if (useCosts) {
2023                int numberColumns = babModel->getNumCols();
2024                int * sort = new int[numberColumns];
2025                double * dsort = new double[numberColumns];
2026                int * priority = new int [numberColumns];
2027                const double * objective = babModel->getObjCoefficients();
2028                const double * lower = babModel->getColLower() ;
2029                const double * upper = babModel->getColUpper() ;
2030                const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
2031                const int * columnLength = matrix->getVectorLengths();
2032                int iColumn;
2033                int n=0;
2034                for (iColumn=0;iColumn<numberColumns;iColumn++) {
2035                  if (babModel->isInteger(iColumn)) {
2036                    sort[n]=n;
2037                    if (useCosts==1)
2038                      dsort[n++]=-fabs(objective[iColumn]);
2039                    else if (useCosts==2)
2040                      dsort[n++]=iColumn;
2041                    else if (useCosts==3)
2042                      dsort[n++]=upper[iColumn]-lower[iColumn];
2043                    else if (useCosts==4)
2044                      dsort[n++]=-(upper[iColumn]-lower[iColumn]);
2045                    else if (useCosts==5)
2046                      dsort[n++]=-columnLength[iColumn];
2047                  }
2048                }
2049                CoinSort_2(dsort,dsort+n,sort);
2050                int level=0;
2051                double last = -1.0e100;
2052                for (int i=0;i<n;i++) {
2053                  int iPut=sort[i];
2054                  if (dsort[i]!=last) {
2055                    level++;
2056                    last=dsort[i];
2057                  }
2058                  priority[iPut]=level;
2059                }
2060                babModel->passInPriorities( priority,false);
2061                delete [] priority;
2062                delete [] sort;
2063                delete [] dsort;
2064              }
2065              // FPump done first as it only works if no solution
2066              CbcHeuristicFPump heuristic4(*babModel);
2067              if (useFpump) {
2068                heuristic4.setMaximumPasses(parameters[whichParam(FPUMPITS,numberParameters,parameters)].intValue());
2069                int pumpTune=parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].intValue();
2070                if (pumpTune>0) {
2071                  /*
2072                    >=1000000 use as accumulate switch
2073                    >=1000 use index+2 as number of large loops
2074                    >=100 use 0.05 objvalue as increment
2075                    >=10 use +0.1 objvalue for cutoff (add)
2076                    1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
2077                    4 and static continuous, 5 as 3 but no internal integers
2078                  */
2079                  double value = babModel->solver()->getObjSense()*babModel->solver()->getObjValue();
2080                  int w = pumpTune/10;
2081                  int c = w % 10;
2082                  w /= 10;
2083                  int i = w % 10;
2084                  w /= 10;
2085                  int r = w;
2086                  int accumulate = r/1000;
2087                  r -= 1000*accumulate;
2088                  // fake cutoff
2089                  printf("Setting ");
2090                  if (c) {
2091                    double cutoff;
2092                    babModel->solver()->getDblParam(OsiDualObjectiveLimit,cutoff);
2093                    cutoff = CoinMin(cutoff,value + 0.1*fabs(value)*c);
2094                    heuristic4.setFakeCutoff(cutoff);
2095                    printf("fake cutoff of %g ",cutoff);
2096                  }
2097                  if (i||r) {
2098                    // also set increment
2099                    heuristic4.setAbsoluteIncrement((0.01*i+0.005)*(fabs(value)+1.0e-12));
2100                    heuristic4.setAccumulate(accumulate);
2101                    heuristic4.setMaximumRetries(r+2);
2102                    if (i)
2103                      printf("increment of %g ",heuristic4.absoluteIncrement());
2104                    if (accumulate)
2105                      printf("accumulate of %d ",accumulate);
2106                    printf("%d retries ",r+2);
2107                  }
2108                  pumpTune = pumpTune%100;
2109                  printf("and setting when to %d\n",pumpTune+10);
2110                  heuristic4.setWhen(pumpTune+10);
2111                }
2112                babModel->addHeuristic(&heuristic4);
2113              }
2114              if (!miplib) {
2115                CbcRounding heuristic1(*babModel);
2116                if (useRounding)
2117                  babModel->addHeuristic(&heuristic1) ;
2118                CbcHeuristicLocal heuristic2(*babModel);
2119                heuristic2.setSearchType(1);
2120                if (useCombine)
2121                  babModel->addHeuristic(&heuristic2);
2122                CbcHeuristicGreedyCover heuristic3(*babModel);
2123                CbcHeuristicGreedyEquality heuristic3a(*babModel);
2124                if (useGreedy) {
2125                  babModel->addHeuristic(&heuristic3);
2126                  babModel->addHeuristic(&heuristic3a);
2127                }
2128                if (useLocalTree) {
2129                  CbcTreeLocal localTree(babModel,NULL,10,0,0,10000,2000);
2130                  babModel->passInTreeHandler(localTree);
2131                }
2132              }
2133              if (type==MIPLIB) {
2134                if (babModel->numberStrong()==5&&babModel->numberBeforeTrust()==5) 
2135                  babModel->setNumberBeforeTrust(50);
2136              }
2137              // add cut generators if wanted
2138              int switches[20];
2139              int numberGenerators=0;
2140              int translate[]={-100,-1,-99,-98,1,1,1,1};
2141              if (probingAction) {
2142                if (probingAction==5||probingAction==7)
2143                  probingGen.setRowCuts(-3); // strengthening etc just at root
2144                if (probingAction==6||probingAction==7) {
2145                  // Number of unsatisfied variables to look at
2146                  probingGen.setMaxProbe(1000);
2147                  probingGen.setMaxProbeRoot(1000);
2148                  // How far to follow the consequences
2149                  probingGen.setMaxLook(50);
2150                  probingGen.setMaxLookRoot(50);
2151                }
2152                babModel->addCutGenerator(&probingGen,translate[probingAction],"Probing");
2153                switches[numberGenerators++]=0;
2154              }
2155              if (gomoryAction) {
2156                babModel->addCutGenerator(&gomoryGen,translate[gomoryAction],"Gomory");
2157                switches[numberGenerators++]=-1;
2158              }
2159              if (knapsackAction) {
2160                babModel->addCutGenerator(&knapsackGen,translate[knapsackAction],"Knapsack");
2161                switches[numberGenerators++]=0;
2162              }
2163              if (redsplitAction) {
2164                babModel->addCutGenerator(&redsplitGen,translate[redsplitAction],"Reduce-and-split");
2165                switches[numberGenerators++]=1;
2166              }
2167              if (cliqueAction) {
2168                babModel->addCutGenerator(&cliqueGen,translate[cliqueAction],"Clique");
2169                switches[numberGenerators++]=0;
2170              }
2171              if (mixedAction) {
2172                babModel->addCutGenerator(&mixedGen,translate[mixedAction],"MixedIntegerRounding2");
2173                switches[numberGenerators++]=-1;
2174              }
2175              if (flowAction) {
2176                babModel->addCutGenerator(&flowGen,translate[flowAction],"FlowCover");
2177                switches[numberGenerators++]=1;
2178              }
2179              if (twomirAction) {
2180                babModel->addCutGenerator(&twomirGen,translate[twomirAction],"TwoMirCuts");
2181                switches[numberGenerators++]=1;
2182              }
2183              if (landpAction) {
2184                babModel->addCutGenerator(&landpGen,translate[landpAction],"LiftAndProject");
2185                switches[numberGenerators++]=1;
2186              }
2187              if (storedCuts) 
2188                babModel->setSpecialOptions(babModel->specialOptions()|64);
2189              // Say we want timings
2190              numberGenerators = babModel->numberCutGenerators();
2191              int iGenerator;
2192              int cutDepth=
2193                parameters[whichParam(CUTDEPTH,numberParameters,parameters)].intValue();
2194              for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
2195                CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
2196                int howOften = generator->howOften();
2197                if (howOften==-98||howOften==-99) 
2198                  generator->setSwitchOffIfLessThan(switches[iGenerator]);
2199                generator->setTiming(true);
2200                if (cutDepth>=0)
2201                  generator->setWhatDepth(cutDepth) ;
2202              }
2203              // Could tune more
2204              if (!miplib) {
2205                babModel->setMinimumDrop(min(5.0e-2,
2206                                             fabs(babModel->getMinimizationObjValue())*1.0e-3+1.0e-4));
2207                if (cutPass==-1234567) {
2208                  if (babModel->getNumCols()<500)
2209                    babModel->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
2210                  else if (babModel->getNumCols()<5000)
2211                    babModel->setMaximumCutPassesAtRoot(100); // use minimum drop
2212                  else
2213                    babModel->setMaximumCutPassesAtRoot(20);
2214                } else {
2215                  babModel->setMaximumCutPassesAtRoot(cutPass);
2216                }
2217                babModel->setMaximumCutPasses(1);
2218              }
2219              // Do more strong branching if small
2220              //if (babModel->getNumCols()<5000)
2221              //babModel->setNumberStrong(20);
2222              // Switch off strong branching if wanted
2223              //if (babModel->getNumCols()>10*babModel->getNumRows())
2224              //babModel->setNumberStrong(0);
2225              if (!noPrinting) {
2226                int iLevel = parameters[log].intValue();
2227                if (iLevel<0) {
2228                  babModel->setPrintingMode(1);
2229                  iLevel = -iLevel;
2230                }
2231                babModel->messageHandler()->setLogLevel(iLevel);
2232                if (babModel->getNumCols()>2000||babModel->getNumRows()>1500||
2233                    babModel->messageHandler()->logLevel()>1)
2234                  babModel->setPrintFrequency(100);
2235              }
2236             
2237              babModel->solver()->setIntParam(OsiMaxNumIterationHotStart,
2238                    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].intValue());
2239              OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2240              // go faster stripes
2241              if (osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
2242                osiclp->setupForRepeatedUse(2,parameters[slog].intValue());
2243              } else {
2244                osiclp->setupForRepeatedUse(0,parameters[slog].intValue());
2245              }
2246              double increment=babModel->getCutoffIncrement();;
2247              int * changed = NULL;
2248              if (!miplib)
2249                changed=analyze( osiclp,numberChanged,increment,false);
2250              if (debugValues) {
2251                if (numberDebugValues==babModel->getNumCols()) {
2252                  // for debug
2253                  babModel->solver()->activateRowCutDebugger(debugValues) ;
2254                } else {
2255                  printf("debug file has incorrect number of columns\n");
2256                }
2257              }
2258              babModel->setCutoffIncrement(CoinMax(babModel->getCutoffIncrement(),increment));
2259              // Turn this off if you get problems
2260              // Used to be automatically set
2261              int mipOptions = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue();
2262              if (mipOptions!=(128|64|1))
2263                printf("mip options %d\n",mipOptions);
2264              osiclp->setSpecialOptions(mipOptions);
2265              if (gapRatio < 1.0e100) {
2266                double value = babModel->solver()->getObjValue() ;
2267                double value2 = gapRatio*(1.0e-5+fabs(value)) ;
2268                babModel->setAllowableGap(value2) ;
2269                std::cout << "Continuous " << value
2270                          << ", so allowable gap set to "
2271                          << value2 << std::endl ;
2272              }
2273              // probably faster to use a basis to get integer solutions
2274              babModel->setSpecialOptions(babModel->specialOptions()|2);
2275              currentBranchModel = babModel;
2276              OsiSolverInterface * strengthenedModel=NULL;
2277              if (type==BAB||type==MIPLIB) {
2278                int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
2279                if (moreMipOptions>=0) {
2280                  printf("more mip options %d\n",moreMipOptions);
2281                  if (((moreMipOptions+1)%1000000)!=0)
2282                    babModel->setSearchStrategy(moreMipOptions%1000000);
2283                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2284                  // go faster stripes
2285                  if( moreMipOptions >=999999) {
2286                    if (osiclp) {
2287                      int save = osiclp->specialOptions();
2288                      osiclp->setupForRepeatedUse(2,0);
2289                      osiclp->setSpecialOptions(save|osiclp->specialOptions());
2290                    }
2291                  } 
2292                }
2293              }
2294              if (type==BAB) {
2295#ifdef COIN_HAS_ASL
2296                if (usingAmpl) {
2297                  priorities=info.priorities;
2298                  branchDirection=info.branchDirection;
2299                  pseudoDown=info.pseudoDown;
2300                  pseudoUp=info.pseudoUp;
2301                  solutionIn=info.primalSolution;
2302                  prioritiesIn = info.priorities;
2303                  if (info.numberSos&&doSOS) {
2304                    // SOS
2305                    numberSOS = info.numberSos;
2306                    sosStart = info.sosStart;
2307                    sosIndices = info.sosIndices;
2308                    sosType = info.sosType;
2309                    sosReference = info.sosReference;
2310                    sosPriority = info.sosPriority;
2311                  }
2312                }
2313#endif               
2314                const int * originalColumns = preProcess ? process.originalColumns() : NULL;
2315                if (solutionIn&&useSolution) {
2316                  if (preProcess) {
2317                    int numberColumns = babModel->getNumCols();
2318                    // extend arrays in case SOS
2319                    int n = originalColumns[numberColumns-1]+1;
2320                    int nSmaller = CoinMin(n,numberOriginalColumns);
2321                    double * solutionIn2 = new double [n];
2322                    int * prioritiesIn2 = new int[n];
2323                    int i;
2324                    for (i=0;i<nSmaller;i++) {
2325                      solutionIn2[i]=solutionIn[i];
2326                      prioritiesIn2[i]=prioritiesIn[i];
2327                    }
2328                    for (;i<n;i++) {
2329                      solutionIn2[i]=0.0;
2330                      prioritiesIn2[i]=1000000;
2331                    }
2332                    int iLast=-1;
2333                    for (i=0;i<numberColumns;i++) {
2334                      int iColumn = originalColumns[i];
2335                      assert (iColumn>iLast);
2336                      iLast=iColumn;
2337                      solutionIn2[i]=solutionIn2[iColumn];
2338                      if (prioritiesIn)
2339                        prioritiesIn2[i]=prioritiesIn2[iColumn];
2340                    }
2341                    babModel->setHotstartSolution(solutionIn2,prioritiesIn2);
2342                    delete [] solutionIn2;
2343                    delete [] prioritiesIn2;
2344                  } else {
2345                    babModel->setHotstartSolution(solutionIn,prioritiesIn);
2346                  }
2347                }
2348                if (preProcess&&process.numberSOS()) {
2349                  int numberSOS = process.numberSOS();
2350                  int numberIntegers = babModel->numberIntegers();
2351                  /* model may not have created objects
2352                     If none then create
2353                  */
2354                  if (!numberIntegers||!babModel->numberObjects()) {
2355                    int type = (pseudoUp) ? 1 : 0;
2356                    babModel->findIntegers(true,type);
2357                    numberIntegers = babModel->numberIntegers();
2358                  }
2359                  OsiObject ** oldObjects = babModel->objects();
2360                  // Do sets and priorities
2361                  OsiObject ** objects = new OsiObject * [numberSOS];
2362                  // set old objects to have low priority
2363                  int numberOldObjects = babModel->numberObjects();
2364                  int numberColumns = babModel->getNumCols();
2365                  for (int iObj = 0;iObj<numberOldObjects;iObj++) {
2366                    oldObjects[iObj]->setPriority(numberColumns+1);
2367                    int iColumn = oldObjects[iObj]->columnNumber();
2368                    assert (iColumn>=0);
2369                    if (iColumn>=numberOriginalColumns)
2370                      continue;
2371                    if (originalColumns)
2372                      iColumn = originalColumns[iColumn];
2373                    if (branchDirection) {
2374                      CbcSimpleInteger * obj =
2375                        dynamic_cast <CbcSimpleInteger *>(oldObjects[iObj]) ;
2376                      if (obj) { 
2377                        obj->setPreferredWay(branchDirection[iColumn]);
2378                      } else {
2379                        CbcObject * obj =
2380                          dynamic_cast <CbcObject *>(oldObjects[iObj]) ;
2381                        assert (obj);
2382                        obj->setPreferredWay(branchDirection[iColumn]);
2383                      }
2384                    }
2385                    if (pseudoUp) {
2386                      CbcSimpleIntegerPseudoCost * obj1a =
2387                        dynamic_cast <CbcSimpleIntegerPseudoCost *>(oldObjects[iObj]) ;
2388                      assert (obj1a);
2389                      if (pseudoDown[iColumn]>0.0)
2390                        obj1a->setDownPseudoCost(pseudoDown[iColumn]);
2391                      if (pseudoUp[iColumn]>0.0)
2392                        obj1a->setUpPseudoCost(pseudoUp[iColumn]);
2393                    }
2394                  }
2395                  const int * starts = process.startSOS();
2396                  const int * which = process.whichSOS();
2397                  const int * type = process.typeSOS();
2398                  const double * weight = process.weightSOS();
2399                  int iSOS;
2400                  for (iSOS =0;iSOS<numberSOS;iSOS++) {
2401                    int iStart = starts[iSOS];
2402                    int n=starts[iSOS+1]-iStart;
2403                    objects[iSOS] = new CbcSOS(babModel,n,which+iStart,weight+iStart,
2404                                               iSOS,type[iSOS]);
2405                    // branch on long sets first
2406                    objects[iSOS]->setPriority(numberColumns-n);
2407                  }
2408                  babModel->addObjects(numberSOS,objects);
2409                  for (iSOS=0;iSOS<numberSOS;iSOS++)
2410                    delete objects[iSOS];
2411                  delete [] objects;
2412                } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
2413                  // do anyway for priorities etc
2414                  int numberIntegers = babModel->numberIntegers();
2415                  /* model may not have created objects
2416                     If none then create
2417                  */
2418                  if (!numberIntegers||!babModel->numberObjects()) {
2419                    int type = (pseudoUp) ? 1 : 0;
2420                    babModel->findIntegers(true,type);
2421                  }
2422                  if (numberSOS) {
2423                    // Do sets and priorities
2424                    OsiObject ** objects = new OsiObject * [numberSOS];
2425                    int iSOS;
2426                    if (originalColumns) {
2427                      // redo sequence numbers
2428                      int numberColumns = babModel->getNumCols();
2429                      int nOld = originalColumns[numberColumns-1]+1;
2430                      int * back = new int[nOld];
2431                      int i;
2432                      for (i=0;i<nOld;i++)
2433                        back[i]=-1;
2434                      for (i=0;i<numberColumns;i++)
2435                        back[originalColumns[i]]=i;
2436                      // Really need better checks
2437                      int nMissing=0;
2438                      int n=sosStart[numberSOS];
2439                      for (i=0;i<n;i++) {
2440                        int iColumn = sosIndices[i];
2441                        int jColumn = back[iColumn];
2442                        if (jColumn>=0) 
2443                          sosIndices[i] = jColumn;
2444                        else 
2445                          nMissing++;
2446                      }
2447                      delete [] back;
2448                      if (nMissing)
2449                        printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
2450                    }
2451                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
2452                      int iStart = sosStart[iSOS];
2453                      int n=sosStart[iSOS+1]-iStart;
2454                      objects[iSOS] = new CbcSOS(babModel,n,sosIndices+iStart,sosReference+iStart,
2455                                                 iSOS,sosType[iSOS]);
2456                      if (sosPriority)
2457                        objects[iSOS]->setPriority(sosPriority[iSOS]);
2458                      else if (!prioritiesIn)
2459                        objects[iSOS]->setPriority(10);  // rather than 1000
2460                    }
2461                    babModel->addObjects(numberSOS,objects);
2462                    for (iSOS=0;iSOS<numberSOS;iSOS++)
2463                      delete objects[iSOS];
2464                    delete [] objects;
2465                  }
2466                  OsiObject ** objects = babModel->objects();
2467                  int numberObjects = babModel->numberObjects();
2468                  for (int iObj = 0;iObj<numberObjects;iObj++) {
2469                    // skip sos
2470                    CbcSOS * objSOS =
2471                      dynamic_cast <CbcSOS *>(objects[iObj]) ;
2472                    if (objSOS)
2473                      continue;
2474                    int iColumn = objects[iObj]->columnNumber();
2475                    assert (iColumn>=0);
2476                    if (originalColumns)
2477                      iColumn = originalColumns[iColumn];
2478                    if (branchDirection) {
2479                      CbcSimpleInteger * obj =
2480                        dynamic_cast <CbcSimpleInteger *>(objects[iObj]) ;
2481                      if (obj) { 
2482                        obj->setPreferredWay(branchDirection[iColumn]);
2483                      } else {
2484                        CbcObject * obj =
2485                          dynamic_cast <CbcObject *>(objects[iObj]) ;
2486                        assert (obj);
2487                        obj->setPreferredWay(branchDirection[iColumn]);
2488                      }
2489                    }
2490                    if (priorities) {
2491                      int iPriority = priorities[iColumn];
2492                      if (iPriority>0)
2493                        objects[iObj]->setPriority(iPriority);
2494                    }
2495                    if (pseudoUp&&pseudoUp[iColumn]) {
2496                      CbcSimpleIntegerPseudoCost * obj1a =
2497                        dynamic_cast <CbcSimpleIntegerPseudoCost *>(objects[iObj]) ;
2498                      assert (obj1a);
2499                      if (pseudoDown[iColumn]>0.0)
2500                        obj1a->setDownPseudoCost(pseudoDown[iColumn]);
2501                      if (pseudoUp[iColumn]>0.0)
2502                        obj1a->setUpPseudoCost(pseudoUp[iColumn]);
2503                    }
2504                  }
2505                }
2506                int statistics = (printOptions>0) ? printOptions: 0;
2507#ifdef COIN_HAS_ASL
2508                if (!usingAmpl) {
2509#endif
2510                  free(priorities);
2511                  priorities=NULL;
2512                  free(branchDirection);
2513                  branchDirection=NULL;
2514                  free(pseudoDown);
2515                  pseudoDown=NULL;
2516                  free(pseudoUp);
2517                  pseudoUp=NULL;
2518                  free(solutionIn);
2519                  solutionIn=NULL;
2520                  free(prioritiesIn);
2521                  prioritiesIn=NULL;
2522                  free(sosStart);
2523                  sosStart=NULL;
2524                  free(sosIndices);
2525                  sosIndices=NULL;
2526                  free(sosType);
2527                  sosType=NULL;
2528                  free(sosReference);
2529                  sosReference=NULL;
2530                  free(sosPriority);
2531                  sosPriority=NULL;
2532#ifdef COIN_HAS_ASL
2533                }
2534#endif               
2535                if (nodeStrategy) {
2536                  // change default
2537                  if (nodeStrategy>1) {
2538                    // up or down
2539                    int way = ((nodeStrategy%1)==1) ? -1 : +1;
2540                    babModel->setPreferredWay(way);
2541#if 0
2542                    OsiObject ** objects = babModel->objects();
2543                    int numberObjects = babModel->numberObjects();
2544                    for (int iObj = 0;iObj<numberObjects;iObj++) {
2545                      CbcObject * obj =
2546                        dynamic_cast <CbcObject *>(objects[iObj]) ;
2547                      assert (obj);
2548                      obj->setPreferredWay(way);
2549                    }
2550#endif
2551                  }
2552                  if (nodeStrategy==1||nodeStrategy>3) {
2553                    // depth
2554                    CbcCompareDefault compare;
2555                    compare.setWeight(-3.0);
2556                    babModel->setNodeComparison(compare);
2557                  }
2558                }
2559                if (cppValue>=0) {
2560                  int prepro = useStrategy ? -1 : preProcess;
2561                  // generate code
2562                  FILE * fp = fopen("user_driver.cpp","w");
2563                  if (fp) {
2564                    // generate enough to do BAB
2565                    babModel->generateCpp(fp,1);
2566                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2567                    // Make general so do factorization
2568                    int factor = osiclp->getModelPtr()->factorizationFrequency();
2569                    osiclp->getModelPtr()->setFactorizationFrequency(200);
2570                    osiclp->generateCpp(fp);
2571                    osiclp->getModelPtr()->setFactorizationFrequency(factor);
2572                    //solveOptions.generateCpp(fp);
2573                    fclose(fp);
2574                    // now call generate code
2575                    generateCode(babModel,"user_driver.cpp",cppValue,prepro);
2576                  } else {
2577                    std::cout<<"Unable to open file user_driver.cpp"<<std::endl;
2578                  }
2579                }
2580                if (useStrategy) {
2581                  CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
2582                  strategy.setupPreProcessing(1);
2583                  babModel->setStrategy(strategy);
2584                }
2585                int testOsiOptions = parameters[whichParam(TESTOSI,numberParameters,parameters)].intValue();
2586                if (testOsiOptions>0) {
2587                  printf("Testing OsiObject options %d\n",testOsiOptions);
2588                  CbcBranchDefaultDecision decision;
2589                  OsiChooseVariable choose(babModel->solver());
2590                  decision.setChooseMethod(choose);
2591                  babModel->setBranchingMethod(decision);
2592                }
2593                checkSOS(babModel, babModel->solver());
2594                babModel->branchAndBound(statistics);
2595                checkSOS(babModel, babModel->solver());
2596              } else if (type==MIPLIB) {
2597                CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
2598                // Set up pre-processing
2599                int translate2[]={9999,1,1,3,2,4,5};
2600                if (preProcess)
2601                  strategy.setupPreProcessing(translate2[preProcess]);
2602                babModel->setStrategy(strategy);
2603                CbcClpUnitTest(*babModel);
2604                goodModel=false;
2605                break;
2606              } else {
2607                strengthenedModel = babModel->strengthenedModel();
2608              }
2609              currentBranchModel = NULL;
2610              osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2611              if (debugFile=="createAfterPre"&&babModel->bestSolution()) {
2612                lpSolver = osiclp->getModelPtr();
2613                //move best solution (should be there -- but ..)
2614                int n = lpSolver->getNumCols();
2615                memcpy(lpSolver->primalColumnSolution(),babModel->bestSolution(),n*sizeof(double));
2616                saveSolution(osiclp->getModelPtr(),"debug.file");
2617              }
2618              if (!noPrinting) {
2619                // Print more statistics
2620                std::cout<<"Cuts at root node changed objective from "<<babModel->getContinuousObjective()
2621                         <<" to "<<babModel->rootObjectiveAfterCuts()<<std::endl;
2622               
2623                numberGenerators = babModel->numberCutGenerators();
2624                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
2625                  CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
2626                  std::cout<<generator->cutGeneratorName()<<" was tried "
2627                           <<generator->numberTimesEntered()<<" times and created "
2628                           <<generator->numberCutsInTotal()<<" cuts of which "
2629                           <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
2630                  if (generator->timing())
2631                    std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
2632                  else
2633                    std::cout<<std::endl;
2634                }
2635              }
2636              time2 = CoinCpuTime();
2637              totalTime += time2-time1;
2638              // For best solution
2639              double * bestSolution = NULL;
2640              if (babModel->getMinimizationObjValue()<1.0e50&&type==BAB) {
2641                // post process
2642                if (preProcess) {
2643                  int n = saveSolver->getNumCols();
2644                  bestSolution = new double [n];
2645                  process.postProcess(*babModel->solver());
2646                  // Solution now back in saveSolver
2647                  babModel->assignSolver(saveSolver);
2648                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
2649                } else {
2650                  int n = babModel->solver()->getNumCols();
2651                  bestSolution = new double [n];
2652                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
2653                }
2654                checkSOS(babModel, babModel->solver());
2655              }
2656              if (type==STRENGTHEN&&strengthenedModel)
2657                clpSolver = dynamic_cast< OsiClpSolverInterface*> (strengthenedModel);
2658              lpSolver = clpSolver->getModelPtr();
2659              if (numberChanged) {
2660                for (int i=0;i<numberChanged;i++) {
2661                  int iColumn=changed[i];
2662                  clpSolver->setContinuous(iColumn);
2663                }
2664                delete [] changed;
2665              }
2666              if (type==BAB) {
2667                //move best solution (should be there -- but ..)
2668                int n = lpSolver->getNumCols();
2669                if (bestSolution)
2670                  memcpy(lpSolver->primalColumnSolution(),bestSolution,n*sizeof(double));
2671                if (debugFile=="create"&&bestSolution) {
2672                  saveSolution(lpSolver,"debug.file");
2673                }
2674                delete [] bestSolution;
2675                std::string statusName[]={"Finished","Stopped on ","Difficulties",
2676                                          "","","User ctrl-c"};
2677                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
2678                int iStat = babModel->status();
2679                int iStat2 = babModel->secondaryStatus();
2680                if (!noPrinting)
2681                  std::cout<<"Result - "<<statusName[iStat]<<minor[iStat2]
2682                           <<" objective "<<babModel->getObjValue()<<
2683                    " after "<<babModel->getNodeCount()<<" nodes and "
2684                           <<babModel->getIterationCount()<<
2685                    " iterations - took "<<time2-time1<<" seconds"<<std::endl;
2686#ifdef COIN_HAS_ASL
2687                if (usingAmpl) {
2688                  double value = babModel->getObjValue()*lpSolver->getObjSense();
2689                  char buf[300];
2690                  int pos=0;
2691                  if (iStat==0) {
2692                    if (babModel->getObjValue()<1.0e40) {
2693                      pos += sprintf(buf+pos,"optimal," );
2694                    } else {
2695                      // infeasible
2696                      iStat=1;
2697                      pos += sprintf(buf+pos,"infeasible,");
2698                    }
2699                  } else if (iStat==1) {
2700                    if (iStat2!=6)
2701                      iStat=3;
2702                    else
2703                      iStat=4;
2704                    pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
2705                  } else if (iStat==2) {
2706                    iStat = 7;
2707                    pos += sprintf(buf+pos,"stopped on difficulties,");
2708                  } else if (iStat==5) {
2709                    iStat = 3;
2710                    pos += sprintf(buf+pos,"stopped on ctrl-c,");
2711                  } else {
2712                    pos += sprintf(buf+pos,"status unknown,");
2713                    iStat=6;
2714                  }
2715                  info.problemStatus=iStat;
2716                  info.objValue = value;
2717                  if (babModel->getObjValue()<1.0e40) 
2718                    pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
2719                                   value);
2720                  sprintf(buf+pos,"\n%d nodes, %d iterations",
2721                          babModel->getNodeCount(),
2722                          babModel->getIterationCount());
2723                  if (bestSolution) {
2724                    free(info.primalSolution);
2725                    info.primalSolution = (double *) malloc(n*sizeof(double));
2726                    CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
2727                    int numberRows = lpSolver->numberRows();
2728                    free(info.dualSolution);
2729                    info.dualSolution = (double *) malloc(numberRows*sizeof(double));
2730                    CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
2731                  } else {
2732                    info.primalSolution=NULL;
2733                    info.dualSolution=NULL;
2734                  }
2735                  // put buffer into info
2736                  strcpy(info.buffer,buf);
2737                }
2738#endif
2739              } else {
2740                std::cout<<"Model strengthened - now has "<<clpSolver->getNumRows()
2741                         <<" rows"<<std::endl;
2742              }
2743              time1 = time2;
2744              delete babModel;
2745              babModel=NULL;
2746            } else {
2747              std::cout << "** Current model not valid" << std::endl ; 
2748            }
2749            break ;
2750          case IMPORT:
2751            {
2752#ifdef COIN_HAS_ASL
2753              if (!usingAmpl) {
2754#endif
2755                free(priorities);
2756                priorities=NULL;
2757                free(branchDirection);
2758                branchDirection=NULL;
2759                free(pseudoDown);
2760                pseudoDown=NULL;
2761                free(pseudoUp);
2762                pseudoUp=NULL;
2763                free(solutionIn);
2764                solutionIn=NULL;
2765                free(prioritiesIn);
2766                prioritiesIn=NULL;
2767                free(sosStart);
2768                sosStart=NULL;
2769                free(sosIndices);
2770                sosIndices=NULL;
2771                free(sosType);
2772                sosType=NULL;
2773                free(sosReference);
2774                sosReference=NULL;
2775                free(sosPriority);
2776                sosPriority=NULL;
2777#ifdef COIN_HAS_ASL
2778              }
2779#endif               
2780              delete babModel;
2781              babModel=NULL;
2782              // get next field
2783              field = CoinReadGetString(argc,argv);
2784              if (field=="$") {
2785                field = parameters[iParam].stringValue();
2786              } else if (field=="EOL") {
2787                parameters[iParam].printString();
2788                break;
2789              } else {
2790                parameters[iParam].setStringValue(field);
2791              }
2792              std::string fileName;
2793              bool canOpen=false;
2794              if (field=="-") {
2795                // stdin
2796                canOpen=true;
2797                fileName = "-";
2798              } else {
2799                bool absolutePath;
2800                if (dirsep=='/') {
2801                  // non Windows (or cygwin)
2802                  absolutePath=(field[0]=='/');
2803                } else {
2804                  //Windows (non cycgwin)
2805                  absolutePath=(field[0]=='\\');
2806                  // but allow for :
2807                  if (strchr(field.c_str(),':'))
2808                    absolutePath=true;
2809                }
2810                if (absolutePath) {
2811                  fileName = field;
2812                } else if (field[0]=='~') {
2813                  char * environVar = getenv("HOME");
2814                  if (environVar) {
2815                    std::string home(environVar);
2816                    field=field.erase(0,1);
2817                    fileName = home+field;
2818                  } else {
2819                    fileName=field;
2820                  }
2821                } else {
2822                  fileName = directory+field;
2823                }
2824                FILE *fp=fopen(fileName.c_str(),"r");
2825                if (fp) {
2826                  // can open - lets go for it
2827                  fclose(fp);
2828                  canOpen=true;
2829                } else {
2830                  std::cout<<"Unable to open file "<<fileName<<std::endl;
2831                }
2832              }
2833              if (canOpen) {
2834                int status =clpSolver->readMps(fileName.c_str(),
2835                                                   keepImportNames!=0,
2836                                                   allowImportErrors!=0);
2837                if (!status||(status>0&&allowImportErrors)) {
2838                  if (keepImportNames) {
2839                    lengthName = lpSolver->lengthNames();
2840                    rowNames = *(lpSolver->rowNames());
2841                    columnNames = *(lpSolver->columnNames());
2842                  } else {
2843                    lengthName=0;
2844                  }
2845                  goodModel=true;
2846                  //Set integers in clpsolver
2847                  const char * info = lpSolver->integerInformation();
2848                  if (info) {
2849                    int numberColumns = lpSolver->numberColumns();
2850                    int i;
2851                    for (i=0;i<numberColumns;i++) {
2852                      if (info[i]) 
2853                        clpSolver->setInteger(i);
2854                    }
2855                  }
2856                  // sets to all slack (not necessary?)
2857                  lpSolver->createStatus();
2858                  time2 = CoinCpuTime();
2859                  totalTime += time2-time1;
2860                  time1=time2;
2861                  // Go to canned file if just input file
2862                  if (CbcOrClpRead_mode==2&&argc==2) {
2863                    // only if ends .mps
2864                    std::string::size_type loc = fileName.find(".mps") ;
2865                    if (loc != std::string::npos &&
2866                        fileName.length() == loc+3)
2867                    { fileName.replace(loc+1,3,"par") ;
2868                      FILE *fp=fopen(fileName.c_str(),"r");
2869                      if (fp) {
2870                        CbcOrClpReadCommand=fp; // Read from that file
2871                        CbcOrClpRead_mode=-1;
2872                      }
2873                    }
2874                  }
2875                } else {
2876                  // errors
2877                  std::cout<<"There were "<<status<<
2878                    " errors on input"<<std::endl;
2879                }
2880              }
2881            }
2882            break;
2883          case EXPORT:
2884            if (goodModel) {
2885              // get next field
2886              field = CoinReadGetString(argc,argv);
2887              if (field=="$") {
2888                field = parameters[iParam].stringValue();
2889              } else if (field=="EOL") {
2890                parameters[iParam].printString();
2891                break;
2892              } else {
2893                parameters[iParam].setStringValue(field);
2894              }
2895              std::string fileName;
2896              bool canOpen=false;
2897              if (field[0]=='/'||field[0]=='\\') {
2898                fileName = field;
2899              } else if (field[0]=='~') {
2900                char * environVar = getenv("HOME");
2901                if (environVar) {
2902                  std::string home(environVar);
2903                  field=field.erase(0,1);
2904                  fileName = home+field;
2905                } else {
2906                  fileName=field;
2907                }
2908              } else {
2909                fileName = directory+field;
2910              }
2911              FILE *fp=fopen(fileName.c_str(),"w");
2912              if (fp) {
2913                // can open - lets go for it
2914                fclose(fp);
2915                canOpen=true;
2916              } else {
2917                std::cout<<"Unable to open file "<<fileName<<std::endl;
2918              }
2919              if (canOpen) {
2920                // If presolve on then save presolved
2921                bool deleteModel2=false;
2922                ClpSimplex * model2 = lpSolver;
2923#ifdef COIN_HAS_ASL
2924                if (info.numberSos&&doSOS&&usingAmpl) {
2925                  // SOS
2926                  numberSOS = info.numberSos;
2927                  sosStart = info.sosStart;
2928                  sosIndices = info.sosIndices;
2929                  sosReference = info.sosReference;
2930                  preSolve=false;
2931                  clpSolver->setSOSData(numberSOS,info.sosType,sosStart,sosIndices,sosReference);
2932                }
2933#endif
2934                if (preSolve) {
2935                  ClpPresolve pinfo;
2936                  int presolveOptions2 = presolveOptions&~0x40000000;
2937                  if ((presolveOptions2&0xffff)!=0)
2938                    pinfo.setPresolveActions(presolveOptions2);
2939                  if ((printOptions&1)!=0)
2940                    pinfo.statistics();
2941                  double presolveTolerance = 
2942                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
2943                  model2 = 
2944                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
2945                                         true,preSolve);
2946                  if (model2) {
2947                    printf("Saving presolved model on %s\n",
2948                           fileName.c_str());
2949                    deleteModel2=true;
2950                  } else {
2951                    printf("Presolved model looks infeasible - saving original on %s\n",
2952                           fileName.c_str());
2953                    deleteModel2=false;
2954                    model2 = lpSolver;
2955
2956                  }
2957                  model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
2958                  if (deleteModel2)
2959                    delete model2;
2960                } else {
2961                  printf("Saving model on %s\n",
2962                           fileName.c_str());
2963                  if (numberSOS) {
2964                    // Convert names
2965                    int iRow;
2966                    int numberRows=model2->numberRows();
2967                    int iColumn;
2968                    int numberColumns=model2->numberColumns();
2969                   
2970                    char ** rowNames = NULL;
2971                    char ** columnNames = NULL;
2972                    if (model2->lengthNames()) {
2973                      rowNames = new char * [numberRows];
2974                      for (iRow=0;iRow<numberRows;iRow++) {
2975                        rowNames[iRow] = 
2976                          strdup(model2->rowName(iRow).c_str());
2977                      }
2978                     
2979                      columnNames = new char * [numberColumns];
2980                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2981                        columnNames[iColumn] = 
2982                          strdup(model2->columnName(iColumn).c_str());
2983                      }
2984                    }
2985                    clpSolver->writeMpsNative(fileName.c_str(),(const char **) rowNames,(const char **) columnNames,
2986                                              (outputFormat-1)/2,1+((outputFormat-1)&1));
2987                    if (rowNames) {
2988                      for (iRow=0;iRow<numberRows;iRow++) {
2989                        free(rowNames[iRow]);
2990                      }
2991                      delete [] rowNames;
2992                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2993                        free(columnNames[iColumn]);
2994                      }
2995                      delete [] columnNames;
2996                    }
2997                  } else {
2998                    model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
2999                  }
3000                }
3001                time2 = CoinCpuTime();
3002                totalTime += time2-time1;
3003                time1=time2;
3004              }
3005            } else {
3006              std::cout<<"** Current model not valid"<<std::endl;
3007            }
3008            break;
3009          case BASISIN:
3010            if (goodModel) {
3011              // get next field
3012              field = CoinReadGetString(argc,argv);
3013              if (field=="$") {
3014                field = parameters[iParam].stringValue();
3015              } else if (field=="EOL") {
3016                parameters[iParam].printString();
3017                break;
3018              } else {
3019                parameters[iParam].setStringValue(field);
3020              }
3021              std::string fileName;
3022              bool canOpen=false;
3023              if (field=="-") {
3024                // stdin
3025                canOpen=true;
3026                fileName = "-";
3027              } else {
3028                if (field[0]=='/'||field[0]=='\\') {
3029                  fileName = field;
3030                } else if (field[0]=='~') {
3031                  char * environVar = getenv("HOME");
3032                  if (environVar) {
3033                    std::string home(environVar);
3034                    field=field.erase(0,1);
3035                    fileName = home+field;
3036                  } else {
3037                    fileName=field;
3038                  }
3039                } else {
3040                  fileName = directory+field;
3041                }
3042                FILE *fp=fopen(fileName.c_str(),"r");
3043                if (fp) {
3044                  // can open - lets go for it
3045                  fclose(fp);
3046                  canOpen=true;
3047                } else {
3048                  std::cout<<"Unable to open file "<<fileName<<std::endl;
3049                }
3050              }
3051              if (canOpen) {
3052                int values = lpSolver->readBasis(fileName.c_str());
3053                if (values==0)
3054                  basisHasValues=-1;
3055                else
3056                  basisHasValues=1;
3057              }
3058            } else {
3059              std::cout<<"** Current model not valid"<<std::endl;
3060            }
3061            break;
3062          case PRIORITYIN:
3063            if (goodModel) {
3064              // get next field
3065              field = CoinReadGetString(argc,argv);
3066              if (field=="$") {
3067                field = parameters[iParam].stringValue();
3068              } else if (field=="EOL") {
3069                parameters[iParam].printString();
3070                break;
3071              } else {
3072                parameters[iParam].setStringValue(field);
3073              }
3074              std::string fileName;
3075              if (field[0]=='/'||field[0]=='\\') {
3076                fileName = field;
3077              } else if (field[0]=='~') {
3078                char * environVar = getenv("HOME");
3079                if (environVar) {
3080                  std::string home(environVar);
3081                  field=field.erase(0,1);
3082                  fileName = home+field;
3083                } else {
3084                  fileName=field;
3085                }
3086              } else {
3087                fileName = directory+field;
3088              }
3089              FILE *fp=fopen(fileName.c_str(),"r");
3090              if (fp) {
3091                // can open - lets go for it
3092                std::string headings[]={"name","number","direction","priority","up","down",
3093                                        "solution","priin"};
3094                int got[]={-1,-1,-1,-1,-1,-1,-1,-1};
3095                int order[8];
3096                assert(sizeof(got)==sizeof(order));
3097                int nAcross=0;
3098                char line[1000];
3099                int numberColumns = lpSolver->numberColumns();
3100                if (!fgets(line,1000,fp)) {
3101                  std::cout<<"Odd file "<<fileName<<std::endl;
3102                } else {
3103                  char * pos = line;
3104                  char * put = line;
3105                  while (*pos>=' '&&*pos!='\n') {
3106                    if (*pos!=' '&&*pos!='\t') {
3107                      *put=tolower(*pos);
3108                      put++;
3109                    }
3110                    pos++;
3111                  }
3112                  *put='\0';
3113                  pos=line;
3114                  int i;
3115                  bool good=true;
3116                  while (pos) {
3117                    char * comma = strchr(pos,',');
3118                    if (comma)
3119                      *comma='\0';
3120                    for (i=0;i<(int) (sizeof(got)/sizeof(int));i++) {
3121                      if (headings[i]==pos) {
3122                        if (got[i]<0) {
3123                          order[nAcross]=i;
3124                          got[i]=nAcross++;
3125                        } else {
3126                          // duplicate
3127                          good=false;
3128                        }
3129                        break;
3130                      }
3131                    }
3132                    if (i==(int) (sizeof(got)/sizeof(int)))
3133                      good=false;
3134                    if (comma) {
3135                      *comma=',';
3136                      pos=comma+1;
3137                    } else {
3138                      break;
3139                    }
3140                  }
3141                  if (got[0]<0&&got[1]<0)
3142                    good=false;
3143                  if (got[0]>=0&&got[1]>=0)
3144                    good=false;
3145                  if (got[0]>=0&&!lpSolver->lengthNames())
3146                    good=false;
3147                  if (good) {
3148                    char ** columnNames = columnNames = new char * [numberColumns];
3149                    pseudoDown= (double *) malloc(numberColumns*sizeof(double));
3150                    pseudoUp = (double *) malloc(numberColumns*sizeof(double));
3151                    branchDirection = (int *) malloc(numberColumns*sizeof(int));
3152                    priorities= (int *) malloc(numberColumns*sizeof(int));
3153                    free(solutionIn);
3154                    solutionIn=NULL;
3155                    free(prioritiesIn);
3156                    prioritiesIn=NULL;
3157                    int iColumn;
3158                    if (got[6]>=0) {
3159                      solutionIn = (double *) malloc(numberColumns*sizeof(double));
3160                      CoinZeroN(solutionIn,numberColumns);
3161                    }
3162                    if (got[7]>=0) {
3163                      prioritiesIn = (int *) malloc(numberColumns*sizeof(int));
3164                      for (iColumn=0;iColumn<numberColumns;iColumn++) 
3165                        prioritiesIn[iColumn]=10000;
3166                    }
3167                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
3168                      columnNames[iColumn] = 
3169                        strdup(lpSolver->columnName(iColumn).c_str());
3170                      pseudoDown[iColumn]=0.0;
3171                      pseudoUp[iColumn]=0.0;
3172                      branchDirection[iColumn]=0;
3173                      priorities[iColumn]=0;
3174                    }
3175                    int nBadPseudo=0;
3176                    int nBadDir=0;
3177                    int nBadPri=0;
3178                    int nBadName=0;
3179                    int nBadLine=0;
3180                    int nLine=0;
3181                    while (fgets(line,1000,fp)) {
3182                      nLine++;
3183                      iColumn = -1;
3184                      double up =0.0;
3185                      double down=0.0;
3186                      int pri=0;
3187                      int dir=0;
3188                      double solValue=COIN_DBL_MAX;
3189                      int priValue=1000000;
3190                      char * pos = line;
3191                      char * put = line;
3192                      while (*pos>=' '&&*pos!='\n') {
3193                        if (*pos!=' '&&*pos!='\t') {
3194                          *put=tolower(*pos);
3195                          put++;
3196                        }
3197                        pos++;
3198                      }
3199                      *put='\0';
3200                      pos=line;
3201                      for (int i=0;i<nAcross;i++) {
3202                        char * comma = strchr(pos,',');
3203                        if (comma) {
3204                          *comma='\0';
3205                        } else if (i<nAcross-1) {
3206                          nBadLine++;
3207                          break;
3208                        }
3209                        switch (order[i]) {
3210                          // name
3211                        case 0:
3212                          for (iColumn=0;iColumn<numberColumns;iColumn++) {
3213                            if (!strcmp(columnNames[iColumn],pos))
3214                              break;
3215                          }
3216                          if (iColumn==numberColumns)
3217                            iColumn=-1;
3218                          break;
3219                          // number
3220                        case 1:
3221                          iColumn = atoi(pos);
3222                          if (iColumn<0||iColumn>=numberColumns)
3223                            iColumn=-1;
3224                          break;
3225                          // direction
3226                        case 2:
3227                          if (*pos=='D')
3228                            dir=-1;
3229                          else if (*pos=='U')
3230                            dir=1;
3231                          else if (*pos=='N')
3232                            dir=0;
3233                          else if (*pos=='1'&&*(pos+1)=='\0')
3234                            dir=1;
3235                          else if (*pos=='0'&&*(pos+1)=='\0')
3236                            dir=0;
3237                          else if (*pos=='1'&&*(pos+1)=='1'&&*(pos+2)=='\0')
3238                            dir=-1;
3239                          else
3240                            dir=-2; // bad
3241                          break;
3242                          // priority
3243                        case 3:
3244                          pri=atoi(pos);
3245                          break;
3246                          // up
3247                        case 4:
3248                          up = atof(pos);
3249                          break;
3250                          // down
3251                        case 5:
3252                          down = atof(pos);
3253                          break;
3254                          // sol value
3255                        case 6:
3256                          solValue = atof(pos);
3257                          break;
3258                          // priority in value
3259                        case 7:
3260                          priValue = atoi(pos);
3261                          break;
3262                        }
3263                        if (comma) {
3264                          *comma=',';
3265                          pos=comma+1;
3266                        }
3267                      }
3268                      if (iColumn>=0) {
3269                        if (down<0.0) {
3270                          nBadPseudo++;
3271                          down=0.0;
3272                        }
3273                        if (up<0.0) {
3274                          nBadPseudo++;
3275                          up=0.0;
3276                        }
3277                        if (!up)
3278                          up=down;
3279                        if (!down)
3280                          down=up;
3281                        if (dir<-1||dir>1) {
3282                          nBadDir++;
3283                          dir=0;
3284                        }
3285                        if (pri<0) {
3286                          nBadPri++;
3287                          pri=0;
3288                        }
3289                        pseudoDown[iColumn]=down;
3290                        pseudoUp[iColumn]=up;
3291                        branchDirection[iColumn]=dir;
3292                        priorities[iColumn]=pri;
3293                        if (solValue!=COIN_DBL_MAX) {
3294                          assert (solutionIn);
3295                          solutionIn[iColumn]=solValue;
3296                        }
3297                        if (priValue!=1000000) {
3298                          assert (prioritiesIn);
3299                          prioritiesIn[iColumn]=priValue;
3300                        }
3301                      } else {
3302                        nBadName++;
3303                      }
3304                    }
3305                    if (!noPrinting) {
3306                      printf("%d fields and %d records",nAcross,nLine);
3307                      if (nBadPseudo)
3308                        printf(" %d bad pseudo costs",nBadPseudo);
3309                      if (nBadDir)
3310                        printf(" %d bad directions",nBadDir);
3311                      if (nBadPri)
3312                        printf(" %d bad priorities",nBadPri);
3313                      if (nBadName)
3314                        printf(" ** %d records did not match on name/sequence",nBadName);
3315                      printf("\n");
3316                    }
3317                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
3318                      free(columnNames[iColumn]);
3319                    }
3320                    delete [] columnNames;
3321                  } else {
3322                    std::cout<<"Duplicate or unknown keyword - or name/number fields wrong"<<line<<std::endl;
3323                  }
3324                }
3325                fclose(fp);
3326              } else {
3327                std::cout<<"Unable to open file "<<fileName<<std::endl;
3328              }
3329            } else {
3330              std::cout<<"** Current model not valid"<<std::endl;
3331            }
3332            break;
3333          case DEBUG:
3334            if (goodModel) {
3335              delete [] debugValues;
3336              debugValues=NULL;
3337              // get next field
3338              field = CoinReadGetString(argc,argv);
3339              if (field=="$") {
3340                field = parameters[iParam].stringValue();
3341              } else if (field=="EOL") {
3342                parameters[iParam].printString();
3343                break;
3344              } else {
3345                parameters[iParam].setStringValue(field);
3346                debugFile=field;
3347                if (debugFile=="create"||
3348                    debugFile=="createAfterPre") {
3349                  printf("Will create a debug file so this run should be a good one\n");
3350                  break;
3351                }
3352              }
3353              std::string fileName;
3354              if (field[0]=='/'||field[0]=='\\') {
3355                fileName = field;
3356              } else if (field[0]=='~') {
3357                char * environVar = getenv("HOME");
3358                if (environVar) {
3359                  std::string home(environVar);
3360                  field=field.erase(0,1);
3361                  fileName = home+field;
3362                } else {
3363                  fileName=field;
3364                }
3365              } else {
3366                fileName = directory+field;
3367              }
3368              FILE *fp=fopen(fileName.c_str(),"rb");
3369              if (fp) {
3370                // can open - lets go for it
3371                int numRows;
3372                double obj;
3373                fread(&numRows,sizeof(int),1,fp);
3374                fread(&numberDebugValues,sizeof(int),1,fp);
3375                fread(&obj,sizeof(double),1,fp);
3376                debugValues = new double[numberDebugValues+numRows];
3377                fread(debugValues,sizeof(double),numRows,fp);
3378                fread(debugValues,sizeof(double),numRows,fp);
3379                fread(debugValues,sizeof(double),numberDebugValues,fp);
3380                printf("%d doubles read into debugValues\n",numberDebugValues);
3381                fclose(fp);
3382              } else {
3383                std::cout<<"Unable to open file "<<fileName<<std::endl;
3384              }
3385            } else {
3386              std::cout<<"** Current model not valid"<<std::endl;
3387            }
3388            break;
3389          case PRINTMASK:
3390            // get next field
3391            {
3392              std::string name = CoinReadGetString(argc,argv);
3393              if (name!="EOL") {
3394                parameters[iParam].setStringValue(name);
3395                printMask = name;
3396              } else {
3397                parameters[iParam].printString();
3398              }
3399            }
3400            break;
3401          case BASISOUT:
3402            if (goodModel) {
3403              // get next field
3404              field = CoinReadGetString(argc,argv);
3405              if (field=="$") {
3406                field = parameters[iParam].stringValue();
3407              } else if (field=="EOL") {
3408                parameters[iParam].printString();
3409                break;
3410              } else {
3411                parameters[iParam].setStringValue(field);
3412              }
3413              std::string fileName;
3414              bool canOpen=false;
3415              if (field[0]=='/'||field[0]=='\\') {
3416                fileName = field;
3417              } else if (field[0]=='~') {
3418                char * environVar = getenv("HOME");
3419                if (environVar) {
3420                  std::string home(environVar);
3421                  field=field.erase(0,1);
3422                  fileName = home+field;
3423                } else {
3424                  fileName=field;
3425                }
3426              } else {
3427                fileName = directory+field;
3428              }
3429              FILE *fp=fopen(fileName.c_str(),"w");
3430              if (fp) {
3431                // can open - lets go for it
3432                fclose(fp);
3433                canOpen=true;
3434              } else {
3435                std::cout<<"Unable to open file "<<fileName<<std::endl;
3436              }
3437              if (canOpen) {
3438                ClpSimplex * model2 = lpSolver;
3439                model2->writeBasis(fileName.c_str(),outputFormat>1,outputFormat-2);
3440                time2 = CoinCpuTime();
3441                totalTime += time2-time1;
3442                time1=time2;
3443              }
3444            } else {
3445              std::cout<<"** Current model not valid"<<std::endl;
3446            }
3447            break;
3448          case SAVE:
3449            {
3450              // get next field
3451              field = CoinReadGetString(argc,argv);
3452              if (field=="$") {
3453                field = parameters[iParam].stringValue();
3454              } else if (field=="EOL") {
3455                parameters[iParam].printString();
3456                break;
3457              } else {
3458                parameters[iParam].setStringValue(field);
3459              }
3460              std::string fileName;
3461              bool canOpen=false;
3462              if (field[0]=='/'||field[0]=='\\') {
3463                fileName = field;
3464              } else if (field[0]=='~') {
3465                char * environVar = getenv("HOME");
3466                if (environVar) {
3467                  std::string home(environVar);
3468                  field=field.erase(0,1);
3469                  fileName = home+field;
3470                } else {
3471                  fileName=field;
3472                }
3473              } else {
3474                fileName = directory+field;
3475              }
3476              FILE *fp=fopen(fileName.c_str(),"wb");
3477              if (fp) {
3478                // can open - lets go for it
3479                fclose(fp);
3480                canOpen=true;
3481              } else {
3482                std::cout<<"Unable to open file "<<fileName<<std::endl;
3483              }
3484              if (canOpen) {
3485                int status;
3486                // If presolve on then save presolved
3487                bool deleteModel2=false;
3488                ClpSimplex * model2 = lpSolver;
3489                if (preSolve) {
3490                  ClpPresolve pinfo;
3491                  double presolveTolerance = 
3492                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
3493                  model2 = 
3494                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
3495                                         false,preSolve);
3496                  if (model2) {
3497                    printf("Saving presolved model on %s\n",
3498                           fileName.c_str());
3499                    deleteModel2=true;
3500                  } else {
3501                    printf("Presolved model looks infeasible - saving original on %s\n",
3502                           fileName.c_str());
3503                    deleteModel2=false;
3504                    model2 = lpSolver;
3505
3506                  }
3507                } else {
3508                  printf("Saving model on %s\n",
3509                           fileName.c_str());
3510                }
3511                status =model2->saveModel(fileName.c_str());
3512                if (deleteModel2)
3513                  delete model2;
3514                if (!status) {
3515                  goodModel=true;
3516                  time2 = CoinCpuTime();
3517                  totalTime += time2-time1;
3518                  time1=time2;
3519                } else {
3520                  // errors
3521                  std::cout<<"There were errors on output"<<std::endl;
3522                }
3523              }
3524            }
3525            break;
3526          case RESTORE:
3527            {
3528              // get next field
3529              field = CoinReadGetString(argc,argv);
3530              if (field=="$") {
3531                field = parameters[iParam].stringValue();
3532              } else if (field=="EOL") {
3533                parameters[iParam].printString();
3534                break;
3535              } else {
3536                parameters[iParam].setStringValue(field);
3537              }
3538              std::string fileName;
3539              bool canOpen=false;
3540              if (field[0]=='/'||field[0]=='\\') {
3541                fileName = field;
3542              } else if (field[0]=='~') {
3543                char * environVar = getenv("HOME");
3544                if (environVar) {
3545                  std::string home(environVar);
3546                  field=field.erase(0,1);
3547                  fileName = home+field;
3548                } else {
3549                  fileName=field;
3550                }
3551              } else {
3552                fileName = directory+field;
3553              }
3554              FILE *fp=fopen(fileName.c_str(),"rb");
3555              if (fp) {
3556                // can open - lets go for it
3557                fclose(fp);
3558                canOpen=true;
3559              } else {
3560                std::cout<<"Unable to open file "<<fileName<<std::endl;
3561              }
3562              if (canOpen) {
3563                int status =lpSolver->restoreModel(fileName.c_str());
3564                if (!status) {
3565                  goodModel=true;
3566                  time2 = CoinCpuTime();
3567                  totalTime += time2-time1;
3568                  time1=time2;
3569                } else {
3570                  // errors
3571                  std::cout<<"There were errors on input"<<std::endl;
3572                }
3573              }
3574            }
3575            break;
3576          case MAXIMIZE:
3577            lpSolver->setOptimizationDirection(-1);
3578            break;
3579          case MINIMIZE:
3580            lpSolver->setOptimizationDirection(1);
3581            break;
3582          case ALLSLACK:
3583            lpSolver->allSlackBasis(true);
3584            break;
3585          case REVERSE:
3586            if (goodModel) {
3587              int iColumn;
3588              int numberColumns=lpSolver->numberColumns();
3589              double * dualColumnSolution = 
3590                lpSolver->dualColumnSolution();
3591              ClpObjective * obj = lpSolver->objectiveAsObject();
3592              assert(dynamic_cast<ClpLinearObjective *> (obj));
3593              double offset;
3594              double * objective = obj->gradient(NULL,NULL,offset,true);
3595              for (iColumn=0;iColumn<numberColumns;iColumn++) {
3596                dualColumnSolution[iColumn] = dualColumnSolution[iColumn];
3597                objective[iColumn] = -objective[iColumn];
3598              }
3599              int iRow;
3600              int numberRows=lpSolver->numberRows();
3601              double * dualRowSolution = 
3602                lpSolver->dualRowSolution();
3603              for (iRow=0;iRow<numberRows;iRow++) 
3604                dualRowSolution[iRow] = dualRowSolution[iRow];
3605            }
3606            break;
3607          case DIRECTORY:
3608            {
3609              std::string name = CoinReadGetString(argc,argv);
3610              if (name!="EOL") {
3611                int length=name.length();
3612                if (name[length-1]=='/'||name[length-1]=='\\')
3613                  directory=name;
3614                else
3615                  directory = name+"/";
3616                parameters[iParam].setStringValue(directory);
3617              } else {
3618                parameters[iParam].printString();
3619              }
3620            }
3621            break;
3622          case STDIN:
3623            CbcOrClpRead_mode=-1;
3624            break;
3625          case NETLIB_DUAL:
3626          case NETLIB_EITHER:
3627          case NETLIB_BARRIER:
3628          case NETLIB_PRIMAL:
3629          case NETLIB_TUNE:
3630            {
3631              // create fields for unitTest
3632              const char * fields[4];
3633              int nFields=2;
3634              fields[0]="fake main from unitTest";
3635              fields[1]="-netlib";
3636              if (directory!=defaultDirectory) {
3637                fields[2]="-netlibDir";
3638                fields[3]=directory.c_str();
3639                nFields=4;
3640              }
3641              int algorithm;
3642              if (type==NETLIB_DUAL) {
3643                std::cerr<<"Doing netlib with dual algorithm"<<std::endl;
3644                algorithm =0;
3645              } else if (type==NETLIB_BARRIER) {
3646                std::cerr<<"Doing netlib with barrier algorithm"<<std::endl;
3647                algorithm =2;
3648              } else if (type==NETLIB_EITHER) {
3649                std::cerr<<"Doing netlib with dual or primal algorithm"<<std::endl;
3650                algorithm =3;
3651              } else if (type==NETLIB_TUNE) {
3652                std::cerr<<"Doing netlib with best algorithm!"<<std::endl;
3653                algorithm =5;
3654                // uncomment next to get active tuning
3655                // algorithm=6;
3656              } else {
3657                std::cerr<<"Doing netlib with primal agorithm"<<std::endl;
3658                algorithm=1;
3659              }
3660              int specialOptions = lpSolver->specialOptions();
3661              lpSolver->setSpecialOptions(0);
3662              mainTest(nFields,fields,algorithm,*lpSolver,
3663                       (preSolve!=0),specialOptions,doVector!=0);
3664            }
3665            break;
3666          case UNITTEST:
3667            {
3668              // create fields for unitTest
3669              const char * fields[3];
3670              int nFields=1;
3671              fields[0]="fake main from unitTest";
3672              if (directory!=defaultDirectory) {
3673                fields[1]="-mpsDir";
3674                fields[2]=directory.c_str();
3675                nFields=3;
3676              }
3677              mainTest(nFields,fields,false,*lpSolver,(preSolve!=0),
3678                       false,doVector!=0);
3679            }
3680            break;
3681          case FAKEBOUND:
3682            if (goodModel) {
3683              // get bound
3684              double value = CoinReadGetDoubleField(argc,argv,&valid);
3685              if (!valid) {
3686                std::cout<<"Setting "<<parameters[iParam].name()<<
3687                  " to DEBUG "<<value<<std::endl;
3688                int iRow;
3689                int numberRows=lpSolver->numberRows();
3690                double * rowLower = lpSolver->rowLower();
3691                double * rowUpper = lpSolver->rowUpper();
3692                for (iRow=0;iRow<numberRows;iRow++) {
3693                  // leave free ones for now
3694                  if (rowLower[iRow]>-1.0e20||rowUpper[iRow]<1.0e20) {
3695                    rowLower[iRow]=CoinMax(rowLower[iRow],-value);
3696                    rowUpper[iRow]=CoinMin(rowUpper[iRow],value);
3697                  }
3698                }
3699                int iColumn;
3700                int numberColumns=lpSolver->numberColumns();
3701                double * columnLower = lpSolver->columnLower();
3702                double * columnUpper = lpSolver->columnUpper();
3703                for (iColumn=0;iColumn<numberColumns;iColumn++) {
3704                  // leave free ones for now
3705                  if (columnLower[iColumn]>-1.0e20||
3706                      columnUpper[iColumn]<1.0e20) {
3707                    columnLower[iColumn]=CoinMax(columnLower[iColumn],-value);
3708                    columnUpper[iColumn]=CoinMin(columnUpper[iColumn],value);
3709                  }
3710                }
3711              } else if (valid==1) {
3712                abort();
3713              } else {
3714                std::cout<<"enter value for "<<parameters[iParam].name()<<
3715                  std::endl;
3716              }
3717            }
3718            break;
3719          case REALLY_SCALE:
3720            if (goodModel) {
3721              ClpSimplex newModel(*lpSolver,
3722                                  lpSolver->scalingFlag());
3723              printf("model really really scaled\n");
3724              *lpSolver=newModel;
3725            }
3726            break;
3727          case USERCLP:
3728            // Replace the sample code by whatever you want
3729            if (goodModel) {
3730              printf("Dummy user clp code - model has %d rows and %d columns\n",
3731                     lpSolver->numberRows(),lpSolver->numberColumns());
3732            }
3733            break;
3734          case USERCBC:
3735            // Replace the sample code by whatever you want
3736            if (goodModel) {
3737#ifndef USER_HAS_FAKE_MAIN
3738              printf("Dummy user cbc code - model has %d rows and %d columns\n",
3739                     model.getNumRows(),model.getNumCols());
3740              // Reduce printout
3741              model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
3742              // Do complete search
3743              model.branchAndBound();
3744#ifdef COIN_HAS_ASL
3745              double objectiveValue=model.getMinimizationObjValue();
3746              int iStat = model.status();
3747              int iStat2 = model.secondaryStatus();
3748#endif
3749#else
3750              // Way of using an existing piece of code
3751              OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model.solver());
3752              ClpSimplex * lpSolver = clpSolver->getModelPtr();
3753              // set time from integer model
3754              double timeToGo = model.getMaximumSeconds();
3755              lpSolver->setMaximumSeconds(timeToGo);
3756              fakeMain(*lpSolver,*clpSolver,model);
3757#ifdef COIN_HAS_ASL
3758              // My actual usage has objective only in clpSolver
3759              double objectiveValue=clpSolver->getObjValue();
3760              int iStat = lpSolver->status();
3761              int iStat2 = lpSolver->secondaryStatus();
3762#endif
3763#endif
3764              // make sure solution back in correct place
3765              clpSolver = dynamic_cast< OsiClpSolverInterface*> (model.solver());
3766              lpSolver = clpSolver->getModelPtr();
3767#ifdef COIN_HAS_ASL
3768              if (usingAmpl) {
3769                int n = clpSolver->getNumCols();
3770                double value = objectiveValue*lpSolver->getObjSense();
3771                char buf[300];
3772                int pos=0;
3773                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
3774                if (iStat==0) {
3775                  if (objectiveValue<1.0e40) {
3776                    pos += sprintf(buf+pos,"optimal," );
3777                  } else {
3778                    // infeasible
3779                    iStat=1;
3780                    pos += sprintf(buf+pos,"infeasible,");
3781                  }
3782                } else if (iStat==1) {
3783                  if (iStat2!=6)
3784                    iStat=3;
3785                  else
3786                    iStat=4;
3787                  pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
3788                } else if (iStat==2) {
3789                  iStat = 7;
3790                  pos += sprintf(buf+pos,"stopped on difficulties,");
3791                } else if (iStat==5) {
3792                  iStat = 3;
3793                  pos += sprintf(buf+pos,"stopped on ctrl-c,");
3794                } else {
3795                  pos += sprintf(buf+pos,"status unknown,");
3796                  iStat=6;
3797                }
3798                info.problemStatus=iStat;
3799                info.objValue = value;
3800                if (objectiveValue<1.0e40) 
3801                  pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
3802                                 value);
3803                sprintf(buf+pos,"\n%d nodes, %d iterations",
3804                        model.getNodeCount(),
3805                        model.getIterationCount());
3806                if (objectiveValue<1.0e50) {
3807                  free(info.primalSolution);
3808                  info.primalSolution = (double *) malloc(n*sizeof(double));
3809                  CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
3810                  int numberRows = lpSolver->numberRows();
3811                  free(info.dualSolution);
3812                  info.dualSolution = (double *) malloc(numberRows*sizeof(double));
3813                  CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
3814                } else {
3815                  info.primalSolution=NULL;
3816                  info.dualSolution=NULL;
3817                }
3818                // put buffer into info
3819                strcpy(info.buffer,buf);
3820              }
3821#endif
3822            }
3823            break;
3824          case HELP:
3825            std::cout<<"Coin Solver version "<<CBCVERSION
3826                     <<", build "<<__DATE__<<std::endl;
3827            std::cout<<"Non default values:-"<<std::endl;
3828            std::cout<<"Perturbation "<<lpSolver->perturbation()<<" (default 100)"
3829                     <<std::endl;
3830            CoinReadPrintit(
3831                    "Presolve being done with 5 passes\n\
3832Dual steepest edge steep/partial on matrix shape and factorization density\n\
3833Clpnnnn taken out of messages\n\
3834If Factorization frequency default then done on size of matrix\n\n\
3835(-)unitTest, (-)netlib or (-)netlibp will do standard tests\n\n\
3836You can switch to interactive mode at any time so\n\
3837clp watson.mps -scaling off -primalsimplex\nis the same as\n\
3838clp watson.mps -\nscaling off\nprimalsimplex"
3839                    );
3840            break;
3841          case SOLUTION:
3842            if (goodModel) {
3843              // get next field
3844              field = CoinReadGetString(argc,argv);
3845              if (field=="$") {
3846                field = parameters[iParam].stringValue();
3847              } else if (field=="EOL") {
3848                parameters[iParam].printString();
3849                break;
3850              } else {
3851                parameters[iParam].setStringValue(field);
3852              }
3853              std::string fileName;
3854              FILE *fp=NULL;
3855              if (field=="-"||field=="EOL"||field=="stdout") {
3856                // stdout
3857                fp=stdout;
3858              } else if (field=="stderr") {
3859                // stderr
3860                fp=stderr;
3861              } else {
3862                if (field[0]=='/'||field[0]=='\\') {
3863                  fileName = field;
3864                } else if (field[0]=='~') {
3865                  char * environVar = getenv("HOME");
3866                  if (environVar) {
3867                    std::string home(environVar);
3868                    field=field.erase(0,1);
3869                    fileName = home+field;
3870                  } else {
3871                    fileName=field;
3872                  }
3873                } else {
3874                  fileName = directory+field;
3875                }
3876                fp=fopen(fileName.c_str(),"w");
3877              }
3878              if (fp) {
3879                // make fancy later on
3880                int iRow;
3881                int numberRows=lpSolver->numberRows();
3882                double * dualRowSolution = lpSolver->dualRowSolution();
3883                double * primalRowSolution = 
3884                  lpSolver->primalRowSolution();
3885                double * rowLower = lpSolver->rowLower();
3886                double * rowUpper = lpSolver->rowUpper();
3887                double primalTolerance = lpSolver->primalTolerance();
3888                char format[6];
3889                sprintf(format,"%%-%ds",CoinMax(lengthName,8));
3890                bool doMask = (printMask!=""&&lengthName);
3891                int * maskStarts=NULL;
3892                int maxMasks=0;
3893                char ** masks =NULL;
3894                if (doMask) {
3895                  int nAst =0;
3896                  const char * pMask2 = printMask.c_str();
3897                  char pMask[100];
3898                  int iChar;
3899                  int lengthMask = strlen(pMask2);
3900                  assert (lengthMask<100);
3901                  if (*pMask2=='"') {
3902                    if (pMask2[lengthMask-1]!='"') {
3903                      printf("mismatched \" in mask %s\n",pMask2);
3904                      break;
3905                    } else {
3906                      strcpy(pMask,pMask2+1);
3907                      *strchr(pMask,'"')='\0';
3908                    }
3909                  } else if (*pMask2=='\'') {
3910                    if (pMask2[lengthMask-1]!='\'') {
3911                      printf("mismatched ' in mask %s\n",pMask2);
3912                      break;
3913                    } else {
3914                      strcpy(pMask,pMask2+1);
3915                      *strchr(pMask,'\'')='\0';
3916                    }
3917                  } else {
3918                    strcpy(pMask,pMask2);
3919                  }
3920                  if (lengthMask>lengthName) {
3921                    printf("mask %s too long - skipping\n",pMask);
3922                    break;
3923                  }
3924                  maxMasks = 1;
3925                  for (iChar=0;iChar<lengthMask;iChar++) {
3926                    if (pMask[iChar]=='*') {
3927                      nAst++;
3928                      maxMasks *= (lengthName+1);
3929                    }
3930                  }
3931                  int nEntries = 1;
3932                  maskStarts = new int[lengthName+2];
3933                  masks = new char * [maxMasks];
3934                  char ** newMasks = new char * [maxMasks];
3935                  int i;
3936                  for (i=0;i<maxMasks;i++) {
3937                    masks[i] = new char[lengthName+1];
3938                    newMasks[i] = new char[lengthName+1];
3939                  }
3940                  strcpy(masks[0],pMask);
3941                  for (int iAst=0;iAst<nAst;iAst++) {
3942                    int nOldEntries = nEntries;
3943                    nEntries=0;
3944                    for (int iEntry = 0;iEntry<nOldEntries;iEntry++) {
3945                      char * oldMask = masks[iEntry];
3946                      char * ast = strchr(oldMask,'*');
3947                      assert (ast);
3948                      int length = strlen(oldMask)-1;
3949                      int nBefore = ast-oldMask;
3950                      int nAfter = length-nBefore;
3951                      // and add null
3952                      nAfter++;
3953                      for (int i=0;i<=lengthName-length;i++) {
3954                        char * maskOut = newMasks[nEntries];
3955                        memcpy(maskOut,oldMask,nBefore);
3956                        for (int k=0;k<i;k++) 
3957                          maskOut[k+nBefore]='?';
3958                        memcpy(maskOut+nBefore+i,ast+1,nAfter);
3959                        nEntries++;
3960                        assert (nEntries<=maxMasks);
3961                      }
3962                    }
3963                    char ** temp = masks;
3964                    masks = newMasks;
3965                    newMasks = temp;
3966                  }
3967                  // Now extend and sort
3968                  int * sort = new int[nEntries];
3969                  for (i=0;i<nEntries;i++) {
3970                    char * maskThis = masks[i];
3971                    int length = strlen(maskThis);
3972                    while (maskThis[length-1]==' ')
3973                      length--;
3974                    maskThis[length]='\0';
3975                    sort[i]=length;
3976                  }
3977                  CoinSort_2(sort,sort+nEntries,masks);
3978                  int lastLength=-1;
3979                  for (i=0;i<nEntries;i++) {
3980                    int length = sort[i];
3981                    while (length>lastLength) 
3982                      maskStarts[++lastLength] = i;
3983                  }
3984                  maskStarts[++lastLength]=nEntries;
3985                  delete [] sort;
3986                  for (i=0;i<maxMasks;i++)
3987                    delete [] newMasks[i];
3988                  delete [] newMasks;
3989                }
3990                if (printMode>2) {
3991                  for (iRow=0;iRow<numberRows;iRow++) {
3992                    int type=printMode-3;
3993                    if (primalRowSolution[iRow]>rowUpper[iRow]+primalTolerance||
3994                        primalRowSolution[iRow]<rowLower[iRow]-primalTolerance) {
3995                      fprintf(fp,"** ");
3996                      type=2;
3997                    } else if (fabs(primalRowSolution[iRow])>1.0e-8) {
3998                      type=1;
3999                    } else if (numberRows<50) {
4000                      type=3;
4001                    }
4002                    if (doMask&&!maskMatches(maskStarts,masks,rowNames[iRow]))
4003                      type =0;
4004                    if (type) {
4005                      fprintf(fp,"%7d ",iRow);
4006                      if (lengthName)
4007                        fprintf(fp,format,rowNames[iRow].c_str());
4008                      fprintf(fp,"%15.8g        %15.8g\n",primalRowSolution[iRow],
4009                              dualRowSolution[iRow]);
4010                    }
4011                  }
4012                }
4013                int iColumn;
4014                int numberColumns=lpSolver->numberColumns();
4015                double * dualColumnSolution = 
4016                  lpSolver->dualColumnSolution();
4017                double * primalColumnSolution = 
4018                  lpSolver->primalColumnSolution();
4019                double * columnLower = lpSolver->columnLower();
4020                double * columnUpper = lpSolver->columnUpper();
4021                if (printMode!=2) {
4022                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4023                    int type=(printMode>3) ? 1 :0;
4024                    if (primalColumnSolution[iColumn]>columnUpper[iColumn]+primalTolerance||
4025                        primalColumnSolution[iColumn]<columnLower[iColumn]-primalTolerance) {
4026                      fprintf(fp,"** ");
4027                      type=2;
4028                    } else if (fabs(primalColumnSolution[iColumn])>1.0e-8) {
4029                      type=1;
4030                    } else if (numberColumns<50) {
4031                      type=3;
4032                    }
4033                    // see if integer
4034                    if ((!lpSolver->isInteger(iColumn)||fabs(primalColumnSolution[iColumn])<1.0e-8)
4035                         &&printMode==1)
4036                      type=0;
4037                    if (doMask&&!maskMatches(maskStarts,masks,
4038                                             columnNames[iColumn]))
4039                      type =0;
4040                    if (type) {
4041                      fprintf(fp,"%7d ",iColumn);
4042                      if (lengthName)
4043                        fprintf(fp,format,columnNames[iColumn].c_str());
4044                      fprintf(fp,"%15.8g        %15.8g\n",
4045                              primalColumnSolution[iColumn],
4046                              dualColumnSolution[iColumn]);
4047                    }
4048                  }
4049                } else {
4050                  // special format suitable for OsiRowCutDebugger
4051                  int n=0;
4052                  bool comma=false;
4053                  bool newLine=false;
4054                  fprintf(fp,"\tint intIndicesV[]={\n");
4055                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4056                    if(primalColumnSolution[iColumn]>0.5&&model.solver()->isInteger(iColumn)) {
4057                      if (comma)
4058                        fprintf(fp,",");
4059                      if (newLine)
4060                        fprintf(fp,"\n");
4061                      fprintf(fp,"%d ",iColumn);
4062                      comma=true;
4063                      newLine=false;
4064                      n++;
4065                      if (n==10) {
4066                        n=0;
4067                        newLine=true;
4068                      }
4069                    }
4070                  }
4071                  fprintf(fp,"};\n");
4072                  n=0;
4073                  comma=false;
4074                  newLine=false;
4075                  fprintf(fp,"\tdouble intSolnV[]={\n");
4076                  for ( iColumn=0;iColumn<numberColumns;iColumn++) {
4077                    if(primalColumnSolution[iColumn]>0.5&&model.solver()->isInteger(iColumn)) {
4078                      if (comma)
4079                        fprintf(fp,",");
4080                      if (newLine)
4081                        fprintf(fp,"\n");
4082                      int value = (int) (primalColumnSolution[iColumn]+0.5);
4083                      fprintf(fp,"%d. ",value);
4084                      comma=true;
4085                      newLine=false;
4086                      n++;
4087                      if (n==10) {
4088                        n=0;
4089                        newLine=true;
4090                      }
4091                    }
4092                  }
4093                  fprintf(fp,"};\n");
4094                }
4095                if (fp!=stdout)
4096                  fclose(fp);
4097                if (masks) {
4098                  delete [] maskStarts;
4099                  for (int i=0;i<maxMasks;i++)
4100                    delete [] masks[i];
4101                  delete [] masks;
4102                }
4103              } else {
4104                std::cout<<"Unable to open file "<<fileName<<std::endl;
4105              }
4106            } else {
4107              std::cout<<"** Current model not valid"<<std::endl;
4108             
4109            }
4110            break;
4111          case SAVESOL:
4112            if (goodModel) {
4113              // get next field
4114              field = CoinReadGetString(argc,argv);
4115              if (field=="$") {
4116                field = parameters[iParam].stringValue();
4117              } else if (field=="EOL") {
4118                parameters[iParam].printString();
4119                break;
4120              } else {
4121                parameters[iParam].setStringValue(field);
4122              }
4123              std::string fileName;
4124              if (field[0]=='/'||field[0]=='\\') {
4125                fileName = field;
4126              } else if (field[0]=='~') {
4127                char * environVar = getenv("HOME");
4128                if (environVar) {
4129                  std::string home(environVar);
4130                  field=field.erase(0,1);
4131                  fileName = home+field;
4132                } else {
4133                  fileName=field;
4134                }
4135              } else {
4136                fileName = directory+field;
4137              }
4138              saveSolution(lpSolver,fileName);
4139            } else {
4140              std::cout<<"** Current model not valid"<<std::endl;
4141             
4142            }
4143            break;
4144          case DUMMY:
4145            break;
4146          default:
4147            abort();
4148          }
4149        } 
4150      } else if (!numberMatches) {
4151        std::cout<<"No match for "<<field<<" - ? for list of commands"
4152                 <<std::endl;
4153      } else if (numberMatches==1) {
4154        if (!numberQuery) {
4155          std::cout<<"Short match for "<<field<<" - completion: ";
4156          std::cout<<parameters[firstMatch].matchName()<<std::endl;
4157        } else if (numberQuery) {
4158          std::cout<<parameters[firstMatch].matchName()<<" : ";
4159          std::cout<<parameters[firstMatch].shortHelp()<<std::endl;
4160          if (numberQuery>=2) 
4161            parameters[firstMatch].printLongHelp();
4162        }
4163      } else {
4164        if (!numberQuery) 
4165          std::cout<<"Multiple matches for "<<field<<" - possible completions:"
4166                   <<std::endl;
4167        else
4168          std::cout<<"Completions of "<<field<<":"<<std::endl;
4169        for ( iParam=0; iParam<numberParameters; iParam++ ) {
4170          int match = parameters[iParam].matches(field);
4171          if (match&&parameters[iParam].displayThis()) {
4172            std::cout<<parameters[iParam].matchName();
4173            if (numberQuery>=2) 
4174              std::cout<<" : "<<parameters[iParam].shortHelp();
4175            std::cout<<std::endl;
4176          }
4177        }
4178      }
4179    }
4180  }
4181  // By now all memory should be freed
4182#ifdef DMALLOC
4183  dmalloc_log_unfreed();
4184  dmalloc_shutdown();
4185#endif
4186  return 0;
4187}   
4188static void breakdown(const char * name, int numberLook, const double * region)
4189{
4190  double range[] = {
4191    -COIN_DBL_MAX,
4192    -1.0e15,-1.0e11,-1.0e8,-1.0e5,-1.0e4,-1.0e3,-1.0e2,-1.0e1,
4193    -1.0,
4194    -1.0e-1,-1.0e-2,-1.0e-3,-1.0e-4,-1.0e-5,-1.0e-8,-1.0e-11,-1.0e-15,
4195    0.0,
4196    1.0e-15,1.0e-11,1.0e-8,1.0e-5,1.0e-4,1.0e-3,1.0e-2,1.0e-1,
4197    1.0,
4198    1.0e1,1.0e2,1.0e3,1.0e4,1.0e5,1.0e8,1.0e11,1.0e15,
4199    COIN_DBL_MAX};
4200  int nRanges = (int) (sizeof(range)/sizeof(double));
4201  int * number = new int[nRanges];
4202  memset(number,0,nRanges*sizeof(int));
4203  int * numberExact = new int[nRanges];
4204  memset(numberExact,0,nRanges*sizeof(int));
4205  int i;
4206  for ( i=0;i<numberLook;i++) {
4207    double value = region[i];
4208    for (int j=0;j<nRanges;j++) {
4209      if (value==range[j]) {
4210        numberExact[j]++;
4211        break;
4212      } else if (value<range[j]) {
4213        number[j]++;
4214        break;
4215      }
4216    }
4217  }
4218  printf("\n%s has %d entries\n",name,numberLook);
4219  for (i=0;i<nRanges;i++) {
4220    if (number[i]) 
4221      printf("%d between %g and %g",number[i],range[i-1],range[i]);
4222    if (numberExact[i]) {
4223      if (number[i])
4224        printf(", ");
4225      printf("%d exactly at %g",numberExact[i],range[i]);
4226    }
4227    if (number[i]+numberExact[i])
4228      printf("\n");
4229  }
4230  delete [] number;
4231  delete [] numberExact;
4232}
4233static void statistics(ClpSimplex * originalModel, ClpSimplex * model)
4234{
4235  int numberColumns = originalModel->numberColumns();
4236  const char * integerInformation  = originalModel->integerInformation(); 
4237  const double * columnLower = originalModel->columnLower();
4238  const double * columnUpper = originalModel->columnUpper();
4239  int numberIntegers=0;
4240  int numberBinary=0;
4241  int iRow,iColumn;
4242  if (integerInformation) {
4243    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4244      if (integerInformation[iColumn]) {
4245        if (columnUpper[iColumn]>columnLower[iColumn]) {
4246          numberIntegers++;
4247          if (columnUpper[iColumn]==0.0&&columnLower[iColumn]==1) 
4248            numberBinary++;
4249        }
4250      }
4251    }
4252  }
4253  numberColumns = model->numberColumns();
4254  int numberRows = model->numberRows();
4255  columnLower = model->columnLower();
4256  columnUpper = model->columnUpper();
4257  const double * rowLower = model->rowLower();
4258  const double * rowUpper = model->rowUpper();
4259  const double * objective = model->objective();
4260  CoinPackedMatrix * matrix = model->matrix();
4261  CoinBigIndex numberElements = matrix->getNumElements();
4262  const int * columnLength = matrix->getVectorLengths();
4263  //const CoinBigIndex * columnStart = matrix->getVectorStarts();
4264  const double * elementByColumn = matrix->getElements();
4265  int * number = new int[numberRows+1];
4266  memset(number,0,(numberRows+1)*sizeof(int));
4267  int numberObjSingletons=0;
4268  /* cType
4269     0 0/inf, 1 0/up, 2 lo/inf, 3 lo/up, 4 free, 5 fix, 6 -inf/0, 7 -inf/up,
4270     8 0/1
4271  */ 
4272  int cType[9];
4273  std::string cName[]={"0.0->inf,","0.0->up,","lo->inf,","lo->up,","free,","fixed,","-inf->0.0,",
4274                       "-inf->up,","0.0->1.0"};
4275  int nObjective=0;
4276  memset(cType,0,sizeof(cType));
4277  for (iColumn=0;iColumn<numberColumns;iColumn++) {
4278    int length=columnLength[iColumn];
4279    if (length==1&&objective[iColumn])
4280      numberObjSingletons++;
4281    number[length]++;
4282    if (objective[iColumn])
4283      nObjective++;
4284    if (columnLower[iColumn]>-1.0e20) {
4285      if (columnLower[iColumn]==0.0) {
4286        if (columnUpper[iColumn]>1.0e20)
4287          cType[0]++;
4288        else if (columnUpper[iColumn]==1.0)
4289          cType[8]++;
4290        else if (columnUpper[iColumn]==0.0)
4291          cType[5]++;
4292        else
4293          cType[1]++;
4294      } else {
4295        if (columnUpper[iColumn]>1.0e20) 
4296          cType[2]++;
4297        else if (columnUpper[iColumn]==columnLower[iColumn])
4298          cType[5]++;
4299        else
4300          cType[3]++;
4301      }
4302    } else {
4303      if (columnUpper[iColumn]>1.0e20) 
4304        cType[4]++;
4305      else if (columnUpper[iColumn]==0.0) 
4306        cType[6]++;
4307      else
4308        cType[7]++;
4309    }
4310  }
4311  /* rType
4312     0 E 0, 1 E 1, 2 E -1, 3 E other, 4 G 0, 5 G 1, 6 G other,
4313     7 L 0,  8 L 1, 9 L other, 10 Range 0/1, 11 Range other, 12 free
4314  */ 
4315  int rType[13];
4316  std::string rName[]={"E 0.0,","E 1.0,","E -1.0,","E other,","G 0.0,","G 1.0,","G other,",
4317                       "L 0.0,","L 1.0,","L other,","Range 0.0->1.0,","Range other,","Free"};
4318  memset(rType,0,sizeof(rType));
4319  for (iRow=0;iRow<numberRows;iRow++) {
4320    if (rowLower[iRow]>-1.0e20) {
4321      if (rowLower[iRow]==0.0) {
4322        if (rowUpper[iRow]>1.0e20)
4323          rType[4]++;
4324        else if (rowUpper[iRow]==1.0)
4325          rType[10]++;
4326        else if (rowUpper[iRow]==0.0)
4327          rType[0]++;
4328        else
4329          rType[11]++;
4330      } else if (rowLower[iRow]==1.0) {
4331        if (rowUpper[iRow]>1.0e20) 
4332          rType[5]++;
4333        else if (rowUpper[iRow]==rowLower[iRow])
4334          rType[1]++;
4335        else
4336          rType[11]++;
4337      } else if (rowLower[iRow]==-1.0) {
4338        if (rowUpper[iRow]>1.0e20) 
4339          rType[6]++;
4340        else if (rowUpper[iRow]==rowLower[iRow])
4341          rType[2]++;
4342        else
4343          rType[11]++;
4344      } else {
4345        if (rowUpper[iRow]>1.0e20) 
4346          rType[6]++;
4347        else if (rowUpper[iRow]==rowLower[iRow])
4348          rType[3]++;
4349        else
4350          rType[11]++;
4351      }
4352    } else {
4353      if (rowUpper[iRow]>1.0e20) 
4354        rType[12]++;
4355      else if (rowUpper[iRow]==0.0) 
4356        rType[7]++;
4357      else if (rowUpper[iRow]==1.0) 
4358        rType[8]++;
4359      else
4360        rType[9]++;
4361    }
4362  }
4363  // Basic statistics
4364  printf("\n\nProblem has %d rows, %d columns (%d with objective) and %d elements\n",
4365         numberRows,numberColumns,nObjective,numberElements);
4366  if (number[0]+number[1]) {
4367    printf("There are ");
4368    if (numberObjSingletons)
4369      printf("%d singletons with objective ",numberObjSingletons);
4370    int numberNoObj = number[1]-numberObjSingletons;
4371    if (numberNoObj)
4372      printf("%d singletons with no objective ",numberNoObj);
4373    if (number[0])
4374      printf("** %d columns have no entries",number[0]);
4375    printf("\n");
4376  }
4377  printf("Column breakdown:\n");
4378  int k;
4379  for (k=0;k<(int) (sizeof(cType)/sizeof(int));k++) {
4380    printf("%d of type %s ",cType[k],cName[k].c_str());
4381    if (((k+1)%3)==0)
4382      printf("\n");
4383  }
4384  if ((k%3)!=0)
4385    printf("\n");
4386  printf("Row breakdown:\n");
4387  for (k=0;k<(int) (sizeof(rType)/sizeof(int));k++) {
4388    printf("%d of type %s ",rType[k],rName[k].c_str());
4389    if (((k+1)%3)==0)
4390      printf("\n");
4391  }
4392  if ((k%3)!=0)
4393    printf("\n");
4394  if (model->logLevel()<2)
4395    return ;
4396  int kMax = model->logLevel()>3 ? 1000000 : 10;
4397  k=0;
4398  for (iRow=1;iRow<=numberRows;iRow++) {
4399    if (number[iRow]) {
4400      k++;
4401      printf("%d columns have %d entries\n",number[iRow],iRow);
4402      if (k==kMax)
4403        break;
4404    }
4405  }
4406  if (k<numberRows) {
4407    int kk=k;
4408    k=0;
4409    for (iRow=numberRows;iRow>=1;iRow--) {
4410      if (number[iRow]) {
4411        k++;
4412        if (k==kMax)
4413          break;
4414      }
4415    }
4416    if (k>kk) {
4417      printf("\n    .........\n\n");
4418      iRow=k;
4419      k=0;
4420      for (;iRow<numberRows;iRow++) {
4421        if (number[iRow]) {
4422          k++;
4423          printf("%d columns have %d entries\n",number[iRow],iRow);
4424          if (k==kMax)
4425            break;
4426        }
4427      }
4428    }
4429  }
4430  delete [] number;
4431  printf("\n\n");
4432  // get row copy
4433  CoinPackedMatrix rowCopy = *matrix;
4434  rowCopy.reverseOrdering();
4435  //const int * column = rowCopy.getIndices();
4436  const int * rowLength = rowCopy.getVectorLengths();
4437  //const CoinBigIndex * rowStart = rowCopy.getVectorStarts();
4438  //const double * element = rowCopy.getElements();
4439  number = new int[numberColumns+1];
4440  memset(number,0,(numberColumns+1)*sizeof(int));
4441  for (iRow=0;iRow<numberRows;iRow++) {
4442    int length=rowLength[iRow];
4443    number[length]++;
4444  }
4445  if (number[0])
4446    printf("** %d rows have no entries\n",number[0]);
4447  k=0;
4448  for (iColumn=1;iColumn<=numberColumns;iColumn++) {
4449    if (number[iColumn]) {
4450      k++;
4451      printf("%d rows have %d entries\n",number[iColumn],iColumn);
4452      if (k==kMax)
4453        break;
4454    }
4455  }
4456  if (k<numberColumns) {
4457    int kk=k;
4458    k=0;
4459    for (iColumn=numberColumns;iColumn>=1;iColumn--) {
4460      if (number[iColumn]) {
4461        k++;
4462        if (k==kMax)
4463          break;
4464      }
4465    }
4466    if (k>kk) {
4467      printf("\n    .........\n\n");
4468      iColumn=k;
4469      k=0;
4470      for (;iColumn<numberColumns;iColumn++) {
4471        if (number[iColumn]) {
4472          k++;
4473          printf("%d rows have %d entries\n",number[iColumn],iColumn);
4474          if (k==kMax)
4475            break;
4476        }
4477      }
4478    }
4479  }
4480  delete [] number;
4481  // Now do breakdown of ranges
4482  breakdown("Elements",numberElements,elementByColumn);
4483  breakdown("RowLower",numberRows,rowLower);
4484  breakdown("RowUpper",numberRows,rowUpper);
4485  breakdown("ColumnLower",numberColumns,columnLower);
4486  breakdown("ColumnUpper",numberColumns,columnUpper);
4487  breakdown("Objective",numberColumns,objective);
4488}
4489static bool maskMatches(const int * starts, char ** masks,
4490                        std::string & check)
4491{
4492  // back to char as I am old fashioned
4493  const char * checkC = check.c_str();
4494  int length = strlen(checkC);
4495  while (checkC[length-1]==' ')
4496    length--;
4497  for (int i=starts[length];i<starts[length+1];i++) {
4498    char * thisMask = masks[i];
4499    int k;
4500    for ( k=0;k<length;k++) {
4501      if (thisMask[k]!='?'&&thisMask[k]!=checkC[k]) 
4502        break;
4503    }
4504    if (k==length)
4505      return true;
4506  }
4507  return false;
4508}
4509static void clean(char * temp)
4510{
4511  char * put = temp;
4512  while (*put>=' ')
4513    put++;
4514  *put='\0';
4515}
4516static void generateCode(CbcModel * model, const char * fileName,int type,int preProcess)
4517{
4518  // options on code generation
4519  bool sizecode = (type&4)!=0;
4520  type &= 3;
4521  FILE * fp = fopen(fileName,"r");
4522  assert (fp);
4523  int numberLines=0;
4524#define MAXLINES 5000
4525#define MAXONELINE 200
4526  char line[MAXLINES][MAXONELINE];
4527  strcpy(line[numberLines++],"0#if defined(_MSC_VER)");
4528  strcpy(line[numberLines++],"0// Turn off compiler warning about long names");
4529  strcpy(line[numberLines++],"0#  pragma warning(disable:4786)");
4530  strcpy(line[numberLines++],"0#endif\n");
4531  strcpy(line[numberLines++],"0#include <cassert>");
4532  strcpy(line[numberLines++],"0#include <iomanip>");
4533  strcpy(line[numberLines++],"0#include \"OsiClpSolverInterface.hpp\"");
4534  strcpy(line[numberLines++],"0#include \"CbcModel.hpp\"");
4535  strcpy(line[numberLines++],"0#include \"CbcCutGenerator.hpp\"");
4536  strcpy(line[numberLines++],"0#include \"CbcStrategy.hpp\"");
4537  strcpy(line[numberLines++],"0#include \"CglPreProcess.hpp\"");
4538  strcpy(line[numberLines++],"0#include \"CoinTime.hpp\"");
4539  if (preProcess>0) 
4540    strcpy(line[numberLines++],"0#include \"CglProbing.hpp\""); // possibly redundant
4541  // To allow generated 5's to be just before branchAndBound - do rest here
4542  strcpy(line[numberLines++],"5  cbcModel->initialSolve();");
4543  strcpy(line[numberLines++],"5  if (clpModel->tightenPrimalBounds()!=0) {");
4544  strcpy(line[numberLines++],"5    std::cout<<\"Problem is infeasible - tightenPrimalBounds!\"<<std::endl;");
4545  strcpy(line[numberLines++],"5    exit(1);");
4546  strcpy(line[numberLines++],"5  }");
4547  strcpy(line[numberLines++],"5  clpModel->dual();  // clean up");
4548  if (sizecode) {
4549    // override some settings
4550    strcpy(line[numberLines++],"5  // compute some things using problem size");
4551    strcpy(line[numberLines++],"5  cbcModel->setMinimumDrop(min(5.0e-2,");
4552    strcpy(line[numberLines++],"5       fabs(cbcModel->getMinimizationObjValue())*1.0e-3+1.0e-4));");
4553    strcpy(line[numberLines++],"5  if (cbcModel->getNumCols()<500)");
4554    strcpy(line[numberLines++],"5    cbcModel->setMaximumCutPassesAtRoot(-100); // always do 100 if possible");
4555    strcpy(line[numberLines++],"5  else if (cbcModel->getNumCols()<5000)");
4556    strcpy(line[numberLines++],"5    cbcModel->setMaximumCutPassesAtRoot(100); // use minimum drop");
4557    strcpy(line[numberLines++],"5  else");
4558    strcpy(line[numberLines++],"5    cbcModel->setMaximumCutPassesAtRoot(20);");
4559    strcpy(line[numberLines++],"5  cbcModel->setMaximumCutPasses(1);");
4560  }
4561  if (preProcess<=0) {
4562    // no preprocessing or strategy
4563    if (preProcess) {
4564      strcpy(line[numberLines++],"5  // Preprocessing using CbcStrategy");
4565      strcpy(line[numberLines++],"5  CbcStrategyDefault strategy(true,5,5);");
4566      strcpy(line[numberLines++],"5  strategy.setupPreProcessing(1);");
4567      strcpy(line[numberLines++],"5  cbcModel->setStrategy(strategy);");
4568    }
4569  } else {
4570    int translate[]={9999,0,0,-1,2,3,-2};
4571    strcpy(line[numberLines++],"5  // Hand coded preprocessing");
4572    strcpy(line[numberLines++],"5  CglPreProcess process;");
4573    strcpy(line[numberLines++],"5  OsiSolverInterface * saveSolver=cbcModel->solver()->clone();");
4574    strcpy(line[numberLines++],"5  // Tell solver we are in Branch and Cut");
4575    strcpy(line[numberLines++],"5  saveSolver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo) ;");
4576    strcpy(line[numberLines++],"5  // Default set of cut generators");
4577    strcpy(line[numberLines++],"5  CglProbing generator1;");
4578    strcpy(line[numberLines++],"5  generator1.setUsingObjective(true);");
4579    strcpy(line[numberLines++],"5  generator1.setMaxPass(3);");
4580    strcpy(line[numberLines++],"5  generator1.setMaxProbeRoot(saveSolver->getNumCols());");
4581    strcpy(line[numberLines++],"5  generator1.setMaxElements(100);");
4582    strcpy(line[numberLines++],"5  generator1.setMaxLookRoot(50);");
4583    strcpy(line[numberLines++],"5  generator1.setRowCuts(3);");
4584    strcpy(line[numberLines++],"5  // Add in generators");
4585    strcpy(line[numberLines++],"5  process.addCutGenerator(&generator1);");
4586    strcpy(line[numberLines++],"5  process.messageHandler()->setLogLevel(cbcModel->logLevel());");
4587    strcpy(line[numberLines++],"5  OsiSolverInterface * solver2 = ");
4588    sprintf(line[numberLines++],"5    process.preProcessNonDefault(*saveSolver,%d,10);",translate[preProcess]);
4589    strcpy(line[numberLines++],"5  // Tell solver we are not in Branch and Cut");
4590    strcpy(line[numberLines++],"5  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;");
4591    strcpy(line[numberLines++],"5  if (solver2)");
4592    strcpy(line[numberLines++],"5    solver2->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;");
4593    strcpy(line[numberLines++],"5  if (!solver2) {");
4594    strcpy(line[numberLines++],"5    std::cout<<\"Pre-processing says infeasible!\"<<std::endl;");
4595    strcpy(line[numberLines++],"5    exit(1);");
4596    strcpy(line[numberLines++],"5  } else {");
4597    strcpy(line[numberLines++],"5    std::cout<<\"processed model has \"<<solver2->getNumRows()");
4598    strcpy(line[numberLines++],"5            <<\" rows, \"<<solver2->getNumCols()");
4599    strcpy(line[numberLines++],"5            <<\" columns and \"<<solver2->getNumElements()");
4600    strcpy(line[numberLines++],"5            <<\" elements\"<<solver2->getNumElements()<<std::endl;");
4601    strcpy(line[numberLines++],"5  }");
4602    strcpy(line[numberLines++],"5  // we have to keep solver2 so pass clone");
4603    strcpy(line[numberLines++],"5  solver2 = solver2->clone();");
4604    strcpy(line[numberLines++],"5  cbcModel->assignSolver(solver2);");
4605    strcpy(line[numberLines++],"5  cbcModel->initialSolve();");
4606  }
4607  while (fgets(line[numberLines],MAXONELINE,fp)) {
4608    assert (numberLines<MAXLINES);
4609    clean(line[numberLines]);
4610    numberLines++;
4611  }
4612  fclose(fp);
4613  strcpy(line[numberLines++],"0\nint main (int argc, const char *argv[])\n{");
4614  strcpy(line[numberLines++],"0  OsiClpSolverInterface solver1;");
4615  strcpy(line[numberLines++],"0  int status=1;");
4616  strcpy(line[numberLines++],"0  if (argc<2)");
4617  strcpy(line[numberLines++],"0    std::cout<<\"Please give file name\"<<std::endl;");
4618  strcpy(line[numberLines++],"0  else");
4619  strcpy(line[numberLines++],"0    status=solver1.readMps(argv[1],\"\");");
4620  strcpy(line[numberLines++],"0  if (status) {");
4621  strcpy(line[numberLines++],"0    std::cout<<\"Bad readMps \"<<argv[1]<<std::endl;");
4622  strcpy(line[numberLines++],"0    exit(1);");
4623  strcpy(line[numberLines++],"0  }\n");
4624  strcpy(line[numberLines++],"0  double time1 = CoinCpuTime();");
4625  strcpy(line[numberLines++],"0  CbcModel model(solver1);");
4626  strcpy(line[numberLines++],"0  // Now do requested saves and modifications");
4627  strcpy(line[numberLines++],"0  CbcModel * cbcModel = & model;");
4628  strcpy(line[numberLines++],"0  OsiSolverInterface * osiModel = model.solver();");
4629  strcpy(line[numberLines++],"0  OsiClpSolverInterface * osiclpModel = dynamic_cast< OsiClpSolverInterface*> (osiModel);");
4630  strcpy(line[numberLines++],"0  ClpSimplex * clpModel = osiclpModel->getModelPtr();");
4631  // add in comments about messages
4632  strcpy(line[numberLines++],"3  // You can save some time by switching off message building");
4633  strcpy(line[numberLines++],"3  // clpModel->messagesPointer()->setDetailMessages(100,10000,(int *) NULL);");
4634  // add in actual solve
4635  strcpy(line[numberLines++],"5  cbcModel->branchAndBound();");
4636  strcpy(line[numberLines++],"8  std::cout<<argv[1]<<\" took \"<<CoinCpuTime()-time1<<\" seconds, \"");
4637  strcpy(line[numberLines++],"8    <<cbcModel->getNodeCount()<<\" nodes with objective \"");
4638  strcpy(line[numberLines++],"8    <<cbcModel->getObjValue()");
4639  strcpy(line[numberLines++],"8    <<(!cbcModel->status() ? \" Finished\" : \" Not finished\")");
4640  strcpy(line[numberLines++],"8    <<std::endl;");
4641  strcpy(line[numberLines++],"5  // For best solution");
4642  strcpy(line[numberLines++],"5  int numberColumns = solver1.getNumCols();");
4643  strcpy(line[numberLines++],"5  if (cbcModel->getMinimizationObjValue()<1.0e50) {");
4644  if (preProcess>0) {
4645    strcpy(line[numberLines++],"5    // post process");
4646    strcpy(line[numberLines++],"5    process.postProcess(*cbcModel->solver());");
4647    strcpy(line[numberLines++],"5    // Solution now back in saveSolver");
4648    strcpy(line[numberLines++],"5    cbcModel->assignSolver(saveSolver);");
4649    strcpy(line[numberLines++],"5    memcpy(cbcModel->bestSolution(),cbcModel->solver()->getColSolution(),");
4650    strcpy(line[numberLines++],"5          numberColumns*sizeof(double));");
4651  }
4652  strcpy(line[numberLines++],"5    // put back in original solver");
4653  strcpy(line[numberLines++],"5    solver1.setColSolution(cbcModel->bestSolution());");
4654  strcpy(line[numberLines++],"5    const double * solution = solver1.getColSolution();");
4655  strcpy(line[numberLines++],"8  \n  // Now you would use solution etc etc\n");
4656  strcpy(line[numberLines++],"5");
4657  strcpy(line[numberLines++],"5    // Get names from solver1 (as OsiSolverInterface may lose)");
4658  strcpy(line[numberLines++],"5    std::vector<std::string> columnNames = *solver1.getModelPtr()->columnNames();");
4659  strcpy(line[numberLines++],"5    ");
4660  strcpy(line[numberLines++],"5    int iColumn;");
4661  strcpy(line[numberLines++],"5    std::cout<<std::setiosflags(std::ios::fixed|std::ios::showpoint)<<std::setw(14);");
4662  strcpy(line[numberLines++],"5    ");
4663  strcpy(line[numberLines++],"5    std::cout<<\"--------------------------------------\"<<std::endl;");
4664  strcpy(line[numberLines++],"5    for (iColumn=0;iColumn<numberColumns;iColumn++) {");
4665  strcpy(line[numberLines++],"5      double value=solution[iColumn];");
4666  strcpy(line[numberLines++],"5      if (fabs(value)>1.0e-7&&solver1.isInteger(iColumn)) ");
4667  strcpy(line[numberLines++],"5 std::cout<<std::setw(6)<<iColumn<<\" \"");
4668  strcpy(line[numberLines++],"5                 <<columnNames[iColumn]<<\" \"");
4669  strcpy(line[numberLines++],"5                 <<value<<std::endl;");
4670  strcpy(line[numberLines++],"5    }");
4671  strcpy(line[numberLines++],"5    std::cout<<\"--------------------------------------\"<<std::endl;");
4672  strcpy(line[numberLines++],"5  ");
4673  strcpy(line[numberLines++],"5    std::cout<<std::resetiosflags(std::ios::fixed|std::ios::showpoint|std::ios::scientific);");
4674  strcpy(line[numberLines++],"5  }");
4675  strcpy(line[numberLines++],"8  return 0;\n}");
4676  fp = fopen(fileName,"w");
4677  assert (fp);
4678
4679  int wanted[9];
4680  memset(wanted,0,sizeof(wanted));
4681  wanted[0]=wanted[3]=wanted[5]=wanted[8]=1;
4682  if (type>0) 
4683    wanted[1]=wanted[6]=1;
4684  if (type>1) 
4685    wanted[2]=wanted[4]=wanted[7]=1;
4686  std::string header[9]=
4687  { "","Save values","Redundant save of default values","Set changed values",
4688    "Redundant set default values","Solve","Restore values","Redundant restore values","Finish up"};
4689  for (int iType=0;iType<9;iType++) {
4690    if (!wanted[iType])
4691      continue;
4692    int n=0;
4693    int iLine;
4694    for (iLine=0;iLine<numberLines;iLine++) {
4695      if (line[iLine][0]=='0'+iType) {
4696        if (!n&&header[iType]!="")
4697          fprintf(fp,"\n  // %s\n\n",header[iType].c_str());
4698        n++;
4699        // skip save and clp as cloned
4700        if (!strstr(line[iLine],"save")||(!strstr(line[iLine],"clpMo")&&
4701                                          !strstr(line[iLine],"_Osi")))
4702          fprintf(fp,"%s\n",line[iLine]+1);
4703      }
4704    }
4705  }
4706  fclose(fp);
4707  printf("C++ file written to %s\n",fileName);
4708}
4709/*
4710  Version 1.00.00 November 16 2005.
4711  This is to stop me (JJF) messing about too much.
4712  Tuning changes should be noted here.
4713  The testing next version may be activated by CBC_NEXT_VERSION
4714  This applies to OsiClp, Clp etc
4715  Version 1.00.01 November 24 2005
4716  Added several classes for advanced users.  This can't affect code (if you don't use it)
4717  Made some tiny changes (for N way branching) which should not change anything.
4718  CbcNWay object class - for N way branching this also allows use of CbcConsequence class.
4719  CbcBranchAllDifferent object class - for branching on general integer variables
4720  to stop them having same value so branches are x >= y+1 and x <= y-1.
4721  Added two new Cgl classes - CglAllDifferent which does column fixing (too slowly)
4722  and CglStored which just has a list of cuts which can be activated.
4723  Modified preprocess option to SOS
4724  Version 1.00.02 December 9 2005
4725  Added use of CbcStrategy to do clean preprocessing
4726  Added use of referenceSolver for cleaner repetition of Cbc
4727  Version 1.01.00 February 2 2006
4728  Added first try at Ampl interface
4729*/
Note: See TracBrowser for help on using the repository browser.