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

Last change on this file since 354 was 354, checked in by ladanyi, 14 years ago

finished switch to subversion

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