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

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

towards common use with other solvers

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