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

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

for local tree search and feasibility pump

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