source: trunk/Test/CoinSolve.cpp @ 237

Last change on this file since 237 was 237, checked in by forrest, 14 years ago

for ampl

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 100.2 KB
Line 
1// Copyright (C) 2004, International Business Machines
2// Corporation and others.  All Rights Reserved.
3   
4#include "CoinPragma.hpp"
5
6#include <cassert>
7#include <cstdio>
8#include <cmath>
9#include <cfloat>
10#include <string>
11#include <iostream>
12
13
14#include "CoinPragma.hpp"
15#include "CoinHelperFunctions.hpp"
16// Same version as CBC
17#define CBCVERSION "1.00.02"
18
19#include "CoinMpsIO.hpp"
20
21#include "ClpFactorization.hpp"
22#include "CoinTime.hpp"
23#include "ClpSimplex.hpp"
24#include "ClpSolve.hpp"
25#include "ClpPackedMatrix.hpp"
26#include "ClpPlusMinusOneMatrix.hpp"
27#include "ClpNetworkMatrix.hpp"
28#include "ClpDualRowSteepest.hpp"
29#include "ClpDualRowDantzig.hpp"
30#include "ClpLinearObjective.hpp"
31#include "ClpPrimalColumnSteepest.hpp"
32#include "ClpPrimalColumnDantzig.hpp"
33#include "ClpPresolve.hpp"
34#include "CbcOrClpParam.hpp"
35#include "OsiRowCutDebugger.hpp"
36#ifdef DMALLOC
37#include "dmalloc.h"
38#endif
39#ifdef WSSMP_BARRIER
40#define FOREIGN_BARRIER
41#endif
42#ifdef UFL_BARRIER
43#define FOREIGN_BARRIER
44#endif
45#ifdef TAUCS_BARRIER
46#define FOREIGN_BARRIER
47#endif
48#include "CoinWarmStartBasis.hpp"
49
50#include "OsiSolverInterface.hpp"
51#include "OsiCuts.hpp"
52#include "OsiRowCut.hpp"
53#include "OsiColCut.hpp"
54
55#include "CglPreProcess.hpp"
56#include "CglCutGenerator.hpp"
57#include "CglGomory.hpp"
58#include "CglProbing.hpp"
59#include "CglKnapsackCover.hpp"
60#include "CglRedSplit.hpp"
61#include "CglClique.hpp"
62#include "CglFlowCover.hpp"
63#include "CglMixedIntegerRounding2.hpp"
64#include "CglTwomir.hpp"
65
66#include "CbcModel.hpp"
67#include "CbcHeuristic.hpp"
68#include "CbcHeuristicLocal.hpp"
69#include "CbcHeuristicGreedy.hpp"
70#include "CbcHeuristicFPump.hpp"
71#include "CbcTreeLocal.hpp"
72#include "CbcCompareActual.hpp"
73#include "CbcBranchActual.hpp"
74#include  "CbcOrClpParam.hpp"
75#include  "CbcCutGenerator.hpp"
76
77#include "OsiClpSolverInterface.hpp"
78#ifdef CBC_AMPL
79#include "Cbc_ampl.h"
80static bool usingAmpl=false;
81#endif
82static double totalTime=0.0;
83
84//#############################################################################
85
86#ifdef NDEBUG
87#undef NDEBUG
88#endif
89// Allow for interrupts
90// But is this threadsafe ? (so switched off by option)
91
92#include "CoinSignal.hpp"
93static CbcModel * currentBranchModel = NULL;
94
95extern "C" {
96   static void signal_handler(int whichSignal)
97   {
98      if (currentBranchModel!=NULL) 
99         currentBranchModel->setMaximumNodes(0); // stop at next node
100      return;
101   }
102}
103
104int mainTest (int argc, const char *argv[],int algorithm,
105              ClpSimplex empty, bool doPresolve,int switchOff);
106int CbcOrClpRead_mode=1;
107FILE * CbcOrClpReadCommand=stdin;
108static bool noPrinting=false;
109static int * analyze(OsiClpSolverInterface * solverMod, int & numberChanged, double & increment,
110                     bool changeInt)
111{
112  OsiSolverInterface * solver = solverMod->clone();
113  if (0) {
114    // just get increment
115    CbcModel model(*solver);
116    model.analyzeObjective();
117    double increment2=model.getCutoffIncrement();
118    printf("initial cutoff increment %g\n",increment2);
119  }
120  const double *objective = solver->getObjCoefficients() ;
121  const double *lower = solver->getColLower() ;
122  const double *upper = solver->getColUpper() ;
123  int numberColumns = solver->getNumCols() ;
124  int numberRows = solver->getNumRows();
125  double direction = solver->getObjSense();
126  int iRow,iColumn;
127
128  // Row copy
129  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
130  const double * elementByRow = matrixByRow.getElements();
131  const int * column = matrixByRow.getIndices();
132  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
133  const int * rowLength = matrixByRow.getVectorLengths();
134
135  // Column copy
136  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
137  const double * element = matrixByCol.getElements();
138  const int * row = matrixByCol.getIndices();
139  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
140  const int * columnLength = matrixByCol.getVectorLengths();
141
142  const double * rowLower = solver->getRowLower();
143  const double * rowUpper = solver->getRowUpper();
144
145  char * ignore = new char [numberRows];
146  int * changed = new int[numberColumns];
147  int * which = new int[numberRows];
148  double * changeRhs = new double[numberRows];
149  memset(changeRhs,0,numberRows*sizeof(double));
150  memset(ignore,0,numberRows);
151  numberChanged=0;
152  int numberInteger=0;
153  for (iColumn=0;iColumn<numberColumns;iColumn++) {
154    if (upper[iColumn] > lower[iColumn]+1.0e-8&&solver->isInteger(iColumn)) 
155      numberInteger++;
156  }
157  bool finished=false;
158  while (!finished) {
159    int saveNumberChanged = numberChanged;
160    for (iRow=0;iRow<numberRows;iRow++) {
161      int numberContinuous=0;
162      double value1=0.0,value2=0.0;
163      bool allIntegerCoeff=true;
164      double sumFixed=0.0;
165      int jColumn1=-1,jColumn2=-1;
166      for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
167        int jColumn = column[j];
168        double value = elementByRow[j];
169        if (upper[jColumn] > lower[jColumn]+1.0e-8) {
170          if (!solver->isInteger(jColumn)) {
171            if (numberContinuous==0) {
172              jColumn1=jColumn;
173              value1=value;
174            } else {
175              jColumn2=jColumn;
176              value2=value;
177            }
178            numberContinuous++;
179          } else {
180            if (fabs(value-floor(value+0.5))>1.0e-12)
181              allIntegerCoeff=false;
182          }
183        } else {
184          sumFixed += lower[jColumn]*value;
185        }
186      }
187      double low = rowLower[iRow];
188      if (low>-1.0e20) {
189        low -= sumFixed;
190        if (fabs(low-floor(low+0.5))>1.0e-12)
191          allIntegerCoeff=false;
192      }
193      double up = rowUpper[iRow];
194      if (up<1.0e20) {
195        up -= sumFixed;
196        if (fabs(up-floor(up+0.5))>1.0e-12)
197          allIntegerCoeff=false;
198      }
199      if (numberContinuous==1) {
200        // see if really integer
201        // This does not allow for complicated cases
202        if (low==up) {
203          if (fabs(value1)>1.0e-3) {
204            value1 = 1.0/value1;
205            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
206              // integer
207              changed[numberChanged++]=jColumn1;
208              solver->setInteger(jColumn1);
209              if (upper[jColumn1]>1.0e20)
210                solver->setColUpper(jColumn1,1.0e20);
211              if (lower[jColumn1]<-1.0e20)
212                solver->setColLower(jColumn1,-1.0e20);
213            }
214          }
215        } else {
216          if (fabs(value1)>1.0e-3) {
217            value1 = 1.0/value1;
218            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
219              // This constraint will not stop it being integer
220              ignore[iRow]=1;
221            }
222          }
223        }
224      } else if (numberContinuous==2) {
225        if (low==up) {
226          /* need general theory - for now just look at 2 cases -
227             1 - +- 1 one in column and just costs i.e. matching objective
228             2 - +- 1 two in column but feeds into G/L row which will try and minimize
229          */
230          if (fabs(value1)==1.0&&value1*value2==-1.0&&!lower[jColumn1]
231              &&!lower[jColumn2]) {
232            int n=0;
233            int i;
234            double objChange=direction*(objective[jColumn1]+objective[jColumn2]);
235            double bound = CoinMin(upper[jColumn1],upper[jColumn2]);
236            bound = CoinMin(bound,1.0e20);
237            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
238              int jRow = row[i];
239              double value = element[i];
240              if (jRow!=iRow) {
241                which[n++]=jRow;
242                changeRhs[jRow]=value;
243              }
244            }
245            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
246              int jRow = row[i];
247              double value = element[i];
248              if (jRow!=iRow) {
249                if (!changeRhs[jRow]) {
250                  which[n++]=jRow;
251                  changeRhs[jRow]=value;
252                } else {
253                  changeRhs[jRow]+=value;
254                }
255              }
256            }
257            if (objChange>=0.0) {
258              // see if all rows OK
259              bool good=true;
260              for (i=0;i<n;i++) {
261                int jRow = which[i];
262                double value = changeRhs[jRow];
263                if (value) {
264                  value *= bound;
265                  if (rowLength[jRow]==1) {
266                    if (value>0.0) {
267                      double rhs = rowLower[jRow];
268                      if (rhs>0.0) {
269                        double ratio =rhs/value;
270                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
271                          good=false;
272                      }
273                    } else {
274                      double rhs = rowUpper[jRow];
275                      if (rhs<0.0) {
276                        double ratio =rhs/value;
277                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
278                          good=false;
279                      }
280                    }
281                  } else if (rowLength[jRow]==2) {
282                    if (value>0.0) {
283                      if (rowLower[jRow]>-1.0e20)
284                        good=false;
285                    } else {
286                      if (rowUpper[jRow]<1.0e20)
287                        good=false;
288                    }
289                  } else {
290                    good=false;
291                  }
292                }
293              }
294              if (good) {
295                // both can be integer
296                changed[numberChanged++]=jColumn1;
297                solver->setInteger(jColumn1);
298                if (upper[jColumn1]>1.0e20)
299                  solver->setColUpper(jColumn1,1.0e20);
300                if (lower[jColumn1]<-1.0e20)
301                  solver->setColLower(jColumn1,-1.0e20);
302                changed[numberChanged++]=jColumn2;
303                solver->setInteger(jColumn2);
304                if (upper[jColumn2]>1.0e20)
305                  solver->setColUpper(jColumn2,1.0e20);
306                if (lower[jColumn2]<-1.0e20)
307                  solver->setColLower(jColumn2,-1.0e20);
308              }
309            }
310            // clear
311            for (i=0;i<n;i++) {
312              changeRhs[which[i]]=0.0;
313            }
314          }
315        }
316      }
317    }
318    for (iColumn=0;iColumn<numberColumns;iColumn++) {
319      if (upper[iColumn] > lower[iColumn]+1.0e-8&&!solver->isInteger(iColumn)) {
320        double value;
321        value = upper[iColumn];
322        if (value<1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
323          continue;
324        value = lower[iColumn];
325        if (value>-1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
326          continue;
327        bool integer=true;
328        for (CoinBigIndex j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
329          int iRow = row[j];
330          if (!ignore[iRow]) {
331            integer=false;
332            break;
333          }
334        }
335        if (integer) {
336          // integer
337          changed[numberChanged++]=iColumn;
338          solver->setInteger(iColumn);
339          if (upper[iColumn]>1.0e20)
340            solver->setColUpper(iColumn,1.0e20);
341          if (lower[iColumn]<-1.0e20)
342            solver->setColLower(iColumn,-1.0e20);
343        }
344      }
345    }
346    finished = numberChanged==saveNumberChanged;
347  }
348  delete [] which;
349  delete [] changeRhs;
350  delete [] ignore;
351  if (numberInteger&&!noPrinting)
352    printf("%d integer variables",numberInteger);
353  if (changeInt) {
354    if (!noPrinting) {
355      if (numberChanged)
356        printf(" and %d variables made integer\n",numberChanged);
357      else
358        printf("\n");
359    }
360    delete [] ignore;
361    //increment=0.0;
362    if (!numberChanged) {
363      delete [] changed;
364      delete solver;
365      return NULL;
366    } else {
367      for (iColumn=0;iColumn<numberColumns;iColumn++) {
368        if (solver->isInteger(iColumn))
369          solverMod->setInteger(iColumn);
370      }
371      delete solver;
372      return changed;
373    }
374  } else {
375    if (!noPrinting) {
376      if (numberChanged)
377        printf(" and %d variables could be made integer\n",numberChanged);
378      else
379        printf("\n");
380    }
381    // just get increment
382    CbcModel model(*solver);
383    if (noPrinting)
384      model.setLogLevel(0);
385    model.analyzeObjective();
386    double increment2=model.getCutoffIncrement();
387    if (increment2>increment) {
388      if (!noPrinting)
389        printf("cutoff increment increased from %g to %g\n",increment,increment2);
390      increment=increment2;
391    }
392    delete solver;
393    numberChanged=0;
394    delete [] changed;
395    return NULL;
396  }
397}
398int main (int argc, const char *argv[])
399{
400  /* Note
401     This is meant as a stand-alone executable to do as much of coin as possible.
402     It should only have one solver known to it.
403  */
404  {
405    double time1 = CoinCpuTime(),time2;
406    bool goodModel=false;
407    CoinSighandler_t saveSignal=SIG_DFL;
408    // register signal handler
409    saveSignal = signal(SIGINT,signal_handler);
410    // Set up all non-standard stuff
411    OsiClpSolverInterface solver1;
412    CbcModel model(solver1);
413    CbcModel * babModel = NULL;
414    model.setNumberBeforeTrust(21);
415    int cutPass=-1234567;
416    OsiSolverInterface * solver = model.solver();
417    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
418    ClpSimplex * lpSolver = clpSolver->getModelPtr();
419    clpSolver->messageHandler()->setLogLevel(0) ;
420    model.messageHandler()->setLogLevel(1);
421#ifdef CBC_AMPL
422    ampl_info info;
423    if (argc>2&&!strcmp(argv[2],"-AMPL")) {
424      usingAmpl=true;
425      int returnCode = readAmpl(&info,argc,const_cast<char **>(argv));
426      if (returnCode)
427        return returnCode;
428      CbcOrClpRead_mode=2; // so will start with parameters
429      // see if log in list
430      noPrinting=true;
431      for (int i=1;i<info.numberArguments;i++) {
432        if (!strcmp(info.arguments[i],"log")) {
433          if (i<info.numberArguments-1&&atoi(info.arguments[i+1])>0)
434            noPrinting=false;
435          break;
436        }
437      }
438      if (noPrinting) {
439        model.messageHandler()->setLogLevel(0);
440        setCbcOrClpPrinting(false);
441      }
442      if (!noPrinting)
443        printf("%d rows, %d columns and %d elements\n",
444               info.numberRows,info.numberColumns,info.numberElements);
445      solver->loadProblem(info.numberColumns,info.numberRows,info.starts,
446                          info.rows,info.elements,
447                          info.columnLower,info.columnUpper,info.objective,
448                          info.rowLower,info.rowUpper);
449      // If we had a solution use it
450      if (info.primalSolution) {
451        solver->setColSolution(info.primalSolution);
452      }
453      // status
454      if (info.rowStatus) {
455        unsigned char * statusArray = lpSolver->statusArray();
456        int i;
457        for (i=0;i<info.numberColumns;i++)
458          statusArray[i]=(char)info.columnStatus[i];
459        statusArray+=info.numberColumns;
460        for (i=0;i<info.numberRows;i++)
461          statusArray[i]=(char)info.rowStatus[i];
462        CoinWarmStartBasis * basis = lpSolver->getBasis();
463        solver->setWarmStart(basis);
464        delete basis;
465      }
466      freeArrays1(&info);
467      // modify objective if necessary
468      solver->setObjSense(info.direction);
469      solver->setDblParam(OsiObjOffset,info.offset);
470      // Set integer variables
471      for (int i=info.numberColumns-info.numberBinary-info.numberIntegers;
472           i<info.numberColumns;i++)
473        solver->setInteger(i);
474      goodModel=true;
475      // change argc etc
476      argc = info.numberArguments;
477      argv = const_cast<const char **>(info.arguments);
478    }
479#endif   
480    // default action on import
481    int allowImportErrors=0;
482    int keepImportNames=1;
483    int doIdiot=-1;
484    int outputFormat=2;
485    int slpValue=-1;
486    int printOptions=0;
487    int printMode=0;
488    int presolveOptions=0;
489    int doCrash=0;
490    int doSprint=-1;
491    int doScaling=1;
492    // set reasonable defaults
493    int preSolve=5;
494    int preProcess=4;
495    bool preSolveFile=false;
496   
497    double djFix=1.0e100;
498    double gapRatio=1.0e100;
499    double tightenFactor=0.0;
500    lpSolver->setPerturbation(50);
501    lpSolver->messageHandler()->setPrefix(false);
502    const char dirsep =  CoinFindDirSeparator();
503    std::string directory = (dirsep == '/' ? "./" : ".\\");
504    std::string defaultDirectory = directory;
505    std::string importFile ="";
506    std::string exportFile ="default.mps";
507    std::string importBasisFile ="";
508    std::string debugFile="";
509    double * debugValues = NULL;
510    int numberDebugValues = -1;
511    int basisHasValues=0;
512    std::string exportBasisFile ="default.bas";
513    std::string saveFile ="default.prob";
514    std::string restoreFile ="default.prob";
515    std::string solutionFile ="stdout";
516    std::string solutionSaveFile ="solution.file";
517#define CBCMAXPARAMETERS 200
518    CbcOrClpParam parameters[CBCMAXPARAMETERS];
519    int numberParameters ;
520    establishParams(numberParameters,parameters) ;
521    parameters[whichParam(BASISIN,numberParameters,parameters)].setStringValue(importBasisFile);
522    parameters[whichParam(BASISOUT,numberParameters,parameters)].setStringValue(exportBasisFile);
523    parameters[whichParam(DEBUG,numberParameters,parameters)].setStringValue(debugFile);
524    parameters[whichParam(DIRECTORY,numberParameters,parameters)].setStringValue(directory);
525    parameters[whichParam(DUALBOUND,numberParameters,parameters)].setDoubleValue(lpSolver->dualBound());
526    parameters[whichParam(DUALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->dualTolerance());
527    parameters[whichParam(EXPORT,numberParameters,parameters)].setStringValue(exportFile);
528    parameters[whichParam(IDIOT,numberParameters,parameters)].setIntValue(doIdiot);
529    parameters[whichParam(IMPORT,numberParameters,parameters)].setStringValue(importFile);
530    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].setDoubleValue(1.0e-8);
531    int slog = whichParam(SOLVERLOGLEVEL,numberParameters,parameters);
532    int log = whichParam(LOGLEVEL,numberParameters,parameters);
533    parameters[slog].setIntValue(0);
534    parameters[log].setIntValue(1);
535    parameters[whichParam(MAXFACTOR,numberParameters,parameters)].setIntValue(lpSolver->factorizationFrequency());
536    parameters[whichParam(MAXITERATION,numberParameters,parameters)].setIntValue(lpSolver->maximumIterations());
537    parameters[whichParam(OUTPUTFORMAT,numberParameters,parameters)].setIntValue(outputFormat);
538    parameters[whichParam(PRESOLVEPASS,numberParameters,parameters)].setIntValue(preSolve);
539    parameters[whichParam(PERTVALUE,numberParameters,parameters)].setIntValue(lpSolver->perturbation());
540    parameters[whichParam(PRIMALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->primalTolerance());
541    parameters[whichParam(PRIMALWEIGHT,numberParameters,parameters)].setDoubleValue(lpSolver->infeasibilityCost());
542    parameters[whichParam(RESTORE,numberParameters,parameters)].setStringValue(restoreFile);
543    parameters[whichParam(SAVE,numberParameters,parameters)].setStringValue(saveFile);
544    //parameters[whichParam(TIMELIMIT,numberParameters,parameters)].setDoubleValue(1.0e8);
545    parameters[whichParam(TIMELIMIT_BAB,numberParameters,parameters)].setDoubleValue(1.0e8);
546    parameters[whichParam(SOLUTION,numberParameters,parameters)].setStringValue(solutionFile);
547    parameters[whichParam(SAVESOL,numberParameters,parameters)].setStringValue(solutionSaveFile);
548    parameters[whichParam(SPRINT,numberParameters,parameters)].setIntValue(doSprint);
549    model.setNumberBeforeTrust(5);
550    parameters[whichParam(NUMBERBEFORE,numberParameters,parameters)].setIntValue(5);
551    parameters[whichParam(MAXNODES,numberParameters,parameters)].setIntValue(model.getMaximumNodes());
552    model.setNumberStrong(5);
553    parameters[whichParam(STRONGBRANCHING,numberParameters,parameters)].setIntValue(model.numberStrong());
554    parameters[whichParam(INFEASIBILITYWEIGHT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcInfeasibilityWeight));
555    parameters[whichParam(INTEGERTOLERANCE,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcIntegerTolerance));
556    parameters[whichParam(INCREMENT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcCutoffIncrement));
557    // Set up likely cut generators and defaults
558    parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("sos");
559    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(128|64|1);
560    parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].setIntValue(-1);
561    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].setIntValue(100);
562
563    CglGomory gomoryGen;
564    // try larger limit
565    gomoryGen.setLimitAtRoot(512);
566    gomoryGen.setLimit(50);
567    // set default action (0=off,1=on,2=root)
568    int gomoryAction=3;
569    parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
570
571    CglProbing probingGen;
572    probingGen.setUsingObjective(true);
573    probingGen.setMaxPass(3);
574    probingGen.setMaxPassRoot(3);
575    // Number of unsatisfied variables to look at
576    probingGen.setMaxProbe(10);
577    probingGen.setMaxProbeRoot(50);
578    // How far to follow the consequences
579    probingGen.setMaxLook(10);
580    probingGen.setMaxLookRoot(50);
581    probingGen.setMaxLookRoot(10);
582    // Only look at rows with fewer than this number of elements
583    probingGen.setMaxElements(200);
584    probingGen.setRowCuts(3);
585    // set default action (0=off,1=on,2=root)
586    int probingAction=3;
587    parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
588
589    CglKnapsackCover knapsackGen;
590    //knapsackGen.switchOnExpensive();
591    // set default action (0=off,1=on,2=root)
592    int knapsackAction=3;
593    parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
594
595    CglRedSplit redsplitGen;
596    //redsplitGen.setLimit(100);
597    // set default action (0=off,1=on,2=root)
598    // Off as seems to give some bad cuts
599    int redsplitAction=2;
600    parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption("root");
601
602    CglClique cliqueGen(false,true);
603    cliqueGen.setStarCliqueReport(false);
604    cliqueGen.setRowCliqueReport(false);
605    cliqueGen.setMinViolation(0.1);
606    // set default action (0=off,1=on,2=root)
607    int cliqueAction=3;
608    parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption("ifmove");
609
610    CglMixedIntegerRounding2 mixedGen;
611    // set default action (0=off,1=on,2=root)
612    int mixedAction=3;
613    parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
614
615    CglFlowCover flowGen;
616    // set default action (0=off,1=on,2=root)
617    int flowAction=3;
618    parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
619
620    CglTwomir twomirGen;
621    twomirGen.setMaxElements(250);
622    // set default action (0=off,1=on,2=root)
623    int twomirAction=2;
624    parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption("root");
625
626    bool useRounding=true;
627    parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption("on");
628    bool useFpump=true;
629    parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption("on");
630    bool useGreedy=true;
631    parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption("on");
632    bool useCombine=true;
633    parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
634    bool useLocalTree=false;
635    parameters[whichParam(COSTSTRATEGY,numberParameters,parameters)].setCurrentOption("off");
636    bool useCosts=false;
637   
638    // total number of commands read
639    int numberGoodCommands=0;
640    // Set false if user does anything advanced
641    bool defaultSettings=true;
642
643    // Hidden stuff for barrier
644    int choleskyType = 0;
645    int gamma=0;
646    int scaleBarrier=0;
647    int doKKT=0;
648    int crossover=2; // do crossover unless quadratic
649    // For names
650    int lengthName = 0;
651    std::vector<std::string> rowNames;
652    std::vector<std::string> columnNames;
653   
654    std::string field;
655    if (!noPrinting) {
656      std::cout<<"Coin Cbc and Clp Solver version "<<CBCVERSION
657               <<", build "<<__DATE__<<std::endl;
658      // Print command line
659      if (argc>1) {
660        printf("command line - ");
661        for (int i=0;i<argc;i++)
662          printf("%s ",argv[i]);
663        printf("\n");
664      }
665    }
666    while (1) {
667      // next command
668      field=CoinReadGetCommand(argc,argv);
669      // exit if null or similar
670      if (!field.length()) {
671        if (numberGoodCommands==1&&goodModel) {
672          // we just had file name - do branch and bound
673          field="branch";
674        } else if (!numberGoodCommands) {
675          // let's give the sucker a hint
676          std::cout
677            <<"CoinSolver takes input from arguments ( - switches to stdin)"
678            <<std::endl
679            <<"Enter ? for list of commands or help"<<std::endl;
680          field="-";
681        } else {
682          break;
683        }
684      }
685     
686      // see if ? at end
687      int numberQuery=0;
688      if (field!="?"&&field!="???") {
689        int length = field.length();
690        int i;
691        for (i=length-1;i>0;i--) {
692          if (field[i]=='?') 
693            numberQuery++;
694          else
695            break;
696        }
697        field=field.substr(0,length-numberQuery);
698      }
699      // find out if valid command
700      int iParam;
701      int numberMatches=0;
702      int firstMatch=-1;
703      for ( iParam=0; iParam<numberParameters; iParam++ ) {
704        int match = parameters[iParam].matches(field);
705        if (match==1) {
706          numberMatches = 1;
707          firstMatch=iParam;
708          break;
709        } else {
710          if (match&&firstMatch<0)
711            firstMatch=iParam;
712          numberMatches += match>>1;
713        }
714      }
715      if (iParam<numberParameters&&!numberQuery) {
716        // found
717        CbcOrClpParam found = parameters[iParam];
718        CbcOrClpParameterType type = found.type();
719        int valid;
720        numberGoodCommands++;
721        if (type==BAB&&goodModel) {
722          // check if any integers
723          if (!lpSolver->integerInformation())
724            type=DUALSIMPLEX;
725        }
726        if (type==GENERALQUERY) {
727          std::cout<<"In argument list keywords have leading - "
728            ", -stdin or just - switches to stdin"<<std::endl;
729          std::cout<<"One command per line (and no -)"<<std::endl;
730          std::cout<<"abcd? gives list of possibilities, if only one + explanation"<<std::endl;
731          std::cout<<"abcd?? adds explanation, if only one fuller help"<<std::endl;
732          std::cout<<"abcd without value (where expected) gives current value"<<std::endl;
733          std::cout<<"abcd value sets value"<<std::endl;
734          std::cout<<"Commands are:"<<std::endl;
735          int maxAcross=5;
736          int limits[]={1,51,101,151,201,251,301,351,401};
737          std::vector<std::string> types;
738          types.push_back("Double parameters:");
739          types.push_back("Branch and Cut double parameters:");
740          types.push_back("Integer parameters:");
741          types.push_back("Branch and Cut integer parameters:");
742          types.push_back("Keyword parameters and others:");
743          types.push_back("Branch and Cut keyword parameters and others:");
744          types.push_back("Actions:");
745          types.push_back("Branch and Cut actions:");
746          int iType;
747          for (iType=0;iType<8;iType++) {
748            int across=0;
749            std::cout<<types[iType]<<"   ";
750            for ( iParam=0; iParam<numberParameters; iParam++ ) {
751              int type = parameters[iParam].type();
752              if (parameters[iParam].displayThis()&&type>=limits[iType]
753                  &&type<limits[iType+1]) {
754                if (!across)
755                  std::cout<<"  ";
756                std::cout<<parameters[iParam].matchName()<<"  ";
757                across++;
758                if (across==maxAcross) {
759                  std::cout<<std::endl;
760                  across=0;
761                }
762              }
763            }
764            if (across)
765              std::cout<<std::endl;
766          }
767        } else if (type==FULLGENERALQUERY) {
768          std::cout<<"Full list of commands is:"<<std::endl;
769          int maxAcross=5;
770          int limits[]={1,51,101,151,201,251,301,351,401};
771          std::vector<std::string> types;
772          types.push_back("Double parameters:");
773          types.push_back("Branch and Cut double parameters:");
774          types.push_back("Integer parameters:");
775          types.push_back("Branch and Cut integer parameters:");
776          types.push_back("Keyword parameters and others:");
777          types.push_back("Branch and Cut keyword parameters and others:");
778          types.push_back("Actions:");
779          types.push_back("Branch and Cut actions:");
780          int iType;
781          for (iType=0;iType<8;iType++) {
782            int across=0;
783            std::cout<<types[iType]<<"  ";
784            for ( iParam=0; iParam<numberParameters; iParam++ ) {
785              int type = parameters[iParam].type();
786              if (type>=limits[iType]
787                  &&type<limits[iType+1]) {
788                if (!across)
789                  std::cout<<"  ";
790                std::cout<<parameters[iParam].matchName()<<"  ";
791                across++;
792                if (across==maxAcross) {
793                  std::cout<<std::endl;
794                  across=0;
795                }
796              }
797            }
798            if (across)
799              std::cout<<std::endl;
800          }
801        } else if (type<101) {
802          // get next field as double
803          double value = CoinReadGetDoubleField(argc,argv,&valid);
804          if (!valid) {
805            if (type<51) {
806              parameters[iParam].setDoubleParameter(lpSolver,value);
807            } else if (type<81) {
808              parameters[iParam].setDoubleParameter(model,value);
809            } else {
810              parameters[iParam].setDoubleParameter(lpSolver,value);
811              switch(type) {
812              case DJFIX:
813                djFix=value;
814                preSolve=5;
815                defaultSettings=false; // user knows what she is doing
816                break;
817              case GAPRATIO:
818                gapRatio=value;
819                break;
820              case TIGHTENFACTOR:
821                tightenFactor=value;
822                defaultSettings=false; // user knows what she is doing
823                break;
824              default:
825                abort();
826              }
827            }
828          } else if (valid==1) {
829            abort();
830          } else {
831            std::cout<<parameters[iParam].name()<<" has value "<<
832              parameters[iParam].doubleValue()<<std::endl;
833          }
834        } else if (type<201) {
835          // get next field as int
836          int value = CoinReadGetIntField(argc,argv,&valid);
837          if (!valid) {
838            if (type<151) {
839              if (parameters[iParam].type()==PRESOLVEPASS)
840                preSolve = value;
841              else if (parameters[iParam].type()==IDIOT)
842                doIdiot = value;
843              else if (parameters[iParam].type()==SPRINT)
844                doSprint = value;
845              else if (parameters[iParam].type()==OUTPUTFORMAT)
846                outputFormat = value;
847              else if (parameters[iParam].type()==SLPVALUE)
848                slpValue = value;
849              else if (parameters[iParam].type()==PRESOLVEOPTIONS)
850                presolveOptions = value;
851              else if (parameters[iParam].type()==PRINTOPTIONS)
852                printOptions = value;
853              else if (parameters[iParam].type()==CUTPASS)
854                cutPass = value;
855              else if (parameters[iParam].type()==FPUMPITS)
856                { useFpump = true;parameters[iParam].setIntValue(value);}
857              parameters[iParam].setIntParameter(lpSolver,value);
858            } else {
859              parameters[iParam].setIntParameter(model,value);
860            }
861          } else if (valid==1) {
862            abort();
863          } else {
864            std::cout<<parameters[iParam].name()<<" has value "<<
865              parameters[iParam].intValue()<<std::endl;
866          }
867        } else if (type<301) {
868          // one of several strings
869          std::string value = CoinReadGetString(argc,argv);
870          int action = parameters[iParam].parameterOption(value);
871          if (action<0) {
872            if (value!="EOL") {
873              // no match
874              parameters[iParam].printOptions();
875            } else {
876              // print current value
877              std::cout<<parameters[iParam].name()<<" has value "<<
878                parameters[iParam].currentOption()<<std::endl;
879            }
880          } else {
881            parameters[iParam].setCurrentOption(action,!noPrinting);
882            // for now hard wired
883            switch (type) {
884            case DIRECTION:
885              if (action==0)
886                lpSolver->setOptimizationDirection(1);
887              else if (action==1)
888                lpSolver->setOptimizationDirection(-1);
889              else
890                lpSolver->setOptimizationDirection(0);
891              break;
892            case DUALPIVOT:
893              if (action==0) {
894                ClpDualRowSteepest steep(3);
895                lpSolver->setDualRowPivotAlgorithm(steep);
896              } else if (action==1) {
897                //ClpDualRowDantzig dantzig;
898                ClpDualRowSteepest dantzig(5);
899                lpSolver->setDualRowPivotAlgorithm(dantzig);
900              } else if (action==2) {
901                // partial steep
902                ClpDualRowSteepest steep(2);
903                lpSolver->setDualRowPivotAlgorithm(steep);
904              } else {
905                ClpDualRowSteepest steep;
906                lpSolver->setDualRowPivotAlgorithm(steep);
907              }
908              break;
909            case PRIMALPIVOT:
910              if (action==0) {
911                ClpPrimalColumnSteepest steep(3);
912                lpSolver->setPrimalColumnPivotAlgorithm(steep);
913              } else if (action==1) {
914                ClpPrimalColumnSteepest steep(0);
915                lpSolver->setPrimalColumnPivotAlgorithm(steep);
916              } else if (action==2) {
917                ClpPrimalColumnDantzig dantzig;
918                lpSolver->setPrimalColumnPivotAlgorithm(dantzig);
919              } else if (action==3) {
920                ClpPrimalColumnSteepest steep(2);
921                lpSolver->setPrimalColumnPivotAlgorithm(steep);
922              } else if (action==4) {
923                ClpPrimalColumnSteepest steep(1);
924                lpSolver->setPrimalColumnPivotAlgorithm(steep);
925              } else if (action==5) {
926                ClpPrimalColumnSteepest steep(4);
927                lpSolver->setPrimalColumnPivotAlgorithm(steep);
928              } else if (action==6) {
929                ClpPrimalColumnSteepest steep(10);
930                lpSolver->setPrimalColumnPivotAlgorithm(steep);
931              }
932              break;
933            case SCALING:
934              lpSolver->scaling(action);
935              solver->setHintParam(OsiDoScale,action!=0,OsiHintTry);
936              doScaling = 1-action;
937              break;
938            case AUTOSCALE:
939              lpSolver->setAutomaticScaling(action!=0);
940              break;
941            case SPARSEFACTOR:
942              lpSolver->setSparseFactorization((1-action)!=0);
943              break;
944            case BIASLU:
945              lpSolver->factorization()->setBiasLU(action);
946              break;
947            case PERTURBATION:
948              if (action==0)
949                lpSolver->setPerturbation(50);
950              else
951                lpSolver->setPerturbation(100);
952              break;
953            case ERRORSALLOWED:
954              allowImportErrors = action;
955              break;
956            case INTPRINT:
957              printMode=action;
958              break;
959              //case ALGORITHM:
960              //algorithm  = action;
961              //defaultSettings=false; // user knows what she is doing
962              //abort();
963              //break;
964            case KEEPNAMES:
965              keepImportNames = 1-action;
966              break;
967            case PRESOLVE:
968              if (action==0)
969                preSolve = 5;
970              else if (action==1)
971                preSolve=0;
972              else if (action==2)
973                preSolve=10;
974              else
975                preSolveFile=true;
976              break;
977            case PFI:
978              lpSolver->factorization()->setForrestTomlin(action==0);
979              break;
980            case CRASH:
981              doCrash=action;
982              break;
983            case MESSAGES:
984              lpSolver->messageHandler()->setPrefix(action!=0);
985              break;
986            case CHOLESKY:
987              choleskyType = action;
988              break;
989            case GAMMA:
990              gamma=action;
991              break;
992            case BARRIERSCALE:
993              scaleBarrier=action;
994              break;
995            case KKT:
996              doKKT=action;
997              break;
998            case CROSSOVER:
999              crossover=action;
1000              break;
1001            case GOMORYCUTS:
1002              defaultSettings=false; // user knows what she is doing
1003              gomoryAction = action;
1004              break;
1005            case PROBINGCUTS:
1006              defaultSettings=false; // user knows what she is doing
1007              probingAction = action;
1008              break;
1009            case KNAPSACKCUTS:
1010              defaultSettings=false; // user knows what she is doing
1011              knapsackAction = action;
1012              break;
1013            case REDSPLITCUTS:
1014              defaultSettings=false; // user knows what she is doing
1015              redsplitAction = action;
1016              break;
1017            case CLIQUECUTS:
1018              defaultSettings=false; // user knows what she is doing
1019              cliqueAction = action;
1020              break;
1021            case FLOWCUTS:
1022              defaultSettings=false; // user knows what she is doing
1023              flowAction = action;
1024              break;
1025            case MIXEDCUTS:
1026              defaultSettings=false; // user knows what she is doing
1027              mixedAction = action;
1028              break;
1029            case TWOMIRCUTS:
1030              defaultSettings=false; // user knows what she is doing
1031              twomirAction = action;
1032              break;
1033            case ROUNDING:
1034              defaultSettings=false; // user knows what she is doing
1035              useRounding = action;
1036              break;
1037            case FPUMP:
1038              defaultSettings=false; // user knows what she is doing
1039              useFpump=action;
1040              break;
1041            case CUTSSTRATEGY:
1042              gomoryAction = action;
1043              probingAction = action;
1044              knapsackAction = action;
1045              cliqueAction = action;
1046              flowAction = action;
1047              mixedAction = action;
1048              twomirAction = action;
1049              parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption(action);
1050              parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption(action);
1051              parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption(action);
1052              if (!action) {
1053                redsplitAction = action;
1054                parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption(action);
1055              }
1056              parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption(action);
1057              parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption(action);
1058              parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption(action);
1059              parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption(action);
1060              break;
1061            case HEURISTICSTRATEGY:
1062              useRounding = action;
1063              useGreedy = action;
1064              useCombine = action;
1065              //useLocalTree = action;
1066              useFpump=action;
1067              parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption(action);
1068              parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption(action);
1069              parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption(action);
1070              //parameters[whichParam(LOCALTREE,numberParameters,parameters)].setCurrentOption(action);
1071              parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption(action);
1072              break;
1073            case GREEDY:
1074              defaultSettings=false; // user knows what she is doing
1075              useGreedy = action;
1076              break;
1077            case COMBINE:
1078              defaultSettings=false; // user knows what she is doing
1079              useCombine = action;
1080              break;
1081            case LOCALTREE:
1082              defaultSettings=false; // user knows what she is doing
1083              useLocalTree = action;
1084              break;
1085            case COSTSTRATEGY:
1086              if (action!=1&&action!=0) {
1087                printf("Pseudo costs not implemented yet\n");
1088              } else {
1089                useCosts=action;
1090              }
1091              break;
1092            case PREPROCESS:
1093              preProcess = action;
1094              break;
1095            default:
1096              abort();
1097            }
1098          }
1099        } else {
1100          // action
1101          if (type==EXIT) {
1102#ifdef CBC_AMPL
1103            if(usingAmpl) {
1104              if (info.numberIntegers||info.numberBinary) {
1105                // integer
1106              } else {
1107                // linear
1108              }
1109              writeAmpl(&info);
1110              freeArrays2(&info);
1111              freeArgs(&info);
1112            }
1113#endif
1114            break; // stop all
1115          }
1116          switch (type) {
1117          case DUALSIMPLEX:
1118          case PRIMALSIMPLEX:
1119          case SOLVECONTINUOUS:
1120          case BARRIER:
1121            if (goodModel) {
1122              ClpSolve::SolveType method;
1123              ClpSolve::PresolveType presolveType;
1124              ClpSimplex * model2 = lpSolver;
1125              if (noPrinting)
1126                lpSolver->setLogLevel(0);
1127              ClpSolve solveOptions;
1128              if (preSolve!=5&&preSolve)
1129                presolveType=ClpSolve::presolveNumber;
1130              else if (preSolve)
1131                presolveType=ClpSolve::presolveOn;
1132              else
1133                presolveType=ClpSolve::presolveOff;
1134              solveOptions.setPresolveType(presolveType,preSolve);
1135              if (type==DUALSIMPLEX||SOLVECONTINUOUS) {
1136                method=ClpSolve::useDual;
1137              } else if (type==PRIMALSIMPLEX) {
1138                method=ClpSolve::usePrimalorSprint;
1139              } else {
1140                method = ClpSolve::useBarrier;
1141                if (crossover==1) {
1142                  method=ClpSolve::useBarrierNoCross;
1143                } else if (crossover==2) {
1144                  ClpObjective * obj = lpSolver->objectiveAsObject();
1145                  if (obj->type()>1) {
1146                    method=ClpSolve::useBarrierNoCross;
1147                    presolveType=ClpSolve::presolveOff;
1148                    solveOptions.setPresolveType(presolveType,preSolve);
1149                  } 
1150                }
1151              }
1152              solveOptions.setSolveType(method);
1153              if(preSolveFile)
1154                presolveOptions |= 0x40000000;
1155              solveOptions.setSpecialOption(4,presolveOptions);
1156              solveOptions.setSpecialOption(5,printOptions);
1157              if (method==ClpSolve::useDual) {
1158                // dual
1159                if (doCrash)
1160                  solveOptions.setSpecialOption(0,1,doCrash); // crash
1161                else if (doIdiot)
1162                  solveOptions.setSpecialOption(0,2,doIdiot); // possible idiot
1163              } else if (method==ClpSolve::usePrimalorSprint) {
1164                // primal
1165                // if slp turn everything off
1166                if (slpValue>0) {
1167                  doCrash=false;
1168                  doSprint=0;
1169                  doIdiot=-1;
1170                  solveOptions.setSpecialOption(1,10,slpValue); // slp
1171                  method=ClpSolve::usePrimal;
1172                }
1173                if (doCrash) {
1174                  solveOptions.setSpecialOption(1,1,doCrash); // crash
1175                } else if (doSprint>0) {
1176                  // sprint overrides idiot
1177                  solveOptions.setSpecialOption(1,3,doSprint); // sprint
1178                } else if (doIdiot>0) {
1179                  solveOptions.setSpecialOption(1,2,doIdiot); // idiot
1180                } else if (slpValue<=0) {
1181                  if (doIdiot==0) {
1182                    if (doSprint==0)
1183                      solveOptions.setSpecialOption(1,4); // all slack
1184                    else
1185                      solveOptions.setSpecialOption(1,9); // all slack or sprint
1186                  } else {
1187                    if (doSprint==0)
1188                      solveOptions.setSpecialOption(1,8); // all slack or idiot
1189                    else
1190                      solveOptions.setSpecialOption(1,7); // initiative
1191                  }
1192                }
1193                if (basisHasValues==-1)
1194                  solveOptions.setSpecialOption(1,11); // switch off values
1195              } else if (method==ClpSolve::useBarrier||method==ClpSolve::useBarrierNoCross) {
1196                int barrierOptions = choleskyType;
1197                if (scaleBarrier)
1198                  barrierOptions |= 8;
1199                if (gamma)
1200                  barrierOptions |= 16;
1201                if (doKKT)
1202                  barrierOptions |= 32;
1203                solveOptions.setSpecialOption(4,barrierOptions);
1204              }
1205              model2->initialSolve(solveOptions);
1206              basisHasValues=1;
1207#ifdef CBC_AMPL
1208              if (usingAmpl) {
1209                double value = model2->getObjValue()*model2->getObjSense();
1210                char buf[300];
1211                int pos=0;
1212                int iStat = model2->status();
1213                if (iStat==0) {
1214                  pos += sprintf(buf+pos,"optimal," );
1215                } else if (iStat==1) {
1216                  // infeasible
1217                  pos += sprintf(buf+pos,"infeasible,");
1218                } else if (iStat==2) {
1219                  // unbounded
1220                  pos += sprintf(buf+pos,"unbounded,");
1221                } else if (iStat==3) {
1222                  pos += sprintf(buf+pos,"stopped on iterations or time,");
1223                } else if (iStat==4) {
1224                  iStat = 7;
1225                  pos += sprintf(buf+pos,"stopped on difficulties,");
1226                } else if (iStat==5) {
1227                  iStat = 3;
1228                  pos += sprintf(buf+pos,"stopped on ctrl-c,");
1229                } else {
1230                  pos += sprintf(buf+pos,"status unknown,");
1231                  iStat=6;
1232                }
1233                info.problemStatus=iStat;
1234                info.objValue = value;
1235                pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
1236                               value);
1237                sprintf(buf+pos,"\n%d iterations",
1238                        model2->getIterationCount());
1239                free(info.primalSolution);
1240                int numberColumns=model2->numberColumns();
1241                info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
1242                CoinCopyN(model2->primalColumnSolution(),numberColumns,info.primalSolution);
1243                int numberRows = model2->numberRows();
1244                free(info.dualSolution);
1245                info.dualSolution = (double *) malloc(numberRows*sizeof(double));
1246                CoinCopyN(model2->dualRowSolution(),numberRows,info.dualSolution);
1247                CoinWarmStartBasis * basis = model2->getBasis();
1248                free(info.rowStatus);
1249                info.rowStatus = (int *) malloc(numberRows*sizeof(int));
1250                free(info.columnStatus);
1251                info.columnStatus = (int *) malloc(numberColumns*sizeof(int));
1252                // Put basis in
1253                int i;
1254                // free,basic,ub,lb are 0,1,2,3
1255                for (i=0;i<numberRows;i++) {
1256                  CoinWarmStartBasis::Status status = basis->getArtifStatus(i);
1257                  info.rowStatus[i]=status;
1258                }
1259                for (i=0;i<numberColumns;i++) {
1260                  CoinWarmStartBasis::Status status = basis->getStructStatus(i);
1261                  info.columnStatus[i]=status;
1262                }
1263                // put buffer into info
1264                strcpy(info.buffer,buf);
1265                delete basis;
1266              }
1267#endif
1268            } else {
1269              std::cout<<"** Current model not valid"<<std::endl;
1270            }
1271            break;
1272          case TIGHTEN:
1273            if (goodModel) {
1274              int numberInfeasibilities = lpSolver->tightenPrimalBounds();
1275              if (numberInfeasibilities)
1276                std::cout<<"** Analysis indicates model infeasible"<<std::endl;
1277            } else {
1278              std::cout<<"** Current model not valid"<<std::endl;
1279            }
1280            break;
1281          case PLUSMINUS:
1282            if (goodModel) {
1283              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
1284              ClpPackedMatrix* clpMatrix =
1285                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
1286              if (clpMatrix) {
1287                ClpPlusMinusOneMatrix * newMatrix = new ClpPlusMinusOneMatrix(*(clpMatrix->matrix()));
1288                if (newMatrix->getIndices()) {
1289                  lpSolver->replaceMatrix(newMatrix);
1290                  delete saveMatrix;
1291                  std::cout<<"Matrix converted to +- one matrix"<<std::endl;
1292                } else {
1293                  std::cout<<"Matrix can not be converted to +- 1 matrix"<<std::endl;
1294                }
1295              } else {
1296                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
1297              }
1298            } else {
1299              std::cout<<"** Current model not valid"<<std::endl;
1300            }
1301            break;
1302          case NETWORK:
1303            if (goodModel) {
1304              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
1305              ClpPackedMatrix* clpMatrix =
1306                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
1307              if (clpMatrix) {
1308                ClpNetworkMatrix * newMatrix = new ClpNetworkMatrix(*(clpMatrix->matrix()));
1309                if (newMatrix->getIndices()) {
1310                  lpSolver->replaceMatrix(newMatrix);
1311                  delete saveMatrix;
1312                  std::cout<<"Matrix converted to network matrix"<<std::endl;
1313                } else {
1314                  std::cout<<"Matrix can not be converted to network matrix"<<std::endl;
1315                }
1316              } else {
1317                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
1318              }
1319            } else {
1320              std::cout<<"** Current model not valid"<<std::endl;
1321            }
1322            break;
1323/*
1324  Run branch-and-cut. First set a few options -- node comparison, scaling. If
1325  the solver is Clp, consider running some presolve code (not yet converted
1326  this to generic OSI) with branch-and-cut. If presolve is disabled, or the
1327  solver is not Clp, simply run branch-and-cut. Print elapsed time at the end.
1328*/
1329          case BAB: // branchAndBound
1330          case STRENGTHEN:
1331            if (goodModel) {
1332              int logLevel = parameters[slog].intValue();
1333              // Reduce printout
1334              if (logLevel<=1)
1335                model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
1336              // Don't switch off all output
1337              {
1338                OsiSolverInterface * solver = model.solver();
1339                OsiClpSolverInterface * si =
1340                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
1341                assert (si != NULL);
1342                si->setSpecialOptions(0x40000000);
1343              }
1344              model.initialSolve();
1345              // If user made settings then use them
1346              if (!defaultSettings) {
1347                OsiSolverInterface * solver = model.solver();
1348                if (!doScaling)
1349                  solver->setHintParam(OsiDoScale,false,OsiHintTry);
1350                OsiClpSolverInterface * si =
1351                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
1352                assert (si != NULL);
1353                // get clp itself
1354                ClpSimplex * modelC = si->getModelPtr();
1355                //if (modelC->tightenPrimalBounds()!=0) {
1356                //std::cout<<"Problem is infeasible!"<<std::endl;
1357                //break;
1358                //}
1359                // bounds based on continuous
1360                if (tightenFactor) {
1361                  if (modelC->tightenPrimalBounds(tightenFactor)!=0) {
1362                    std::cout<<"Problem is infeasible!"<<std::endl;
1363                    break;
1364                  }
1365                }
1366                if (djFix<1.0e20) {
1367                  // do some fixing
1368                  int numberColumns = modelC->numberColumns();
1369                  int i;
1370                  const char * type = modelC->integerInformation();
1371                  double * lower = modelC->columnLower();
1372                  double * upper = modelC->columnUpper();
1373                  double * solution = modelC->primalColumnSolution();
1374                  double * dj = modelC->dualColumnSolution();
1375                  int numberFixed=0;
1376                  for (i=0;i<numberColumns;i++) {
1377                    if (type[i]) {
1378                      double value = solution[i];
1379                      if (value<lower[i]+1.0e-5&&dj[i]>djFix) {
1380                        solution[i]=lower[i];
1381                        upper[i]=lower[i];
1382                        numberFixed++;
1383                      } else if (value>upper[i]-1.0e-5&&dj[i]<-djFix) {
1384                        solution[i]=upper[i];
1385                        lower[i]=upper[i];
1386                        numberFixed++;
1387                      }
1388                    }
1389                  }
1390                  printf("%d columns fixed\n",numberFixed);
1391                }
1392              }
1393              // See if we want preprocessing
1394              OsiSolverInterface * saveSolver=NULL;
1395              CglPreProcess process;
1396              delete babModel;
1397              babModel = new CbcModel(model);
1398              OsiSolverInterface * solver3 = clpSolver->clone();
1399              babModel->assignSolver(solver3);
1400              OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
1401              int numberChanged=0;
1402              if (clpSolver2->messageHandler()->logLevel())
1403                clpSolver2->messageHandler()->setLogLevel(1);
1404              if (logLevel>-1)
1405                clpSolver2->messageHandler()->setLogLevel(logLevel);
1406              lpSolver = clpSolver2->getModelPtr();
1407              if (lpSolver->factorizationFrequency()==200) {
1408                // User did not touch preset
1409                int numberRows = lpSolver->numberRows();
1410                const int cutoff1=10000;
1411                const int cutoff2=100000;
1412                const int base=75;
1413                const int freq0 = 50;
1414                const int freq1=200;
1415                const int freq2=400;
1416                const int maximum=1000;
1417                int frequency;
1418                if (numberRows<cutoff1)
1419                  frequency=base+numberRows/freq0;
1420                else if (numberRows<cutoff2)
1421                  frequency=base+cutoff1/freq0 + (numberRows-cutoff1)/freq1;
1422                else
1423                  frequency=base+cutoff1/freq0 + (cutoff2-cutoff1)/freq1 + (numberRows-cutoff2)/freq2;
1424                lpSolver->setFactorizationFrequency(CoinMin(maximum,frequency));
1425              }
1426              time2 = CoinCpuTime();
1427              totalTime += time2-time1;
1428              time1 = time2;
1429              double timeLeft = babModel->getMaximumSeconds();
1430              if (preProcess&&type==BAB) {
1431                saveSolver=babModel->solver()->clone();
1432                /* Do not try and produce equality cliques and
1433                   do up to 10 passes */
1434                OsiSolverInterface * solver2;
1435                {
1436                  // Tell solver we are in Branch and Cut
1437                  saveSolver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo) ;
1438                  // Default set of cut generators
1439                  CglProbing generator1;
1440                  generator1.setUsingObjective(true);
1441                  generator1.setMaxPass(3);
1442                  generator1.setMaxProbeRoot(saveSolver->getNumCols());
1443                  generator1.setMaxElements(100);
1444                  generator1.setMaxLookRoot(50);
1445                  generator1.setRowCuts(3);
1446                  // Add in generators
1447                  process.addCutGenerator(&generator1);
1448                  int translate[]={9999,0,0,-1,2};
1449                  process.messageHandler()->setLogLevel(babModel->logLevel());
1450                  solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],10);
1451                  // Tell solver we are not in Branch and Cut
1452                  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
1453                  if (solver2)
1454                    solver2->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
1455                }
1456                if (!noPrinting) {
1457                  if (!solver2) {
1458                    printf("Pre-processing says infeasible\n");
1459                    break;
1460                  } else {
1461                    printf("processed model has %d rows, %d columns and %d elements\n",
1462                           solver2->getNumRows(),solver2->getNumCols(),solver2->getNumElements());
1463                  }
1464                }
1465                //solver2->resolve();
1466                if (preProcess==2) {
1467                  OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (solver2);
1468                  ClpSimplex * lpSolver = clpSolver2->getModelPtr();
1469                  lpSolver->writeMps("presolved.mps",0,1,lpSolver->optimizationDirection());
1470                  printf("Preprocessed model (minimization) on presolved.mps\n");
1471                }
1472                // we have to keep solver2 so pass clone
1473                solver2 = solver2->clone();
1474                babModel->assignSolver(solver2);
1475                babModel->initialSolve();
1476                babModel->setMaximumSeconds(timeLeft-(CoinCpuTime()-time1));
1477              }
1478              // now tighten bounds
1479              {
1480                OsiClpSolverInterface * si =
1481                  dynamic_cast<OsiClpSolverInterface *>(babModel->solver()) ;
1482                assert (si != NULL);
1483                // get clp itself
1484                ClpSimplex * modelC = si->getModelPtr();
1485                if (noPrinting)
1486                  modelC->setLogLevel(0);
1487                if (modelC->tightenPrimalBounds()!=0) {
1488                  std::cout<<"Problem is infeasible!"<<std::endl;
1489                  break;
1490                }
1491                modelC->dual();
1492              }
1493              if (debugValues) {
1494                // for debug
1495                std::string problemName ;
1496                babModel->solver()->getStrParam(OsiProbName,problemName) ;
1497                //babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
1498                twomirGen.probname_=strdup(problemName.c_str());
1499                // checking seems odd
1500                //redsplitGen.set_given_optsol(babModel->solver()->getRowCutDebuggerAlways()->optimalSolution(),
1501                //                         babModel->getNumCols());
1502              }
1503              if (useCosts) {
1504                int numberColumns = babModel->getNumCols();
1505                int * sort = new int[numberColumns];
1506                double * dsort = new double[numberColumns];
1507                int * priority = new int [numberColumns];
1508                const double * objective = babModel->getObjCoefficients();
1509                int iColumn;
1510                int n=0;
1511                for (iColumn=0;iColumn<numberColumns;iColumn++) {
1512                  if (babModel->isInteger(iColumn)) {
1513                    sort[n]=n;
1514                    dsort[n++]=-objective[iColumn];
1515                  }
1516                }
1517                CoinSort_2(dsort,dsort+n,sort);
1518                int level=0;
1519                double last = -1.0e100;
1520                for (int i=0;i<n;i++) {
1521                  int iPut=sort[i];
1522                  if (dsort[i]!=last) {
1523                    level++;
1524                    last=dsort[i];
1525                  }
1526                  priority[iPut]=level;
1527                }
1528                babModel->passInPriorities( priority,false);
1529                delete [] priority;
1530                delete [] sort;
1531                delete [] dsort;
1532              }
1533              // FPump done first as it only works if no solution
1534              CbcHeuristicFPump heuristic4(*babModel);
1535              if (useFpump) {
1536                heuristic4.setMaximumPasses(parameters[whichParam(FPUMPITS,numberParameters,parameters)].intValue());
1537                babModel->addHeuristic(&heuristic4);
1538              }
1539              CbcRounding heuristic1(*babModel);
1540              if (useRounding)
1541                babModel->addHeuristic(&heuristic1) ;
1542              CbcHeuristicLocal heuristic2(*babModel);
1543              heuristic2.setSearchType(1);
1544              if (useCombine)
1545                babModel->addHeuristic(&heuristic2);
1546              CbcHeuristicGreedyCover heuristic3(*babModel);
1547              CbcHeuristicGreedyEquality heuristic3a(*babModel);
1548              if (useGreedy) {
1549                babModel->addHeuristic(&heuristic3);
1550                babModel->addHeuristic(&heuristic3a);
1551              }
1552              if (useLocalTree) {
1553                CbcTreeLocal localTree(babModel,NULL,10,0,0,10000,2000);
1554                babModel->passInTreeHandler(localTree);
1555              }
1556              // add cut generators if wanted
1557              int switches[20];
1558              int numberGenerators=0;
1559              if (probingAction==1) {
1560                babModel->addCutGenerator(&probingGen,-1,"Probing");
1561                switches[numberGenerators++]=0;
1562              } else if (probingAction>=2) {
1563                babModel->addCutGenerator(&probingGen,-101+probingAction,"Probing");
1564                switches[numberGenerators++]=0;
1565              }
1566              if (gomoryAction==1) {
1567                babModel->addCutGenerator(&gomoryGen,-1,"Gomory");
1568                switches[numberGenerators++]=1;
1569              } else if (gomoryAction>=2) {
1570                babModel->addCutGenerator(&gomoryGen,-101+gomoryAction,"Gomory");
1571                switches[numberGenerators++]=-1;
1572              }
1573              if (knapsackAction==1) {
1574                babModel->addCutGenerator(&knapsackGen,-1,"Knapsack");
1575                switches[numberGenerators++]=0;
1576              } else if (knapsackAction>=2) {
1577                babModel->addCutGenerator(&knapsackGen,-101+knapsackAction,"Knapsack");
1578                switches[numberGenerators++]=0;
1579              }
1580              if (redsplitAction==1) {
1581                babModel->addCutGenerator(&redsplitGen,-1,"Reduce-and-split");
1582                switches[numberGenerators++]=1;
1583              } else if (redsplitAction>=2) {
1584                babModel->addCutGenerator(&redsplitGen,-101+redsplitAction,"Reduce-and-split");
1585                switches[numberGenerators++]=1;
1586              }
1587              if (cliqueAction==1) {
1588                babModel->addCutGenerator(&cliqueGen,-1,"Clique");
1589                switches[numberGenerators++]=1;
1590              } else if (cliqueAction>=2) {
1591                babModel->addCutGenerator(&cliqueGen,-101+cliqueAction,"Clique");
1592                switches[numberGenerators++]=-1;
1593              }
1594              if (mixedAction==1) {
1595                babModel->addCutGenerator(&mixedGen,-1,"MixedIntegerRounding2");
1596                switches[numberGenerators++]=1;
1597              } else if (mixedAction>=2) {
1598                babModel->addCutGenerator(&mixedGen,-101+mixedAction,"MixedIntegerRounding2");
1599                switches[numberGenerators++]=-1;
1600              }
1601              if (flowAction==1) {
1602                babModel->addCutGenerator(&flowGen,-1,"FlowCover");
1603                switches[numberGenerators++]=1;
1604              } else if (flowAction>=2) {
1605                babModel->addCutGenerator(&flowGen,-101+flowAction,"FlowCover");
1606                switches[numberGenerators++]=1;
1607              }
1608              if (twomirAction==1) {
1609                babModel->addCutGenerator(&twomirGen,-1,"TwoMirCuts");
1610                switches[numberGenerators++]=1;
1611              } else if (twomirAction>=2) {
1612                babModel->addCutGenerator(&twomirGen,-101+twomirAction,"TwoMirCuts");
1613                switches[numberGenerators++]=1;
1614              }
1615              // Say we want timings
1616              numberGenerators = babModel->numberCutGenerators();
1617              int iGenerator;
1618              int cutDepth=
1619                parameters[whichParam(CUTDEPTH,numberParameters,parameters)].intValue();
1620              for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
1621                CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
1622                int howOften = generator->howOften();
1623                if (howOften==-98||howOften==-99) 
1624                  generator->setSwitchOffIfLessThan(switches[iGenerator]);
1625                generator->setTiming(true);
1626                if (cutDepth>=0)
1627                  generator->setWhatDepth(cutDepth) ;
1628              }
1629              // Could tune more
1630              babModel->setMinimumDrop(min(5.0e-2,
1631                                        fabs(babModel->getMinimizationObjValue())*1.0e-3+1.0e-4));
1632              if (cutPass==-1234567) {
1633                if (babModel->getNumCols()<500)
1634                  babModel->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
1635                else if (babModel->getNumCols()<5000)
1636                  babModel->setMaximumCutPassesAtRoot(100); // use minimum drop
1637                else
1638                  babModel->setMaximumCutPassesAtRoot(20);
1639              } else {
1640                babModel->setMaximumCutPassesAtRoot(cutPass);
1641              }
1642              babModel->setMaximumCutPasses(1);
1643             
1644              // Do more strong branching if small
1645              //if (babModel->getNumCols()<5000)
1646              //babModel->setNumberStrong(20);
1647              // Switch off strong branching if wanted
1648              //if (babModel->getNumCols()>10*babModel->getNumRows())
1649              //babModel->setNumberStrong(0);
1650              if (!noPrinting) {
1651                babModel->messageHandler()->setLogLevel(parameters[log].intValue());
1652                if (babModel->getNumCols()>2000||babModel->getNumRows()>1500||
1653                    babModel->messageHandler()->logLevel()>1)
1654                  babModel->setPrintFrequency(100);
1655              }
1656             
1657              babModel->solver()->setIntParam(OsiMaxNumIterationHotStart,
1658                    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].intValue());
1659              OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
1660              // go faster stripes
1661              if (osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
1662                osiclp->setupForRepeatedUse(2,parameters[slog].intValue());
1663              } else {
1664                osiclp->setupForRepeatedUse(0,parameters[slog].intValue());
1665              }
1666              double increment=babModel->getCutoffIncrement();;
1667              int * changed = analyze( osiclp,numberChanged,increment,false);
1668              if (debugValues) {
1669                if (numberDebugValues==babModel->getNumCols()) {
1670                  // for debug
1671                  babModel->solver()->activateRowCutDebugger(debugValues) ;
1672                } else {
1673                  printf("debug file has incorrect number of columns\n");
1674                }
1675              }
1676              babModel->setCutoffIncrement(CoinMax(babModel->getCutoffIncrement(),increment));
1677              // Turn this off if you get problems
1678              // Used to be automatically set
1679              int mipOptions = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue();
1680              if (mipOptions!=(128|64|1))
1681                printf("mip options %d\n",mipOptions);
1682              osiclp->setSpecialOptions(mipOptions);
1683              if (gapRatio < 1.0e100) {
1684                double value = babModel->solver()->getObjValue() ;
1685                double value2 = gapRatio*(1.0e-5+fabs(value)) ;
1686                babModel->setAllowableGap(value2) ;
1687                std::cout << "Continuous " << value
1688                          << ", so allowable gap set to "
1689                          << value2 << std::endl ;
1690              }
1691              // probably faster to use a basis to get integer solutions
1692              babModel->setSpecialOptions(2);
1693              currentBranchModel = babModel;
1694              OsiSolverInterface * strengthenedModel=NULL;
1695              if (type==BAB) {
1696                int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
1697                if (moreMipOptions>=0) {
1698                  printf("more mip options %d\n",moreMipOptions);
1699                  babModel->setSearchStrategy(moreMipOptions);
1700                }
1701                if (preProcess&&process.numberSOS()) {
1702                  int numberSOS = process.numberSOS();
1703                  int numberIntegers = babModel->numberIntegers();
1704                  /* model may not have created objects
1705                     If none then create
1706                  */
1707                  if (!numberIntegers||!babModel->numberObjects()) {
1708                    babModel->findIntegers(true);
1709                    numberIntegers = babModel->numberIntegers();
1710                  }
1711                  CbcObject ** oldObjects = babModel->objects();
1712                  // Do sets and priorities
1713                  CbcObject ** objects = new CbcObject * [numberSOS];
1714                  // set old objects to have low priority
1715                  int numberOldObjects = babModel->numberObjects();
1716                  int numberColumns = babModel->getNumCols();
1717                  for (int iObj = 0;iObj<numberOldObjects;iObj++)
1718                    oldObjects[iObj]->setPriority(numberColumns+1);
1719                  const int * starts = process.startSOS();
1720                  const int * which = process.whichSOS();
1721                  const int * type = process.typeSOS();
1722                  const double * weight = process.weightSOS();
1723                  int iSOS;
1724                  for (iSOS =0;iSOS<numberSOS;iSOS++) {
1725                    int iStart = starts[iSOS];
1726                    int n=starts[iSOS+1]-iStart;
1727                    objects[iSOS] = new CbcSOS(babModel,n,which+iStart,weight+iStart,
1728                                               iSOS,type[iSOS]);
1729                    // branch on long sets first
1730                    objects[iSOS]->setPriority(numberColumns-n);
1731                  }
1732                  babModel->addObjects(numberSOS,objects);
1733                  for (iSOS=0;iSOS<numberSOS;iSOS++)
1734                    delete objects[iSOS];
1735                  delete [] objects;
1736                }
1737                int statistics = (printOptions>0) ? printOptions: 0;
1738                babModel->branchAndBound(statistics);
1739              } else {
1740                strengthenedModel = babModel->strengthenedModel();
1741              }
1742              currentBranchModel = NULL;
1743              osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
1744              if (debugFile=="createAfterPre"&&babModel->bestSolution()) {
1745                lpSolver = osiclp->getModelPtr();
1746                //move best solution (should be there -- but ..)
1747                int n = lpSolver->getNumCols();
1748                memcpy(lpSolver->primalColumnSolution(),babModel->bestSolution(),n*sizeof(double));
1749                saveSolution(osiclp->getModelPtr(),"debug.file");
1750              }
1751              if (!noPrinting) {
1752                // Print more statistics
1753                std::cout<<"Cuts at root node changed objective from "<<babModel->getContinuousObjective()
1754                         <<" to "<<babModel->rootObjectiveAfterCuts()<<std::endl;
1755               
1756                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
1757                  CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
1758                  std::cout<<generator->cutGeneratorName()<<" was tried "
1759                           <<generator->numberTimesEntered()<<" times and created "
1760                           <<generator->numberCutsInTotal()<<" cuts of which "
1761                           <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
1762                  if (generator->timing())
1763                    std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
1764                  else
1765                    std::cout<<std::endl;
1766                }
1767              }
1768              time2 = CoinCpuTime();
1769              totalTime += time2-time1;
1770              // For best solution
1771              double * bestSolution = NULL;
1772              if (babModel->getMinimizationObjValue()<1.0e50&&type==BAB) {
1773                // post process
1774                if (preProcess) {
1775                  int n = saveSolver->getNumCols();
1776                  bestSolution = new double [n];
1777                  process.postProcess(*babModel->solver());
1778                  // Solution now back in saveSolver
1779                  babModel->assignSolver(saveSolver);
1780                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
1781                } else {
1782                  int n = babModel->solver()->getNumCols();
1783                  bestSolution = new double [n];
1784                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
1785                }
1786              }
1787              if (type==STRENGTHEN&&strengthenedModel)
1788                clpSolver = dynamic_cast< OsiClpSolverInterface*> (strengthenedModel);
1789              lpSolver = clpSolver->getModelPtr();
1790              if (debugFile=="create"&&bestSolution) {
1791                saveSolution(lpSolver,"debug.file");
1792              }
1793              if (numberChanged) {
1794                for (int i=0;i<numberChanged;i++) {
1795                  int iColumn=changed[i];
1796                  clpSolver->setContinuous(iColumn);
1797                }
1798                delete [] changed;
1799              }
1800              if (type==BAB) {
1801                //move best solution (should be there -- but ..)
1802                int n = lpSolver->getNumCols();
1803                if (bestSolution)
1804                  memcpy(lpSolver->primalColumnSolution(),bestSolution,n*sizeof(double));
1805                delete [] bestSolution;
1806                std::string statusName[]={"Finished","Stopped on ","Difficulties",
1807                                          "","","User ctrl-c"};
1808                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
1809                int iStat = babModel->status();
1810                int iStat2 = babModel->secondaryStatus();
1811                if (!noPrinting)
1812                  std::cout<<"Result - "<<statusName[iStat]<<minor[iStat2]
1813                           <<" objective "<<babModel->getObjValue()<<
1814                    " after "<<babModel->getNodeCount()<<" nodes and "
1815                           <<babModel->getIterationCount()<<
1816                    " iterations - took "<<time2-time1<<" seconds"<<std::endl;
1817#ifdef CBC_AMPL
1818                if (usingAmpl) {
1819                  double value = babModel->getObjValue()*lpSolver->getObjSense();
1820                  char buf[300];
1821                  int pos=0;
1822                  if (iStat==0) {
1823                    if (babModel->getObjValue()<1.0e40) {
1824                      pos += sprintf(buf+pos,"optimal," );
1825                    } else {
1826                      // infeasible
1827                      iStat=1;
1828                      pos += sprintf(buf+pos,"infeasible,");
1829                    }
1830                  } else if (iStat==1) {
1831                    if (iStat2!=6)
1832                      iStat=3;
1833                    else
1834                      iStat=4;
1835                    pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
1836                  } else if (iStat==2) {
1837                    iStat = 7;
1838                    pos += sprintf(buf+pos,"stopped on difficulties,");
1839                  } else if (iStat==5) {
1840                    iStat = 3;
1841                    pos += sprintf(buf+pos,"stopped on ctrl-c,");
1842                  } else {
1843                    pos += sprintf(buf+pos,"status unknown,");
1844                    iStat=6;
1845                  }
1846                  info.problemStatus=iStat;
1847                  info.objValue = value;
1848                  if (babModel->getObjValue()<1.0e40) 
1849                    pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
1850                                   value);
1851                  sprintf(buf+pos,"\n%d nodes, %d iterations",
1852                          babModel->getNodeCount(),
1853                          babModel->getIterationCount());
1854                  if (bestSolution) {
1855                    free(info.primalSolution);
1856                    info.primalSolution = (double *) malloc(n*sizeof(double));
1857                    CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
1858                    int numberRows = lpSolver->numberRows();
1859                    free(info.dualSolution);
1860                    info.dualSolution = (double *) malloc(numberRows*sizeof(double));
1861                    CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
1862                  } else {
1863                    info.primalSolution=NULL;
1864                    info.dualSolution=NULL;
1865                  }
1866                  // put buffer into info
1867                  strcpy(info.buffer,buf);
1868                }
1869#endif
1870              } else {
1871                std::cout<<"Model strengthend - now has "<<clpSolver->getNumRows()
1872                         <<" rows"<<std::endl;
1873              }
1874              time1 = time2;
1875              delete babModel;
1876              babModel=NULL;
1877            } else {
1878              std::cout << "** Current model not valid" << std::endl ; 
1879            }
1880            break ;
1881          case IMPORT:
1882            {
1883              delete babModel;
1884              babModel=NULL;
1885              // get next field
1886              field = CoinReadGetString(argc,argv);
1887              if (field=="$") {
1888                field = parameters[iParam].stringValue();
1889              } else if (field=="EOL") {
1890                parameters[iParam].printString();
1891                break;
1892              } else {
1893                parameters[iParam].setStringValue(field);
1894              }
1895              std::string fileName;
1896              bool canOpen=false;
1897              if (field=="-") {
1898                // stdin
1899                canOpen=true;
1900                fileName = "-";
1901              } else {
1902                bool absolutePath;
1903                if (dirsep=='/') {
1904                  // non Windows (or cygwin)
1905                  absolutePath=(field[0]=='/');
1906                } else {
1907                  //Windows (non cycgwin)
1908                  absolutePath=(field[0]=='\\');
1909                  // but allow for :
1910                  if (strchr(field.c_str(),':'))
1911                    absolutePath=true;
1912                }
1913                if (absolutePath) {
1914                  fileName = field;
1915                } else if (field[0]=='~') {
1916                  char * environVar = getenv("HOME");
1917                  if (environVar) {
1918                    std::string home(environVar);
1919                    field=field.erase(0,1);
1920                    fileName = home+field;
1921                  } else {
1922                    fileName=field;
1923                  }
1924                } else {
1925                  fileName = directory+field;
1926                }
1927                FILE *fp=fopen(fileName.c_str(),"r");
1928                if (fp) {
1929                  // can open - lets go for it
1930                  fclose(fp);
1931                  canOpen=true;
1932                } else {
1933                  std::cout<<"Unable to open file "<<fileName<<std::endl;
1934                }
1935              }
1936              if (canOpen) {
1937                int status =lpSolver->readMps(fileName.c_str(),
1938                                                   keepImportNames!=0,
1939                                                   allowImportErrors!=0);
1940                if (!status||(status>0&&allowImportErrors)) {
1941                  if (keepImportNames) {
1942                    lengthName = lpSolver->lengthNames();
1943                    rowNames = *(lpSolver->rowNames());
1944                    columnNames = *(lpSolver->columnNames());
1945                  } else {
1946                    lengthName=0;
1947                  }
1948                  goodModel=true;
1949                  //Set integers in clpsolver
1950                  const char * info = lpSolver->integerInformation();
1951                  if (info) {
1952                    int numberColumns = lpSolver->numberColumns();
1953                    int i;
1954                    for (i=0;i<numberColumns;i++) {
1955                      if (info[i]) 
1956                        clpSolver->setInteger(i);
1957                    }
1958                  }
1959                  // sets to all slack (not necessary?)
1960                  lpSolver->createStatus();
1961                  time2 = CoinCpuTime();
1962                  totalTime += time2-time1;
1963                  time1=time2;
1964                  // Go to canned file if just input file
1965                  if (CbcOrClpRead_mode==2&&argc==2) {
1966                    // only if ends .mps
1967                    char * find = strstr(fileName.c_str(),".mps");
1968                    if (find&&find[4]=='\0') {
1969                      find[1]='p'; find[2]='a';find[3]='r';
1970                      FILE *fp=fopen(fileName.c_str(),"r");
1971                      if (fp) {
1972                        CbcOrClpReadCommand=fp; // Read from that file
1973                        CbcOrClpRead_mode=-1;
1974                      }
1975                    }
1976                  }
1977                } else {
1978                  // errors
1979                  std::cout<<"There were "<<status<<
1980                    " errors on input"<<std::endl;
1981                }
1982              }
1983            }
1984            break;
1985          case EXPORT:
1986            if (goodModel) {
1987              // get next field
1988              field = CoinReadGetString(argc,argv);
1989              if (field=="$") {
1990                field = parameters[iParam].stringValue();
1991              } else if (field=="EOL") {
1992                parameters[iParam].printString();
1993                break;
1994              } else {
1995                parameters[iParam].setStringValue(field);
1996              }
1997              std::string fileName;
1998              bool canOpen=false;
1999              if (field[0]=='/'||field[0]=='\\') {
2000                fileName = field;
2001              } else if (field[0]=='~') {
2002                char * environVar = getenv("HOME");
2003                if (environVar) {
2004                  std::string home(environVar);
2005                  field=field.erase(0,1);
2006                  fileName = home+field;
2007                } else {
2008                  fileName=field;
2009                }
2010              } else {
2011                fileName = directory+field;
2012              }
2013              FILE *fp=fopen(fileName.c_str(),"w");
2014              if (fp) {
2015                // can open - lets go for it
2016                fclose(fp);
2017                canOpen=true;
2018              } else {
2019                std::cout<<"Unable to open file "<<fileName<<std::endl;
2020              }
2021              if (canOpen) {
2022                // If presolve on then save presolved
2023                bool deleteModel2=false;
2024                ClpSimplex * model2 = lpSolver;
2025                if (preSolve) {
2026                  ClpPresolve pinfo;
2027                  int presolveOptions2 = presolveOptions&~0x40000000;
2028                  if ((presolveOptions2&0xffff)!=0)
2029                    pinfo.setPresolveActions(presolveOptions2);
2030                  if ((printOptions&1)!=0)
2031                    pinfo.statistics();
2032                  double presolveTolerance = 
2033                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
2034                  model2 = 
2035                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
2036                                         true,preSolve);
2037                  if (model2) {
2038                    printf("Saving presolved model on %s\n",
2039                           fileName.c_str());
2040                    deleteModel2=true;
2041                  } else {
2042                    printf("Presolved model looks infeasible - saving original on %s\n",
2043                           fileName.c_str());
2044                    deleteModel2=false;
2045                    model2 = lpSolver;
2046
2047                  }
2048                } else {
2049                  printf("Saving model on %s\n",
2050                           fileName.c_str());
2051                }
2052#if 0
2053                // Convert names
2054                int iRow;
2055                int numberRows=model2->numberRows();
2056                int iColumn;
2057                int numberColumns=model2->numberColumns();
2058
2059                char ** rowNames = NULL;
2060                char ** columnNames = NULL;
2061                if (model2->lengthNames()) {
2062                  rowNames = new char * [numberRows];
2063                  for (iRow=0;iRow<numberRows;iRow++) {
2064                    rowNames[iRow] =
2065                      strdup(model2->rowName(iRow).c_str());
2066#ifdef STRIPBLANKS
2067                    char * xx = rowNames[iRow];
2068                    int i;
2069                    int length = strlen(xx);
2070                    int n=0;
2071                    for (i=0;i<length;i++) {
2072                      if (xx[i]!=' ')
2073                        xx[n++]=xx[i];
2074                    }
2075                    xx[n]='\0';
2076#endif
2077                  }
2078                 
2079                  columnNames = new char * [numberColumns];
2080                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2081                    columnNames[iColumn] =
2082                      strdup(model2->columnName(iColumn).c_str());
2083#ifdef STRIPBLANKS
2084                    char * xx = columnNames[iColumn];
2085                    int i;
2086                    int length = strlen(xx);
2087                    int n=0;
2088                    for (i=0;i<length;i++) {
2089                      if (xx[i]!=' ')
2090                        xx[n++]=xx[i];
2091                    }
2092                    xx[n]='\0';
2093#endif
2094                  }
2095                }
2096                CoinMpsIO writer;
2097                writer.setMpsData(*model2->matrix(), COIN_DBL_MAX,
2098                                  model2->getColLower(), model2->getColUpper(),
2099                                  model2->getObjCoefficients(),
2100                                  (const char*) 0 /*integrality*/,
2101                                  model2->getRowLower(), model2->getRowUpper(),
2102                                  columnNames, rowNames);
2103                // Pass in array saying if each variable integer
2104                writer.copyInIntegerInformation(model2->integerInformation());
2105                writer.setObjectiveOffset(model2->objectiveOffset());
2106                writer.writeMps(fileName.c_str(),0,1,1);
2107                if (rowNames) {
2108                  for (iRow=0;iRow<numberRows;iRow++) {
2109                    free(rowNames[iRow]);
2110                  }
2111                  delete [] rowNames;
2112                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2113                    free(columnNames[iColumn]);
2114                  }
2115                  delete [] columnNames;
2116                }
2117#else
2118                model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
2119#endif
2120                if (deleteModel2)
2121                  delete model2;
2122                time2 = CoinCpuTime();
2123                totalTime += time2-time1;
2124                time1=time2;
2125              }
2126            } else {
2127              std::cout<<"** Current model not valid"<<std::endl;
2128            }
2129            break;
2130          case BASISIN:
2131            if (goodModel) {
2132              // get next field
2133              field = CoinReadGetString(argc,argv);
2134              if (field=="$") {
2135                field = parameters[iParam].stringValue();
2136              } else if (field=="EOL") {
2137                parameters[iParam].printString();
2138                break;
2139              } else {
2140                parameters[iParam].setStringValue(field);
2141              }
2142              std::string fileName;
2143              bool canOpen=false;
2144              if (field=="-") {
2145                // stdin
2146                canOpen=true;
2147                fileName = "-";
2148              } else {
2149                if (field[0]=='/'||field[0]=='\\') {
2150                  fileName = field;
2151                } else if (field[0]=='~') {
2152                  char * environVar = getenv("HOME");
2153                  if (environVar) {
2154                    std::string home(environVar);
2155                    field=field.erase(0,1);
2156                    fileName = home+field;
2157                  } else {
2158                    fileName=field;
2159                  }
2160                } else {
2161                  fileName = directory+field;
2162                }
2163                FILE *fp=fopen(fileName.c_str(),"r");
2164                if (fp) {
2165                  // can open - lets go for it
2166                  fclose(fp);
2167                  canOpen=true;
2168                } else {
2169                  std::cout<<"Unable to open file "<<fileName<<std::endl;
2170                }
2171              }
2172              if (canOpen) {
2173                int values = lpSolver->readBasis(fileName.c_str());
2174                if (values==0)
2175                  basisHasValues=-1;
2176                else
2177                  basisHasValues=1;
2178              }
2179            } else {
2180              std::cout<<"** Current model not valid"<<std::endl;
2181            }
2182            break;
2183          case DEBUG:
2184            if (goodModel) {
2185              delete [] debugValues;
2186              debugValues=NULL;
2187              // get next field
2188              field = CoinReadGetString(argc,argv);
2189              if (field=="$") {
2190                field = parameters[iParam].stringValue();
2191              } else if (field=="EOL") {
2192                parameters[iParam].printString();
2193                break;
2194              } else {
2195                parameters[iParam].setStringValue(field);
2196                debugFile=field;
2197                if (debugFile=="create"||
2198                    debugFile=="createAfterPre") {
2199                  printf("Will create a debug file so this run should be a good one\n");
2200                  break;
2201                }
2202              }
2203              std::string fileName;
2204              if (field[0]=='/'||field[0]=='\\') {
2205                fileName = field;
2206              } else if (field[0]=='~') {
2207                char * environVar = getenv("HOME");
2208                if (environVar) {
2209                  std::string home(environVar);
2210                  field=field.erase(0,1);
2211                  fileName = home+field;
2212                } else {
2213                  fileName=field;
2214                }
2215              } else {
2216                fileName = directory+field;
2217              }
2218              FILE *fp=fopen(fileName.c_str(),"rb");
2219              if (fp) {
2220                // can open - lets go for it
2221                int numRows;
2222                double obj;
2223                fread(&numRows,sizeof(int),1,fp);
2224                fread(&numberDebugValues,sizeof(int),1,fp);
2225                fread(&obj,sizeof(double),1,fp);
2226                debugValues = new double[numberDebugValues+numRows];
2227                fread(debugValues,sizeof(double),numRows,fp);
2228                fread(debugValues,sizeof(double),numRows,fp);
2229                fread(debugValues,sizeof(double),numberDebugValues,fp);
2230                printf("%d doubles read into debugValues\n",numberDebugValues);
2231                fclose(fp);
2232              } else {
2233                std::cout<<"Unable to open file "<<fileName<<std::endl;
2234              }
2235            } else {
2236              std::cout<<"** Current model not valid"<<std::endl;
2237            }
2238            break;
2239          case BASISOUT:
2240            if (goodModel) {
2241              // get next field
2242              field = CoinReadGetString(argc,argv);
2243              if (field=="$") {
2244                field = parameters[iParam].stringValue();
2245              } else if (field=="EOL") {
2246                parameters[iParam].printString();
2247                break;
2248              } else {
2249                parameters[iParam].setStringValue(field);
2250              }
2251              std::string fileName;
2252              bool canOpen=false;
2253              if (field[0]=='/'||field[0]=='\\') {
2254                fileName = field;
2255              } else if (field[0]=='~') {
2256                char * environVar = getenv("HOME");
2257                if (environVar) {
2258                  std::string home(environVar);
2259                  field=field.erase(0,1);
2260                  fileName = home+field;
2261                } else {
2262                  fileName=field;
2263                }
2264              } else {
2265                fileName = directory+field;
2266              }
2267              FILE *fp=fopen(fileName.c_str(),"w");
2268              if (fp) {
2269                // can open - lets go for it
2270                fclose(fp);
2271                canOpen=true;
2272              } else {
2273                std::cout<<"Unable to open file "<<fileName<<std::endl;
2274              }
2275              if (canOpen) {
2276                ClpSimplex * model2 = lpSolver;
2277                model2->writeBasis(fileName.c_str(),outputFormat>1,outputFormat-2);
2278                time2 = CoinCpuTime();
2279                totalTime += time2-time1;
2280                time1=time2;
2281              }
2282            } else {
2283              std::cout<<"** Current model not valid"<<std::endl;
2284            }
2285            break;
2286          case SAVE:
2287            {
2288              // get next field
2289              field = CoinReadGetString(argc,argv);
2290              if (field=="$") {
2291                field = parameters[iParam].stringValue();
2292              } else if (field=="EOL") {
2293                parameters[iParam].printString();
2294                break;
2295              } else {
2296                parameters[iParam].setStringValue(field);
2297              }
2298              std::string fileName;
2299              bool canOpen=false;
2300              if (field[0]=='/'||field[0]=='\\') {
2301                fileName = field;
2302              } else if (field[0]=='~') {
2303                char * environVar = getenv("HOME");
2304                if (environVar) {
2305                  std::string home(environVar);
2306                  field=field.erase(0,1);
2307                  fileName = home+field;
2308                } else {
2309                  fileName=field;
2310                }
2311              } else {
2312                fileName = directory+field;
2313              }
2314              FILE *fp=fopen(fileName.c_str(),"wb");
2315              if (fp) {
2316                // can open - lets go for it
2317                fclose(fp);
2318                canOpen=true;
2319              } else {
2320                std::cout<<"Unable to open file "<<fileName<<std::endl;
2321              }
2322              if (canOpen) {
2323                int status;
2324                // If presolve on then save presolved
2325                bool deleteModel2=false;
2326                ClpSimplex * model2 = lpSolver;
2327                if (preSolve) {
2328                  ClpPresolve pinfo;
2329                  double presolveTolerance = 
2330                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
2331                  model2 = 
2332                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
2333                                         false,preSolve);
2334                  if (model2) {
2335                    printf("Saving presolved model on %s\n",
2336                           fileName.c_str());
2337                    deleteModel2=true;
2338                  } else {
2339                    printf("Presolved model looks infeasible - saving original on %s\n",
2340                           fileName.c_str());
2341                    deleteModel2=false;
2342                    model2 = lpSolver;
2343
2344                  }
2345                } else {
2346                  printf("Saving model on %s\n",
2347                           fileName.c_str());
2348                }
2349                status =model2->saveModel(fileName.c_str());
2350                if (deleteModel2)
2351                  delete model2;
2352                if (!status) {
2353                  goodModel=true;
2354                  time2 = CoinCpuTime();
2355                  totalTime += time2-time1;
2356                  time1=time2;
2357                } else {
2358                  // errors
2359                  std::cout<<"There were errors on output"<<std::endl;
2360                }
2361              }
2362            }
2363            break;
2364          case RESTORE:
2365            {
2366              // get next field
2367              field = CoinReadGetString(argc,argv);
2368              if (field=="$") {
2369                field = parameters[iParam].stringValue();
2370              } else if (field=="EOL") {
2371                parameters[iParam].printString();
2372                break;
2373              } else {
2374                parameters[iParam].setStringValue(field);
2375              }
2376              std::string fileName;
2377              bool canOpen=false;
2378              if (field[0]=='/'||field[0]=='\\') {
2379                fileName = field;
2380              } else if (field[0]=='~') {
2381                char * environVar = getenv("HOME");
2382                if (environVar) {
2383                  std::string home(environVar);
2384                  field=field.erase(0,1);
2385                  fileName = home+field;
2386                } else {
2387                  fileName=field;
2388                }
2389              } else {
2390                fileName = directory+field;
2391              }
2392              FILE *fp=fopen(fileName.c_str(),"rb");
2393              if (fp) {
2394                // can open - lets go for it
2395                fclose(fp);
2396                canOpen=true;
2397              } else {
2398                std::cout<<"Unable to open file "<<fileName<<std::endl;
2399              }
2400              if (canOpen) {
2401                int status =lpSolver->restoreModel(fileName.c_str());
2402                if (!status) {
2403                  goodModel=true;
2404                  time2 = CoinCpuTime();
2405                  totalTime += time2-time1;
2406                  time1=time2;
2407                } else {
2408                  // errors
2409                  std::cout<<"There were errors on input"<<std::endl;
2410                }
2411              }
2412            }
2413            break;
2414          case MAXIMIZE:
2415            lpSolver->setOptimizationDirection(-1);
2416            break;
2417          case MINIMIZE:
2418            lpSolver->setOptimizationDirection(1);
2419            break;
2420          case ALLSLACK:
2421            lpSolver->allSlackBasis(true);
2422            break;
2423          case REVERSE:
2424            if (goodModel) {
2425              int iColumn;
2426              int numberColumns=lpSolver->numberColumns();
2427              double * dualColumnSolution = 
2428                lpSolver->dualColumnSolution();
2429              ClpObjective * obj = lpSolver->objectiveAsObject();
2430              assert(dynamic_cast<ClpLinearObjective *> (obj));
2431              double offset;
2432              double * objective = obj->gradient(NULL,NULL,offset,true);
2433              for (iColumn=0;iColumn<numberColumns;iColumn++) {
2434                dualColumnSolution[iColumn] = dualColumnSolution[iColumn];
2435                objective[iColumn] = -objective[iColumn];
2436              }
2437              int iRow;
2438              int numberRows=lpSolver->numberRows();
2439              double * dualRowSolution = 
2440                lpSolver->dualRowSolution();
2441              for (iRow=0;iRow<numberRows;iRow++) 
2442                dualRowSolution[iRow] = dualRowSolution[iRow];
2443            }
2444            break;
2445          case DIRECTORY:
2446            {
2447              std::string name = CoinReadGetString(argc,argv);
2448              if (name!="EOL") {
2449                int length=name.length();
2450                if (name[length-1]=='/'||name[length-1]=='\\')
2451                  directory=name;
2452                else
2453                  directory = name+"/";
2454                parameters[iParam].setStringValue(directory);
2455              } else {
2456                parameters[iParam].printString();
2457              }
2458            }
2459            break;
2460          case STDIN:
2461            CbcOrClpRead_mode=-1;
2462            break;
2463          case NETLIB_DUAL:
2464          case NETLIB_EITHER:
2465          case NETLIB_BARRIER:
2466          case NETLIB_PRIMAL:
2467          case NETLIB_TUNE:
2468            {
2469              // create fields for unitTest
2470              const char * fields[4];
2471              int nFields=2;
2472              fields[0]="fake main from unitTest";
2473              fields[1]="-netlib";
2474              if (directory!=defaultDirectory) {
2475                fields[2]="-netlibDir";
2476                fields[3]=directory.c_str();
2477                nFields=4;
2478              }
2479              int algorithm;
2480              if (type==NETLIB_DUAL) {
2481                std::cerr<<"Doing netlib with dual algorithm"<<std::endl;
2482                algorithm =0;
2483              } else if (type==NETLIB_BARRIER) {
2484                std::cerr<<"Doing netlib with barrier algorithm"<<std::endl;
2485                algorithm =2;
2486              } else if (type==NETLIB_EITHER) {
2487                std::cerr<<"Doing netlib with dual or primal algorithm"<<std::endl;
2488                algorithm =3;
2489              } else if (type==NETLIB_TUNE) {
2490                std::cerr<<"Doing netlib with best algorithm!"<<std::endl;
2491                algorithm =5;
2492                // uncomment next to get active tuning
2493                // algorithm=6;
2494              } else {
2495                std::cerr<<"Doing netlib with primal agorithm"<<std::endl;
2496                algorithm=1;
2497              }
2498              int specialOptions = lpSolver->specialOptions();
2499              lpSolver->setSpecialOptions(0);
2500              mainTest(nFields,fields,algorithm,*lpSolver,
2501                       (preSolve!=0),specialOptions);
2502            }
2503            break;
2504          case UNITTEST:
2505            {
2506              // create fields for unitTest
2507              const char * fields[3];
2508              int nFields=1;
2509              fields[0]="fake main from unitTest";
2510              if (directory!=defaultDirectory) {
2511                fields[1]="-mpsDir";
2512                fields[2]=directory.c_str();
2513                nFields=3;
2514              }
2515              mainTest(nFields,fields,false,*lpSolver,(preSolve!=0),
2516                       false);
2517            }
2518            break;
2519          case FAKEBOUND:
2520            if (goodModel) {
2521              // get bound
2522              double value = CoinReadGetDoubleField(argc,argv,&valid);
2523              if (!valid) {
2524                std::cout<<"Setting "<<parameters[iParam].name()<<
2525                  " to DEBUG "<<value<<std::endl;
2526                int iRow;
2527                int numberRows=lpSolver->numberRows();
2528                double * rowLower = lpSolver->rowLower();
2529                double * rowUpper = lpSolver->rowUpper();
2530                for (iRow=0;iRow<numberRows;iRow++) {
2531                  // leave free ones for now
2532                  if (rowLower[iRow]>-1.0e20||rowUpper[iRow]<1.0e20) {
2533                    rowLower[iRow]=CoinMax(rowLower[iRow],-value);
2534                    rowUpper[iRow]=CoinMin(rowUpper[iRow],value);
2535                  }
2536                }
2537                int iColumn;
2538                int numberColumns=lpSolver->numberColumns();
2539                double * columnLower = lpSolver->columnLower();
2540                double * columnUpper = lpSolver->columnUpper();
2541                for (iColumn=0;iColumn<numberColumns;iColumn++) {
2542                  // leave free ones for now
2543                  if (columnLower[iColumn]>-1.0e20||
2544                      columnUpper[iColumn]<1.0e20) {
2545                    columnLower[iColumn]=CoinMax(columnLower[iColumn],-value);
2546                    columnUpper[iColumn]=CoinMin(columnUpper[iColumn],value);
2547                  }
2548                }
2549              } else if (valid==1) {
2550                abort();
2551              } else {
2552                std::cout<<"enter value for "<<parameters[iParam].name()<<
2553                  std::endl;
2554              }
2555            }
2556            break;
2557          case REALLY_SCALE:
2558            if (goodModel) {
2559              ClpSimplex newModel(*lpSolver,
2560                                  lpSolver->scalingFlag());
2561              printf("model really really scaled\n");
2562              *lpSolver=newModel;
2563            }
2564            break;
2565          case HELP:
2566            std::cout<<"Coin Solver version "<<CBCVERSION
2567                     <<", build "<<__DATE__<<std::endl;
2568            std::cout<<"Non default values:-"<<std::endl;
2569            std::cout<<"Perturbation "<<lpSolver->perturbation()<<" (default 100)"
2570                     <<std::endl;
2571            CoinReadPrintit(
2572                    "Presolve being done with 5 passes\n\
2573Dual steepest edge steep/partial on matrix shape and factorization density\n\
2574Clpnnnn taken out of messages\n\
2575If Factorization frequency default then done on size of matrix\n\n\
2576(-)unitTest, (-)netlib or (-)netlibp will do standard tests\n\n\
2577You can switch to interactive mode at any time so\n\
2578clp watson.mps -scaling off -primalsimplex\nis the same as\n\
2579clp watson.mps -\nscaling off\nprimalsimplex"
2580                    );
2581            break;
2582          case SOLUTION:
2583            if (goodModel) {
2584              // get next field
2585              field = CoinReadGetString(argc,argv);
2586              if (field=="$") {
2587                field = parameters[iParam].stringValue();
2588              } else if (field=="EOL") {
2589                parameters[iParam].printString();
2590                break;
2591              } else {
2592                parameters[iParam].setStringValue(field);
2593              }
2594              std::string fileName;
2595              FILE *fp=NULL;
2596              if (field=="-"||field=="EOL"||field=="stdout") {
2597                // stdout
2598                fp=stdout;
2599              } else if (field=="stderr") {
2600                // stderr
2601                fp=stderr;
2602              } else {
2603                if (field[0]=='/'||field[0]=='\\') {
2604                  fileName = field;
2605                } else if (field[0]=='~') {
2606                  char * environVar = getenv("HOME");
2607                  if (environVar) {
2608                    std::string home(environVar);
2609                    field=field.erase(0,1);
2610                    fileName = home+field;
2611                  } else {
2612                    fileName=field;
2613                  }
2614                } else {
2615                  fileName = directory+field;
2616                }
2617                fp=fopen(fileName.c_str(),"w");
2618              }
2619              if (fp) {
2620                // make fancy later on
2621                int iRow;
2622                int numberRows=lpSolver->numberRows();
2623                double * dualRowSolution = lpSolver->dualRowSolution();
2624                double * primalRowSolution = 
2625                  lpSolver->primalRowSolution();
2626                double * rowLower = lpSolver->rowLower();
2627                double * rowUpper = lpSolver->rowUpper();
2628                double primalTolerance = lpSolver->primalTolerance();
2629                char format[6];
2630                sprintf(format,"%%-%ds",CoinMax(lengthName,8));
2631                if (printMode>2) {
2632                  for (iRow=0;iRow<numberRows;iRow++) {
2633                    int type=printMode-3;
2634                    if (primalRowSolution[iRow]>rowUpper[iRow]+primalTolerance||
2635                        primalRowSolution[iRow]<rowLower[iRow]-primalTolerance) {
2636                      fprintf(fp,"** ");
2637                      type=2;
2638                    } else if (fabs(primalRowSolution[iRow])>1.0e-8) {
2639                      type=1;
2640                    } else if (numberRows<50) {
2641                      type=3;
2642                    }
2643                    if (type) {
2644                      fprintf(fp,"%7d ",iRow);
2645                      if (lengthName)
2646                        fprintf(fp,format,rowNames[iRow].c_str());
2647                      fprintf(fp,"%15.8g        %15.8g\n",primalRowSolution[iRow],
2648                              dualRowSolution[iRow]);
2649                    }
2650                  }
2651                }
2652                int iColumn;
2653                int numberColumns=lpSolver->numberColumns();
2654                double * dualColumnSolution = 
2655                  lpSolver->dualColumnSolution();
2656                double * primalColumnSolution = 
2657                  lpSolver->primalColumnSolution();
2658                double * columnLower = lpSolver->columnLower();
2659                double * columnUpper = lpSolver->columnUpper();
2660                if (printMode!=2) {
2661                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2662                    int type=0;
2663                    if (primalColumnSolution[iColumn]>columnUpper[iColumn]+primalTolerance||
2664                        primalColumnSolution[iColumn]<columnLower[iColumn]-primalTolerance) {
2665                      fprintf(fp,"** ");
2666                      type=2;
2667                    } else if (fabs(primalColumnSolution[iColumn])>1.0e-8) {
2668                      type=1;
2669                    } else if (numberColumns<50) {
2670                      type=3;
2671                    }
2672                    // see if integer
2673                    if ((!lpSolver->isInteger(iColumn)||fabs(primalColumnSolution[iColumn])<1.0e-8)
2674                         &&printMode==1)
2675                      type=0;
2676                    if (type) {
2677                      fprintf(fp,"%7d ",iColumn);
2678                      if (lengthName)
2679                        fprintf(fp,format,columnNames[iColumn].c_str());
2680                      fprintf(fp,"%15.8g        %15.8g\n",
2681                              primalColumnSolution[iColumn],
2682                              dualColumnSolution[iColumn]);
2683                    }
2684                  }
2685                } else {
2686                  // special format suitable for OsiRowCutDebugger
2687                  int n=0;
2688                  bool comma=false;
2689                  bool newLine=false;
2690                  fprintf(fp,"\tint intIndicesV[]={\n");
2691                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2692                    if(primalColumnSolution[iColumn]>0.5&&model.solver()->isInteger(iColumn)) {
2693                      if (comma)
2694                        fprintf(fp,",");
2695                      if (newLine)
2696                        fprintf(fp,"\n");
2697                      fprintf(fp,"%d ",iColumn);
2698                      comma=true;
2699                      newLine=false;
2700                      n++;
2701                      if (n==10) {
2702                        n=0;
2703                        newLine=true;
2704                      }
2705                    }
2706                  }
2707                  fprintf(fp,"};\n");
2708                  n=0;
2709                  comma=false;
2710                  newLine=false;
2711                  fprintf(fp,"\tdouble intSolnV[]={\n");
2712                  for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2713                    if(primalColumnSolution[iColumn]>0.5&&model.solver()->isInteger(iColumn)) {
2714                      if (comma)
2715                        fprintf(fp,",");
2716                      if (newLine)
2717                        fprintf(fp,"\n");
2718                      int value = (int) (primalColumnSolution[iColumn]+0.5);
2719                      fprintf(fp,"%d. ",value);
2720                      comma=true;
2721                      newLine=false;
2722                      n++;
2723                      if (n==10) {
2724                        n=0;
2725                        newLine=true;
2726                      }
2727                    }
2728                  }
2729                  fprintf(fp,"};\n");
2730                }
2731                if (fp!=stdout)
2732                  fclose(fp);
2733              } else {
2734                std::cout<<"Unable to open file "<<fileName<<std::endl;
2735              }
2736            } else {
2737              std::cout<<"** Current model not valid"<<std::endl;
2738             
2739            }
2740            break;
2741          case SAVESOL:
2742            if (goodModel) {
2743              // get next field
2744              field = CoinReadGetString(argc,argv);
2745              if (field=="$") {
2746                field = parameters[iParam].stringValue();
2747              } else if (field=="EOL") {
2748                parameters[iParam].printString();
2749                break;
2750              } else {
2751                parameters[iParam].setStringValue(field);
2752              }
2753              std::string fileName;
2754              if (field[0]=='/'||field[0]=='\\') {
2755                fileName = field;
2756              } else if (field[0]=='~') {
2757                char * environVar = getenv("HOME");
2758                if (environVar) {
2759                  std::string home(environVar);
2760                  field=field.erase(0,1);
2761                  fileName = home+field;
2762                } else {
2763                  fileName=field;
2764                }
2765              } else {
2766                fileName = directory+field;
2767              }
2768              saveSolution(lpSolver,fileName);
2769            } else {
2770              std::cout<<"** Current model not valid"<<std::endl;
2771             
2772            }
2773            break;
2774          default:
2775            abort();
2776          }
2777        } 
2778      } else if (!numberMatches) {
2779        std::cout<<"No match for "<<field<<" - ? for list of commands"
2780                 <<std::endl;
2781      } else if (numberMatches==1) {
2782        if (!numberQuery) {
2783          std::cout<<"Short match for "<<field<<" - completion: ";
2784          std::cout<<parameters[firstMatch].matchName()<<std::endl;
2785        } else if (numberQuery) {
2786          std::cout<<parameters[firstMatch].matchName()<<" : ";
2787          std::cout<<parameters[firstMatch].shortHelp()<<std::endl;
2788          if (numberQuery>=2) 
2789            parameters[firstMatch].printLongHelp();
2790        }
2791      } else {
2792        if (!numberQuery) 
2793          std::cout<<"Multiple matches for "<<field<<" - possible completions:"
2794                   <<std::endl;
2795        else
2796          std::cout<<"Completions of "<<field<<":"<<std::endl;
2797        for ( iParam=0; iParam<numberParameters; iParam++ ) {
2798          int match = parameters[iParam].matches(field);
2799          if (match&&parameters[iParam].displayThis()) {
2800            std::cout<<parameters[iParam].matchName();
2801            if (numberQuery>=2) 
2802              std::cout<<" : "<<parameters[iParam].shortHelp();
2803            std::cout<<std::endl;
2804          }
2805        }
2806      }
2807    }
2808  }
2809  // By now all memory should be freed
2810#ifdef DMALLOC
2811  dmalloc_log_unfreed();
2812  dmalloc_shutdown();
2813#endif
2814  return 0;
2815}   
2816/*
2817  Version 1.00.00 November 16 2005.
2818  This is to stop me (JJF) messing about too much.
2819  Tuning changes should be noted here.
2820  The testing next version may be activated by CBC_NEXT_VERSION
2821  This applies to OsiClp, Clp etc
2822  Version 1.00.01 November 24 2005
2823  Added several classes for advanced users.  This can't affect code (if you don't use it)
2824  Made some tiny changes (for N way branching) which should not change anything.
2825  CbcNWay object class - for N way branching this also allows use of CbcConsequence class.
2826  CbcBranchAllDifferent object class - for branching on general integer variables
2827  to stop them having same value so branches are x >= y+1 and x <= y-1.
2828  Added two new Cgl classes - CglAllDifferent which does column fixing (too slowly)
2829  and CglStored which just has a list of cuts which can be activated.
2830  Modified preprocess option to SOS
2831  Version 1.00.02 December 9 2005
2832  Added use of CbcStrategy to do clean preprocessing
2833  Added use of referenceSolver for cleaner repetition of Cbc
2834*/
Note: See TracBrowser for help on using the repository browser.