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

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

failed without ASL

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