source: trunk/Test/CoinSolve.cpp @ 212

Last change on this file since 212 was 212, checked in by forrest, 16 years ago

1.00

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