source: trunk/Cbc/src/CoinSolve.cpp @ 364

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

for ampl preprocess sos

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