source: trunk/Cbc/src/CbcSolver.cpp @ 1212

Last change on this file since 1212 was 1212, checked in by forrest, 10 years ago

fixes

  • Property svn:keywords set to Id
File size: 386.6 KB
</
Line 
1/* $Id: CbcSolver.cpp 1212 2009-08-21 16:19:13Z forrest $ */
2// Copyright (C) 2007, International Business Machines
3// Corporation and others.  All Rights Reserved.
4   
5#include "CbcConfig.h"
6#include "CoinPragma.hpp"
7
8#include <cassert>
9#include <cstdio>
10#include <cstdlib>
11#include <cmath>
12#include <cfloat>
13#include <cstring>
14#include <iostream>
15
16#include "CoinPragma.hpp"
17#include "CoinHelperFunctions.hpp"
18
19#include "CoinMpsIO.hpp"
20#include "CoinModel.hpp"
21
22#include "ClpFactorization.hpp"
23#include "ClpQuadraticObjective.hpp"
24#include "CoinTime.hpp"
25#include "ClpSimplex.hpp"
26#include "ClpSimplexOther.hpp"
27#include "ClpSolve.hpp"
28#include "ClpMessage.hpp"
29#include "ClpPackedMatrix.hpp"
30#include "ClpPlusMinusOneMatrix.hpp"
31#include "ClpNetworkMatrix.hpp"
32#include "ClpDualRowSteepest.hpp"
33#include "ClpDualRowDantzig.hpp"
34#include "ClpLinearObjective.hpp"
35#include "ClpPrimalColumnSteepest.hpp"
36#include "ClpPrimalColumnDantzig.hpp"
37
38#include "ClpPresolve.hpp"
39#include "CbcOrClpParam.hpp"
40#include "OsiRowCutDebugger.hpp"
41#include "OsiChooseVariable.hpp"
42#include "OsiAuxInfo.hpp"
43// Version
44#ifndef CBCVERSION
45#define CBCVERSION "2.31.00"
46#endif
47//#define ORBITAL
48#ifdef ORBITAL
49#include "CbcOrbital.hpp"
50#endif
51//#define USER_HAS_FAKE_CLP
52//#define USER_HAS_FAKE_CBC
53//#define CLP_MALLOC_STATISTICS
54#ifdef CLP_MALLOC_STATISTICS
55#include <malloc.h>
56#include <exception>
57#include <new>
58#include "stolen_from_ekk_malloc.cpp"
59static double malloc_times=0.0;
60static double malloc_total=0.0;
61static int malloc_amount[]={0,32,128,256,1024,4096,16384,65536,262144,INT_MAX};
62static int malloc_n=10;
63double malloc_counts[10]={0,0,0,0,0,0,0,0,0,0};
64bool malloc_counts_on=true;
65void * operator new (size_t size) throw (std::bad_alloc)
66{
67  malloc_times ++;
68  malloc_total += size;
69  int i;
70  for (i=0;i<malloc_n;i++) {
71    if ((int) size<=malloc_amount[i]) {
72      malloc_counts[i]++;
73      break;
74    }
75  }
76#ifdef DEBUG_MALLOC
77  void *p;
78  if (malloc_counts_on)
79    p =stolen_from_ekk_mallocBase(size);
80  else
81    p =malloc(size);
82#else
83  void * p =malloc(size);
84#endif
85  //char * xx = (char *) p;
86  //memset(xx,0,size);
87  // Initialize random seed
88  //CoinSeedRandom(987654321);
89  return p;
90}
91void operator delete (void *p) throw()
92{
93#ifdef DEBUG_MALLOC
94  if (malloc_counts_on)
95    stolen_from_ekk_freeBase(p);
96  else
97    free(p);
98#else
99  free(p);
100#endif
101}
102static void malloc_stats2()
103{
104  double average = malloc_total/malloc_times;
105  printf("count %g bytes %g - average %g\n",malloc_times,malloc_total,average);
106  for (int i=0;i<malloc_n;i++) 
107    printf("%g ",malloc_counts[i]);
108  printf("\n");
109  malloc_times=0.0;
110  malloc_total=0.0;
111  memset(malloc_counts,0,sizeof(malloc_counts));
112  // print results
113}
114#else
115//void stolen_from_ekk_memory(void * dummy,int type)
116//{
117//}
118//bool malloc_counts_on=false;
119#endif
120//#define DMALLOC
121#ifdef DMALLOC
122#include "dmalloc.h"
123#endif
124#ifdef WSSMP_BARRIER
125#define FOREIGN_BARRIER
126#endif
127#ifdef UFL_BARRIER
128#define FOREIGN_BARRIER
129#endif
130#ifdef TAUCS_BARRIER
131#define FOREIGN_BARRIER
132#endif
133static int initialPumpTune=-1;
134#include "CoinWarmStartBasis.hpp"
135
136#include "OsiSolverInterface.hpp"
137#include "OsiCuts.hpp"
138#include "OsiRowCut.hpp"
139#include "OsiColCut.hpp"
140#ifndef COIN_HAS_LINK
141#define COIN_HAS_LINK
142#endif
143#ifdef COIN_HAS_LINK
144#include "CbcLinked.hpp"
145#endif
146#include "CglPreProcess.hpp"
147#include "CglCutGenerator.hpp"
148#include "CglGomory.hpp"
149#include "CglProbing.hpp"
150#include "CglKnapsackCover.hpp"
151#include "CglRedSplit.hpp"
152#include "CglClique.hpp"
153#include "CglFlowCover.hpp"
154#include "CglMixedIntegerRounding2.hpp"
155#include "CglTwomir.hpp"
156#include "CglDuplicateRow.hpp"
157#include "CglStored.hpp"
158#include "CglLandP.hpp"
159#include "CglResidualCapacity.hpp"
160#ifdef ZERO_HALF_CUTS
161#include "CglZeroHalf.hpp"
162#endif
163
164#include "CbcModel.hpp"
165#include "CbcHeuristic.hpp"
166#include "CbcHeuristicLocal.hpp"
167#include "CbcHeuristicPivotAndFix.hpp"
168#include "CbcHeuristicRandRound.hpp"
169#include "CbcHeuristicGreedy.hpp"
170#include "CbcHeuristicFPump.hpp"
171#include "CbcHeuristicRINS.hpp"
172#include "CbcHeuristicDiveCoefficient.hpp"
173#include "CbcHeuristicDiveFractional.hpp"
174#include "CbcHeuristicDiveGuided.hpp"
175#include "CbcHeuristicDiveVectorLength.hpp"
176#include "CbcHeuristicDivePseudoCost.hpp"
177#include "CbcHeuristicDiveLineSearch.hpp"
178#include "CbcTreeLocal.hpp"
179#include "CbcCompareActual.hpp"
180#include "CbcBranchActual.hpp"
181#include "CbcBranchLotsize.hpp"
182#include  "CbcOrClpParam.hpp"
183#include  "CbcCutGenerator.hpp"
184#include  "CbcStrategy.hpp"
185#include "CbcBranchCut.hpp"
186
187#include "OsiClpSolverInterface.hpp"
188#include "CbcSolver.hpp"
189//#define IN_BRANCH_AND_BOUND (0x01000000|262144)
190#define IN_BRANCH_AND_BOUND (0x01000000|262144|128|1024|2048)
191//#define IN_BRANCH_AND_BOUND (0x01000000|262144|128)
192CbcSolver::CbcSolver()
193  : babModel_(NULL),
194    userFunction_(NULL),
195    statusUserFunction_(NULL),
196    originalSolver_(NULL),
197    originalCoinModel_(NULL),
198    cutGenerator_(NULL),
199    numberUserFunctions_(0),
200    numberCutGenerators_(0),
201    startTime_(CoinCpuTime()),
202    parameters_(NULL),
203    numberParameters_(0),
204    doMiplib_(false),
205    noPrinting_(false),
206    readMode_(1)
207{
208  callBack_ = new CbcStopNow();
209  fillParameters();
210}
211CbcSolver::CbcSolver(const OsiClpSolverInterface & solver)
212  : babModel_(NULL),
213    userFunction_(NULL),
214    statusUserFunction_(NULL),
215    originalSolver_(NULL),
216    originalCoinModel_(NULL),
217    cutGenerator_(NULL),
218    numberUserFunctions_(0),
219    numberCutGenerators_(0),
220    startTime_(CoinCpuTime()),
221    parameters_(NULL),
222    numberParameters_(0),
223    doMiplib_(false),
224    noPrinting_(false),
225    readMode_(1)
226{
227  callBack_ = new CbcStopNow();
228  model_ = CbcModel(solver);
229  fillParameters();
230}
231CbcSolver::CbcSolver(const CbcModel & solver)
232  : babModel_(NULL),
233    userFunction_(NULL),
234    statusUserFunction_(NULL),
235    originalSolver_(NULL),
236    originalCoinModel_(NULL),
237    cutGenerator_(NULL),
238    numberUserFunctions_(0),
239    numberCutGenerators_(0),
240    startTime_(CoinCpuTime()),
241    parameters_(NULL),
242    numberParameters_(0),
243    doMiplib_(false),
244    noPrinting_(false),
245    readMode_(1)
246{
247  callBack_ = new CbcStopNow();
248  model_ = solver;
249  fillParameters();
250}
251CbcSolver::~CbcSolver()
252{
253  int i;
254  for (i=0;i<numberUserFunctions_;i++)
255    delete userFunction_[i];
256  delete [] userFunction_;
257  for (i=0;i<numberCutGenerators_;i++)
258    delete cutGenerator_[i];
259  delete [] cutGenerator_;
260  delete [] statusUserFunction_;
261  delete originalSolver_;
262  delete originalCoinModel_;
263  delete babModel_;
264  delete [] parameters_;
265  delete callBack_;
266}
267// Copy constructor
268CbcSolver::CbcSolver ( const CbcSolver & rhs)
269  : model_(rhs.model_),
270    babModel_(NULL),
271    userFunction_(NULL),
272    statusUserFunction_(NULL),
273    numberUserFunctions_(rhs.numberUserFunctions_),
274    startTime_(CoinCpuTime()),
275    parameters_(NULL),
276    numberParameters_(rhs.numberParameters_),
277    doMiplib_(rhs.doMiplib_),
278    noPrinting_(rhs.noPrinting_),
279    readMode_(rhs.readMode_)
280{
281  fillParameters();
282  if (rhs.babModel_)
283    babModel_ = new CbcModel(*rhs.babModel_);
284  userFunction_ = new CbcUser * [numberUserFunctions_];
285  int i;
286  for (i=0;i<numberUserFunctions_;i++)
287    userFunction_[i] = rhs.userFunction_[i]->clone();
288  for (i=0;i<numberParameters_;i++)
289    parameters_[i]=rhs.parameters_[i];
290  for (i=0;i<numberCutGenerators_;i++)
291    cutGenerator_[i] = rhs.cutGenerator_[i]->clone();
292  callBack_ = rhs.callBack_->clone();
293  originalSolver_ = NULL;
294  if (rhs.originalSolver_) {
295    OsiSolverInterface * temp = rhs.originalSolver_->clone();
296    originalSolver_ = dynamic_cast<OsiClpSolverInterface *> (temp);
297    assert (originalSolver_);
298  }
299  originalCoinModel_ = NULL;
300  if (rhs.originalCoinModel_)
301    originalCoinModel_ = new CoinModel(*rhs.originalCoinModel_);
302}
303// Assignment operator
304CbcSolver &
305CbcSolver::operator=(const CbcSolver & rhs)
306{
307  if (this != &rhs) {
308    int i;
309    for (i=0;i<numberUserFunctions_;i++)
310      delete userFunction_[i];
311    delete [] userFunction_;
312    for (i=0;i<numberCutGenerators_;i++)
313      delete cutGenerator_[i];
314    delete [] cutGenerator_;
315    delete [] statusUserFunction_;
316    delete originalSolver_;
317    delete originalCoinModel_;
318    statusUserFunction_ = NULL;
319    delete babModel_;
320    delete [] parameters_;
321    delete callBack_;
322    numberUserFunctions_ = rhs.numberUserFunctions_;
323    startTime_ = rhs.startTime_;
324    numberParameters_= rhs.numberParameters_;
325    for (i=0;i<numberParameters_;i++)
326      parameters_[i]=rhs.parameters_[i];
327    for (i=0;i<numberCutGenerators_;i++)
328      cutGenerator_[i] = rhs.cutGenerator_[i]->clone();
329    noPrinting_ = rhs.noPrinting_;
330    readMode_ = rhs.readMode_;
331    doMiplib_ = rhs.doMiplib_;
332    model_ = rhs.model_;
333    if (rhs.babModel_)
334      babModel_ = new CbcModel(*rhs.babModel_);
335    else
336      babModel_ = NULL;
337    userFunction_ = new CbcUser * [numberUserFunctions_];
338    for (i=0;i<numberUserFunctions_;i++)
339      userFunction_[i] = rhs.userFunction_[i]->clone();
340    callBack_ = rhs.callBack_->clone();
341    originalSolver_ = NULL;
342    if (rhs.originalSolver_) {
343      OsiSolverInterface * temp = rhs.originalSolver_->clone();
344      originalSolver_ = dynamic_cast<OsiClpSolverInterface *> (temp);
345    assert (originalSolver_);
346    }
347    originalCoinModel_ = NULL;
348    if (rhs.originalCoinModel_)
349      originalCoinModel_ = new CoinModel(*rhs.originalCoinModel_);
350}
351  return *this;
352}
353// Get int value
354int CbcSolver::intValue(CbcOrClpParameterType type) const
355{
356  return parameters_[whichParam(type,numberParameters_,parameters_)].intValue();
357}
358// Set int value
359void CbcSolver::setIntValue(CbcOrClpParameterType type,int value)
360{
361  parameters_[whichParam(type,numberParameters_,parameters_)].setIntValue(value);
362}
363// Get double value
364double CbcSolver::doubleValue(CbcOrClpParameterType type) const
365{
366  return parameters_[whichParam(type,numberParameters_,parameters_)].doubleValue();
367}
368// Set double value
369void CbcSolver::setDoubleValue(CbcOrClpParameterType type,double value)
370{
371  parameters_[whichParam(type,numberParameters_,parameters_)].setDoubleValue(value);
372}
373// User function (NULL if no match)
374CbcUser * CbcSolver::userFunction(const char * name) const
375{
376  int i;
377  for (i=0;i<numberUserFunctions_;i++) {
378    if (!strcmp(name,userFunction_[i]->name().c_str()))
379      break;
380  }
381  if (i<numberUserFunctions_)
382    return userFunction_[i];
383  else
384    return NULL;
385}
386void CbcSolver::fillParameters()
387{
388  int maxParam = 200;
389  CbcOrClpParam * parameters = new CbcOrClpParam [maxParam];
390  numberParameters_=0 ;
391  establishParams(numberParameters_,parameters) ;
392  assert (numberParameters_<=maxParam);
393  parameters_ = new CbcOrClpParam [numberParameters_];
394  int i;
395  for (i=0;i<numberParameters_;i++)
396    parameters_[i]=parameters[i];
397  delete [] parameters;
398  const char dirsep =  CoinFindDirSeparator();
399  std::string directory;
400  std::string dirSample;
401  std::string dirNetlib;
402  std::string dirMiplib;
403  if (dirsep == '/') {
404    directory = "./";
405    dirSample = "../../Data/Sample/";
406    dirNetlib = "../../Data/Netlib/";
407    dirMiplib = "../../Data/miplib3/";
408  } else {
409    directory = ".\\";
410    dirSample = "..\\..\\Data\\Sample\\";
411    dirNetlib = "..\\..\\Data\\Netlib\\";
412    dirMiplib = "..\\..\\Data\\miplib3\\";
413  }
414  std::string defaultDirectory = directory;
415  std::string importFile ="";
416  std::string exportFile ="default.mps";
417  std::string importBasisFile ="";
418  std::string importPriorityFile ="";
419  std::string debugFile="";
420  std::string printMask="";
421  std::string exportBasisFile ="default.bas";
422  std::string saveFile ="default.prob";
423  std::string restoreFile ="default.prob";
424  std::string solutionFile ="stdout";
425  std::string solutionSaveFile ="solution.file";
426  int doIdiot=-1;
427  int outputFormat=2;
428  int substitution=3;
429  int dualize=3;
430  int preSolve=5;
431  int doSprint=-1;
432  int testOsiParameters=-1;
433  int createSolver=0;
434  ClpSimplex * lpSolver;
435  OsiClpSolverInterface * clpSolver;
436  if (model_.solver()) {
437    clpSolver = dynamic_cast<OsiClpSolverInterface *> (model_.solver());
438    assert (clpSolver);
439    lpSolver = clpSolver->getModelPtr();
440    assert (lpSolver);
441  } else {
442    lpSolver = new ClpSimplex();
443    clpSolver = new OsiClpSolverInterface(lpSolver,true);
444    createSolver =1 ;
445  }
446  parameters_[whichParam(BASISIN,numberParameters_,parameters_)].setStringValue(importBasisFile);
447  parameters_[whichParam(PRIORITYIN,numberParameters_,parameters_)].setStringValue(importPriorityFile);
448  parameters_[whichParam(BASISOUT,numberParameters_,parameters_)].setStringValue(exportBasisFile);
449  parameters_[whichParam(DEBUG,numberParameters_,parameters_)].setStringValue(debugFile);
450  parameters_[whichParam(PRINTMASK,numberParameters_,parameters_)].setStringValue(printMask);
451  parameters_[whichParam(DIRECTORY,numberParameters_,parameters_)].setStringValue(directory);
452  parameters_[whichParam(DIRSAMPLE,numberParameters_,parameters_)].setStringValue(dirSample);
453  parameters_[whichParam(DIRNETLIB,numberParameters_,parameters_)].setStringValue(dirNetlib);
454  parameters_[whichParam(DIRMIPLIB,numberParameters_,parameters_)].setStringValue(dirMiplib);
455  parameters_[whichParam(DUALBOUND,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualBound());
456  parameters_[whichParam(DUALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualTolerance());
457  parameters_[whichParam(EXPORT,numberParameters_,parameters_)].setStringValue(exportFile);
458  parameters_[whichParam(IDIOT,numberParameters_,parameters_)].setIntValue(doIdiot);
459  parameters_[whichParam(IMPORT,numberParameters_,parameters_)].setStringValue(importFile);
460  parameters_[whichParam(PRESOLVETOLERANCE,numberParameters_,parameters_)].setDoubleValue(1.0e-8);
461  int iParam = whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_);
462  int value=1;
463  clpSolver->messageHandler()->setLogLevel(1) ;
464  lpSolver->setLogLevel(1);
465  parameters_[iParam].setIntValue(value);
466  iParam = whichParam(LOGLEVEL,numberParameters_,parameters_);
467  model_.messageHandler()->setLogLevel(value);
468  parameters_[iParam].setIntValue(value);
469  parameters_[whichParam(MAXFACTOR,numberParameters_,parameters_)].setIntValue(lpSolver->factorizationFrequency());
470  parameters_[whichParam(MAXITERATION,numberParameters_,parameters_)].setIntValue(lpSolver->maximumIterations());
471  parameters_[whichParam(OUTPUTFORMAT,numberParameters_,parameters_)].setIntValue(outputFormat);
472  parameters_[whichParam(PRESOLVEPASS,numberParameters_,parameters_)].setIntValue(preSolve);
473  parameters_[whichParam(PERTVALUE,numberParameters_,parameters_)].setIntValue(lpSolver->perturbation());
474  parameters_[whichParam(PRIMALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->primalTolerance());
475  parameters_[whichParam(PRIMALWEIGHT,numberParameters_,parameters_)].setDoubleValue(lpSolver->infeasibilityCost());
476  parameters_[whichParam(RESTORE,numberParameters_,parameters_)].setStringValue(restoreFile);
477  parameters_[whichParam(SAVE,numberParameters_,parameters_)].setStringValue(saveFile);
478  //parameters_[whichParam(TIMELIMIT,numberParameters_,parameters_)].setDoubleValue(1.0e8);
479  parameters_[whichParam(TIMELIMIT_BAB,numberParameters_,parameters_)].setDoubleValue(1.0e8);
480  parameters_[whichParam(SOLUTION,numberParameters_,parameters_)].setStringValue(solutionFile);
481  parameters_[whichParam(SAVESOL,numberParameters_,parameters_)].setStringValue(solutionSaveFile);
482  parameters_[whichParam(SPRINT,numberParameters_,parameters_)].setIntValue(doSprint);
483  parameters_[whichParam(SUBSTITUTION,numberParameters_,parameters_)].setIntValue(substitution);
484  parameters_[whichParam(DUALIZE,numberParameters_,parameters_)].setIntValue(dualize);
485  parameters_[whichParam(NUMBERBEFORE,numberParameters_,parameters_)].setIntValue(model_.numberBeforeTrust());
486  parameters_[whichParam(MAXNODES,numberParameters_,parameters_)].setIntValue(model_.getMaximumNodes());
487  parameters_[whichParam(STRONGBRANCHING,numberParameters_,parameters_)].setIntValue(model_.numberStrong());
488  parameters_[whichParam(INFEASIBILITYWEIGHT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcInfeasibilityWeight));
489  parameters_[whichParam(INTEGERTOLERANCE,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcIntegerTolerance));
490  parameters_[whichParam(INCREMENT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcCutoffIncrement));
491  parameters_[whichParam(TESTOSI,numberParameters_,parameters_)].setIntValue(testOsiParameters);
492  parameters_[whichParam(FPUMPTUNE,numberParameters_,parameters_)].setIntValue(1003);
493  initialPumpTune=1003;
494#ifdef CBC_THREAD
495  parameters_[whichParam(THREADS,numberParameters_,parameters_)].setIntValue(0);
496#endif
497  // Set up likely cut generators and defaults
498  parameters_[whichParam(PREPROCESS,numberParameters_,parameters_)].setCurrentOption("sos");
499  parameters_[whichParam(MIPOPTIONS,numberParameters_,parameters_)].setIntValue(1057);
500  parameters_[whichParam(CUTPASSINTREE,numberParameters_,parameters_)].setIntValue(1);
501  parameters_[whichParam(MOREMIPOPTIONS,numberParameters_,parameters_)].setIntValue(-1);
502  parameters_[whichParam(MAXHOTITS,numberParameters_,parameters_)].setIntValue(100);
503  parameters_[whichParam(CUTSSTRATEGY,numberParameters_,parameters_)].setCurrentOption("on");
504  parameters_[whichParam(HEURISTICSTRATEGY,numberParameters_,parameters_)].setCurrentOption("on");
505  parameters_[whichParam(NODESTRATEGY,numberParameters_,parameters_)].setCurrentOption("fewest");
506  parameters_[whichParam(GOMORYCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove");
507  parameters_[whichParam(PROBINGCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove");
508  parameters_[whichParam(KNAPSACKCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove");
509#ifdef ZERO_HALF_CUTS
510  parameters_[whichParam(ZEROHALFCUTS,numberParameters_,parameters_)].setCurrentOption("off");
511#endif
512  parameters_[whichParam(REDSPLITCUTS,numberParameters_,parameters_)].setCurrentOption("off");
513  parameters_[whichParam(CLIQUECUTS,numberParameters_,parameters_)].setCurrentOption("ifmove");
514  parameters_[whichParam(MIXEDCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove");
515  parameters_[whichParam(FLOWCUTS,numberParameters_,parameters_)].setCurrentOption("ifmove");
516  parameters_[whichParam(TWOMIRCUTS,numberParameters_,parameters_)].setCurrentOption("root");
517  parameters_[whichParam(LANDPCUTS,numberParameters_,parameters_)].setCurrentOption("off");
518  parameters_[whichParam(RESIDCUTS,numberParameters_,parameters_)].setCurrentOption("off");
519  parameters_[whichParam(ROUNDING,numberParameters_,parameters_)].setCurrentOption("on");
520  parameters_[whichParam(FPUMP,numberParameters_,parameters_)].setCurrentOption("on");
521  parameters_[whichParam(GREEDY,numberParameters_,parameters_)].setCurrentOption("on");
522  parameters_[whichParam(COMBINE,numberParameters_,parameters_)].setCurrentOption("on");
523  parameters_[whichParam(CROSSOVER2,numberParameters_,parameters_)].setCurrentOption("off");
524  parameters_[whichParam(PIVOTANDFIX,numberParameters_,parameters_)].setCurrentOption("off");
525  parameters_[whichParam(RANDROUND,numberParameters_,parameters_)].setCurrentOption("off");
526  parameters_[whichParam(NAIVE,numberParameters_,parameters_)].setCurrentOption("off");
527  parameters_[whichParam(RINS,numberParameters_,parameters_)].setCurrentOption("off");
528  parameters_[whichParam(DINS,numberParameters_,parameters_)].setCurrentOption("off");
529  parameters_[whichParam(RENS,numberParameters_,parameters_)].setCurrentOption("off");
530  parameters_[whichParam(LOCALTREE,numberParameters_,parameters_)].setCurrentOption("off");
531  parameters_[whichParam(COSTSTRATEGY,numberParameters_,parameters_)].setCurrentOption("off");
532  if (createSolver)
533    delete clpSolver;
534}
535void CbcSolver::fillValuesInSolver()
536{
537  OsiSolverInterface * solver = model_.solver();
538  OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
539  assert (clpSolver);
540  noPrinting_ = (clpSolver->getModelPtr()->logLevel()==0);
541  CoinMessageHandler * generalMessageHandler = clpSolver->messageHandler();
542  generalMessageHandler->setPrefix(true);
543  ClpSimplex * lpSolver = clpSolver->getModelPtr();
544  lpSolver->setPerturbation(50);
545  lpSolver->messageHandler()->setPrefix(false);
546  parameters_[whichParam(DUALBOUND,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualBound());
547  parameters_[whichParam(DUALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->dualTolerance());
548  int iParam = whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_);
549  int value=parameters_[iParam].intValue();
550  clpSolver->messageHandler()->setLogLevel(value) ;
551  lpSolver->setLogLevel(value);
552  iParam = whichParam(LOGLEVEL,numberParameters_,parameters_);
553  value=parameters_[iParam].intValue();
554  model_.messageHandler()->setLogLevel(value);
555  parameters_[whichParam(LOGLEVEL,numberParameters_,parameters_)].setIntValue(model_.logLevel());
556  parameters_[whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_)].setIntValue(lpSolver->logLevel());
557  parameters_[whichParam(MAXFACTOR,numberParameters_,parameters_)].setIntValue(lpSolver->factorizationFrequency());
558  parameters_[whichParam(MAXITERATION,numberParameters_,parameters_)].setIntValue(lpSolver->maximumIterations());
559  parameters_[whichParam(PERTVALUE,numberParameters_,parameters_)].setIntValue(lpSolver->perturbation());
560  parameters_[whichParam(PRIMALTOLERANCE,numberParameters_,parameters_)].setDoubleValue(lpSolver->primalTolerance());
561  parameters_[whichParam(PRIMALWEIGHT,numberParameters_,parameters_)].setDoubleValue(lpSolver->infeasibilityCost());
562  parameters_[whichParam(NUMBERBEFORE,numberParameters_,parameters_)].setIntValue(model_.numberBeforeTrust());
563  parameters_[whichParam(MAXNODES,numberParameters_,parameters_)].setIntValue(model_.getMaximumNodes());
564  parameters_[whichParam(STRONGBRANCHING,numberParameters_,parameters_)].setIntValue(model_.numberStrong());
565  parameters_[whichParam(INFEASIBILITYWEIGHT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcInfeasibilityWeight));
566  parameters_[whichParam(INTEGERTOLERANCE,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcIntegerTolerance));
567  parameters_[whichParam(INCREMENT,numberParameters_,parameters_)].setDoubleValue(model_.getDblParam(CbcModel::CbcCutoffIncrement));
568}
569// Add user function
570void 
571CbcSolver::addUserFunction(CbcUser * function)
572{
573  CbcUser ** temp = new CbcUser * [numberUserFunctions_+1];
574  int i;
575  for (i=0;i<numberUserFunctions_;i++)
576    temp[i]=userFunction_[i];
577  delete [] userFunction_;
578  userFunction_ = temp;
579  userFunction_[numberUserFunctions_++] = function->clone();
580  delete [] statusUserFunction_;
581  statusUserFunction_ = NULL;
582}
583// Set user call back
584void 
585CbcSolver::setUserCallBack(CbcStopNow * function)
586{
587  delete callBack_;
588  callBack_ = function->clone();
589}
590// Copy of model on initial load (will contain output solutions)
591void 
592CbcSolver::setOriginalSolver(OsiClpSolverInterface * originalSolver)
593{
594  delete originalSolver_;
595  OsiSolverInterface * temp = originalSolver->clone();
596  originalSolver_ = dynamic_cast<OsiClpSolverInterface *> (temp);
597  assert (originalSolver_);
598 
599}
600// Copy of model on initial load
601void 
602CbcSolver::setOriginalCoinModel(CoinModel * originalCoinModel)
603{
604  delete originalCoinModel_;
605  originalCoinModel_ = new CoinModel(*originalCoinModel);
606}
607// Add cut generator
608void 
609CbcSolver::addCutGenerator(CglCutGenerator * generator)
610{
611  CglCutGenerator ** temp = new CglCutGenerator * [numberCutGenerators_+1];
612  int i;
613  for (i=0;i<numberCutGenerators_;i++)
614    temp[i]=cutGenerator_[i];
615  delete [] cutGenerator_;
616  cutGenerator_ = temp;
617  cutGenerator_[numberCutGenerators_++] = generator->clone();
618}
619// User stuff (base class)
620CbcUser::CbcUser()
621  : coinModel_(NULL),
622    userName_("null")
623{
624}
625CbcUser::~CbcUser()
626{
627  delete coinModel_;
628}
629// Copy constructor
630CbcUser::CbcUser ( const CbcUser & rhs)
631{
632  if (rhs.coinModel_)
633    coinModel_ = new CoinModel(*rhs.coinModel_);
634  else
635    coinModel_ = NULL;
636  userName_ = rhs.userName_;
637}
638// Assignment operator
639CbcUser &
640CbcUser::operator=(const CbcUser & rhs)
641{
642  if (this != &rhs) {
643    if (rhs.coinModel_)
644      coinModel_ = new CoinModel(*rhs.coinModel_);
645    else
646      coinModel_ = NULL;
647    userName_ = rhs.userName_;
648  }
649  return *this;
650}
651/* Updates model_ from babModel_ according to returnMode
652   returnMode -
653   0 model and solver untouched - babModel updated
654   1 model updated - just with solution basis etc
655   2 model updated i.e. as babModel (babModel NULL) (only use without preprocessing!)
656*/
657void
658CbcSolver::updateModel(ClpSimplex * model2, int returnMode)
659{
660  if (!returnMode)
661    return;
662  if (returnMode==2&&babModel_)
663    model_ = *babModel_;
664  if (model2) {
665    // Only continuous valid
666    // update with basis etc
667    // map states
668    /* clp status
669       -1 - unknown e.g. before solve or if postSolve says not optimal
670       0 - optimal
671       1 - primal infeasible
672       2 - dual infeasible
673       3 - stopped on iterations or time
674       4 - stopped due to errors
675       5 - stopped by event handler (virtual int ClpEventHandler::event()) */
676    /* cbc status
677       -1 before branchAndBound
678       0 finished - check isProvenOptimal or isProvenInfeasible to see if solution found
679       (or check value of best solution)
680       1 stopped - on maxnodes, maxsols, maxtime
681       2 difficulties so run was abandoned
682       (5 event user programmed event occurred) */
683    /* clp secondary status of problem - may get extended
684       0 - none
685       1 - primal infeasible because dual limit reached OR probably primal
686       infeasible but can't prove it (main status 4)
687       2 - scaled problem optimal - unscaled problem has primal infeasibilities
688       3 - scaled problem optimal - unscaled problem has dual infeasibilities
689       4 - scaled problem optimal - unscaled problem has primal and dual infeasibilities
690       5 - giving up in primal with flagged variables
691       6 - failed due to empty problem check
692       7 - postSolve says not optimal
693       8 - failed due to bad element check
694       9 - status was 3 and stopped on time
695       100 up - translation of enum from ClpEventHandler
696    */
697    /* cbc secondary status of problem
698       -1 unset (status_ will also be -1)
699       0 search completed with solution
700       1 linear relaxation not feasible (or worse than cutoff)
701       2 stopped on gap
702       3 stopped on nodes
703       4 stopped on time
704       5 stopped on user event
705       6 stopped on solutions
706       7 linear relaxation unbounded
707    */
708    int iStatus = model2->status();
709    int iStatus2 = model2->secondaryStatus();
710    if (iStatus==0) {
711      iStatus2=0;
712    } else if (iStatus==1) {
713      iStatus=0;
714      iStatus2=1; // say infeasible
715    } else if (iStatus==2) {
716      iStatus=0;
717      iStatus2=7; // say unbounded
718    } else if (iStatus==3) {
719      iStatus=1;
720      if (iStatus2==9)
721        iStatus2=4;
722      else
723        iStatus2=3; // Use nodes - as closer than solutions
724    } else if (iStatus==4) {
725      iStatus=2; // difficulties
726      iStatus2=0; 
727    }
728    model_.setProblemStatus(iStatus);
729    model_.setSecondaryStatus(iStatus2);
730    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model_.solver());
731    ClpSimplex * lpSolver = clpSolver->getModelPtr();
732    if (model2!=lpSolver) {
733      lpSolver->moveInfo(*model2);
734    }
735    clpSolver->setWarmStart(NULL); // synchronize bases
736    if (originalSolver_) {
737      ClpSimplex * lpSolver2 = originalSolver_->getModelPtr();
738      assert (model2!=lpSolver2);
739      lpSolver2->moveInfo(*model2);
740      originalSolver_->setWarmStart(NULL); // synchronize bases
741    }
742  } else if (returnMode==1) {
743    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (model_.solver());
744    ClpSimplex * lpSolver = clpSolver->getModelPtr();
745    if (babModel_) {
746      model_.moveInfo(*babModel_);
747      int numberColumns = babModel_->getNumCols();
748      if (babModel_->bestSolution())
749        model_.setBestSolution(babModel_->bestSolution(),numberColumns,babModel_->getMinimizationObjValue());
750      OsiClpSolverInterface * clpSolver1 = dynamic_cast< OsiClpSolverInterface*> (babModel_->solver());
751      ClpSimplex * lpSolver1 = clpSolver1->getModelPtr();
752      if (lpSolver1!=lpSolver&&model_.bestSolution()) {
753        lpSolver->moveInfo(*lpSolver1);
754      }
755    }
756    clpSolver->setWarmStart(NULL); // synchronize bases
757  }
758  if (returnMode==2) {
759    delete babModel_;
760    babModel_=NULL;
761  }
762}
763// Stop now stuff (base class)
764CbcStopNow::CbcStopNow()
765{
766}
767CbcStopNow::~CbcStopNow()
768{
769}
770// Copy constructor
771CbcStopNow::CbcStopNow ( const CbcStopNow & )
772{
773}
774// Assignment operator
775CbcStopNow &
776CbcStopNow::operator=(const CbcStopNow & rhs)
777{
778  if (this != &rhs) {
779  }
780  return *this;
781}
782// Clone
783CbcStopNow *
784CbcStopNow::clone() const
785{
786  return new CbcStopNow(*this);
787}
788#ifdef CLP_FAST_CODE
789// force new style solver
790#ifndef NEW_STYLE_SOLVER
791#define NEW_STYLE_SOLVER 1
792#endif
793#else
794// Not new style solver
795#ifndef NEW_STYLE_SOLVER
796#define NEW_STYLE_SOLVER 0
797#endif
798#endif
799#ifdef CPX_KEEP_RESULTS
800#define CBC_OTHER_SOLVER 1
801#endif
802#ifdef COIN_HAS_CPX
803#include "OsiCpxSolverInterface.hpp"
804#endif
805#ifdef CBC_OTHER_SOLVER
806#if CBC_OTHER_SOLVER==1
807#include "OsiCpxSolverInterface.hpp"
808#endif
809#undef NEW_STYLE_SOLVER
810#define NEW_STYLE_SOLVER 0
811#endif
812//#undef COIN_HAS_ASL
813#ifdef COIN_HAS_ASL
814#include "Cbc_ampl.h"
815#endif
816static double totalTime=0.0;
817static void statistics(ClpSimplex * originalModel, ClpSimplex * model);
818static bool maskMatches(const int * starts, char ** masks,
819                        std::string & check);
820static void generateCode(CbcModel * model, const char * fileName,int type,int preProcess);
821//#ifdef NDEBUG
822//#undef NDEBUG
823//#endif
824// define (probably dummy fake main programs for UserClp and UserCbc
825void fakeMain (ClpSimplex & model,OsiSolverInterface & osiSolver, CbcModel & babSolver);
826void fakeMain2 (ClpSimplex & model,OsiClpSolverInterface & osiSolver,int options);
827
828// Allow for interrupts
829// But is this threadsafe ? (so switched off by option)
830
831#include "CoinSignal.hpp"
832static CbcModel * currentBranchModel = NULL;
833
834extern "C" {
835  static void signal_handler(int /*whichSignal*/)
836   {
837     if (currentBranchModel!=NULL) {
838       currentBranchModel->setMaximumNodes(0); // stop at next node
839       currentBranchModel->setMaximumSeconds(0.0); // stop
840     }
841     return;
842   }
843}
844//#define CBC_SIG_TRAP
845#ifdef CBC_SIG_TRAP
846#include <setjmp.h>
847static sigjmp_buf cbc_seg_buffer;
848extern "C" {
849   static void signal_handler_error(int whichSignal)
850   {
851     siglongjmp(cbc_seg_buffer,1);
852   }
853}
854#endif
855#if 0
856/* Updates model_ from babModel_ according to returnMode
857   returnMode -
858   0 model and solver untouched
859   1 model updated - just with solution basis etc
860*/
861static void updateModel(CbcModel & model_,int returnMode)
862{
863  if (!returnMode)
864    return;
865  assert (returnMode==1);
866}
867#endif
868int CbcOrClpRead_mode=1;
869FILE * CbcOrClpReadCommand=stdin;
870extern int CbcOrClpEnvironmentIndex;
871static bool noPrinting=false;
872#ifndef CBC_OTHER_SOLVER
873#if NEW_STYLE_SOLVER
874int * CbcSolver::analyze(OsiClpSolverInterface * solverMod, int & numberChanged, double & increment,
875                     bool changeInt,  CoinMessageHandler * generalMessageHandler)
876#else
877static int * analyze(OsiClpSolverInterface * solverMod, int & numberChanged, double & increment,
878                     bool changeInt,  CoinMessageHandler * generalMessageHandler)
879#endif
880{
881#if NEW_STYLE_SOLVER==0
882  bool noPrinting_ = noPrinting;
883#endif
884  OsiSolverInterface * solver = solverMod->clone();
885  char generalPrint[200];
886  if (0) {
887    // just get increment
888    CbcModel model(*solver);
889    model.analyzeObjective();
890    double increment2=model.getCutoffIncrement();
891    printf("initial cutoff increment %g\n",increment2);
892  }
893  const double *objective = solver->getObjCoefficients() ;
894  const double *lower = solver->getColLower() ;
895  const double *upper = solver->getColUpper() ;
896  int numberColumns = solver->getNumCols() ;
897  int numberRows = solver->getNumRows();
898  double direction = solver->getObjSense();
899  int iRow,iColumn;
900
901  // Row copy
902  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
903  const double * elementByRow = matrixByRow.getElements();
904  const int * column = matrixByRow.getIndices();
905  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
906  const int * rowLength = matrixByRow.getVectorLengths();
907
908  // Column copy
909  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
910  const double * element = matrixByCol.getElements();
911  const int * row = matrixByCol.getIndices();
912  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
913  const int * columnLength = matrixByCol.getVectorLengths();
914
915  const double * rowLower = solver->getRowLower();
916  const double * rowUpper = solver->getRowUpper();
917
918  char * ignore = new char [numberRows];
919  int * changed = new int[numberColumns];
920  int * which = new int[numberRows];
921  double * changeRhs = new double[numberRows];
922  memset(changeRhs,0,numberRows*sizeof(double));
923  memset(ignore,0,numberRows);
924  numberChanged=0;
925  int numberInteger=0;
926  for (iColumn=0;iColumn<numberColumns;iColumn++) {
927    if (upper[iColumn] > lower[iColumn]+1.0e-8&&solver->isInteger(iColumn)) 
928      numberInteger++;
929  }
930  bool finished=false;
931  while (!finished) {
932    int saveNumberChanged = numberChanged;
933    for (iRow=0;iRow<numberRows;iRow++) {
934      int numberContinuous=0;
935      double value1=0.0,value2=0.0;
936      bool allIntegerCoeff=true;
937      double sumFixed=0.0;
938      int jColumn1=-1,jColumn2=-1;
939      for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
940        int jColumn = column[j];
941        double value = elementByRow[j];
942        if (upper[jColumn] > lower[jColumn]+1.0e-8) {
943          if (!solver->isInteger(jColumn)) {
944            if (numberContinuous==0) {
945              jColumn1=jColumn;
946              value1=value;
947            } else {
948              jColumn2=jColumn;
949              value2=value;
950            }
951            numberContinuous++;
952          } else {
953            if (fabs(value-floor(value+0.5))>1.0e-12)
954              allIntegerCoeff=false;
955          }
956        } else {
957          sumFixed += lower[jColumn]*value;
958        }
959      }
960      double low = rowLower[iRow];
961      if (low>-1.0e20) {
962        low -= sumFixed;
963        if (fabs(low-floor(low+0.5))>1.0e-12)
964          allIntegerCoeff=false;
965      }
966      double up = rowUpper[iRow];
967      if (up<1.0e20) {
968        up -= sumFixed;
969        if (fabs(up-floor(up+0.5))>1.0e-12)
970          allIntegerCoeff=false;
971      }
972      if (!allIntegerCoeff)
973        continue; // can't do
974      if (numberContinuous==1) {
975        // see if really integer
976        // This does not allow for complicated cases
977        if (low==up) {
978          if (fabs(value1)>1.0e-3) {
979            value1 = 1.0/value1;
980            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
981              // integer
982              changed[numberChanged++]=jColumn1;
983              solver->setInteger(jColumn1);
984              if (upper[jColumn1]>1.0e20)
985                solver->setColUpper(jColumn1,1.0e20);
986              if (lower[jColumn1]<-1.0e20)
987                solver->setColLower(jColumn1,-1.0e20);
988            }
989          }
990        } else {
991          if (fabs(value1)>1.0e-3) {
992            value1 = 1.0/value1;
993            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
994              // This constraint will not stop it being integer
995              ignore[iRow]=1;
996            }
997          }
998        }
999      } else if (numberContinuous==2) {
1000        if (low==up) {
1001          /* need general theory - for now just look at 2 cases -
1002             1 - +- 1 one in column and just costs i.e. matching objective
1003             2 - +- 1 two in column but feeds into G/L row which will try and minimize
1004          */
1005          if (fabs(value1)==1.0&&value1*value2==-1.0&&!lower[jColumn1]
1006              &&!lower[jColumn2]) {
1007            int n=0;
1008            int i;
1009            double objChange=direction*(objective[jColumn1]+objective[jColumn2]);
1010            double bound = CoinMin(upper[jColumn1],upper[jColumn2]);
1011            bound = CoinMin(bound,1.0e20);
1012            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
1013              int jRow = row[i];
1014              double value = element[i];
1015              if (jRow!=iRow) {
1016                which[n++]=jRow;
1017                changeRhs[jRow]=value;
1018              }
1019            }
1020            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
1021              int jRow = row[i];
1022              double value = element[i];
1023              if (jRow!=iRow) {
1024                if (!changeRhs[jRow]) {
1025                  which[n++]=jRow;
1026                  changeRhs[jRow]=value;
1027                } else {
1028                  changeRhs[jRow]+=value;
1029                }
1030              }
1031            }
1032            if (objChange>=0.0) {
1033              // see if all rows OK
1034              bool good=true;
1035              for (i=0;i<n;i++) {
1036                int jRow = which[i];
1037                double value = changeRhs[jRow];
1038                if (value) {
1039                  value *= bound;
1040                  if (rowLength[jRow]==1) {
1041                    if (value>0.0) {
1042                      double rhs = rowLower[jRow];
1043                      if (rhs>0.0) {
1044                        double ratio =rhs/value;
1045                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
1046                          good=false;
1047                      }
1048                    } else {
1049                      double rhs = rowUpper[jRow];
1050                      if (rhs<0.0) {
1051                        double ratio =rhs/value;
1052                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
1053                          good=false;
1054                      }
1055                    }
1056                  } else if (rowLength[jRow]==2) {
1057                    if (value>0.0) {
1058                      if (rowLower[jRow]>-1.0e20)
1059                        good=false;
1060                    } else {
1061                      if (rowUpper[jRow]<1.0e20)
1062                        good=false;
1063                    }
1064                  } else {
1065                    good=false;
1066                  }
1067                }
1068              }
1069              if (good) {
1070                // both can be integer
1071                changed[numberChanged++]=jColumn1;
1072                solver->setInteger(jColumn1);
1073                if (upper[jColumn1]>1.0e20)
1074                  solver->setColUpper(jColumn1,1.0e20);
1075                if (lower[jColumn1]<-1.0e20)
1076                  solver->setColLower(jColumn1,-1.0e20);
1077                changed[numberChanged++]=jColumn2;
1078                solver->setInteger(jColumn2);
1079                if (upper[jColumn2]>1.0e20)
1080                  solver->setColUpper(jColumn2,1.0e20);
1081                if (lower[jColumn2]<-1.0e20)
1082                  solver->setColLower(jColumn2,-1.0e20);
1083              }
1084            }
1085            // clear
1086            for (i=0;i<n;i++) {
1087              changeRhs[which[i]]=0.0;
1088            }
1089          }
1090        }
1091      }
1092    }
1093    for (iColumn=0;iColumn<numberColumns;iColumn++) {
1094      if (upper[iColumn] > lower[iColumn]+1.0e-8&&!solver->isInteger(iColumn)) {
1095        double value;
1096        value = upper[iColumn];
1097        if (value<1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
1098          continue;
1099        value = lower[iColumn];
1100        if (value>-1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
1101          continue;
1102        bool integer=true;
1103        for (CoinBigIndex j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
1104          int iRow = row[j];
1105          if (!ignore[iRow]) {
1106            integer=false;
1107            break;
1108          }
1109        }
1110        if (integer) {
1111          // integer
1112          changed[numberChanged++]=iColumn;
1113          solver->setInteger(iColumn);
1114          if (upper[iColumn]>1.0e20)
1115            solver->setColUpper(iColumn,1.0e20);
1116          if (lower[iColumn]<-1.0e20)
1117            solver->setColLower(iColumn,-1.0e20);
1118        }
1119      }
1120    }
1121    finished = numberChanged==saveNumberChanged;
1122  }
1123  delete [] which;
1124  delete [] changeRhs;
1125  delete [] ignore;
1126  //if (numberInteger&&!noPrinting_)
1127  //printf("%d integer variables",numberInteger);
1128  if (changeInt) {
1129    //if (!noPrinting_) {
1130    //if (numberChanged)
1131    //  printf(" and %d variables made integer\n",numberChanged);
1132    //else
1133    //  printf("\n");
1134    //}
1135    delete [] ignore;
1136    //increment=0.0;
1137    if (!numberChanged) {
1138      delete [] changed;
1139      delete solver;
1140      return NULL;
1141    } else {
1142      for (iColumn=0;iColumn<numberColumns;iColumn++) {
1143        if (solver->isInteger(iColumn))
1144          solverMod->setInteger(iColumn);
1145      }
1146      delete solver;
1147      return changed;
1148    }
1149  } else {
1150    //if (!noPrinting_) {
1151    //if (numberChanged)
1152    //  printf(" and %d variables could be made integer\n",numberChanged);
1153    //else
1154    //  printf("\n");
1155    //}
1156    // just get increment
1157    int logLevel=generalMessageHandler->logLevel();
1158    CbcModel model(*solver);
1159    model.passInMessageHandler(generalMessageHandler);
1160    if (noPrinting_)
1161      model.setLogLevel(0);
1162    model.analyzeObjective();
1163    generalMessageHandler->setLogLevel(logLevel);
1164    double increment2=model.getCutoffIncrement();
1165    if (increment2>increment&&increment2>0.0) {
1166      if (!noPrinting_) {
1167        sprintf(generalPrint,"Cutoff increment increased from %g to %g",increment,increment2);
1168        CoinMessages generalMessages = solverMod->getModelPtr()->messages();
1169        generalMessageHandler->message(CLP_GENERAL,generalMessages)
1170          << generalPrint
1171          <<CoinMessageEol;
1172      }
1173      increment=increment2;
1174    }
1175    delete solver;
1176    numberChanged=0;
1177    delete [] changed;
1178    return NULL;
1179  }
1180}
1181#endif
1182#if 1
1183#include "ClpSimplexOther.hpp"
1184
1185// Crunch down model
1186static void 
1187crunchIt(ClpSimplex * model)
1188{
1189#if 0
1190  model->dual();
1191#else
1192  int numberColumns = model->numberColumns();
1193  int numberRows = model->numberRows();
1194  // Use dual region
1195  double * rhs = model->dualRowSolution();
1196  int * whichRow = new int[3*numberRows];
1197  int * whichColumn = new int[2*numberColumns];
1198  int nBound;
1199  ClpSimplex * small = static_cast<ClpSimplexOther *> (model)->crunch(rhs,whichRow,whichColumn,
1200                                                               nBound,false,false);
1201  if (small) {
1202    small->dual();
1203    if (small->problemStatus()==0) {
1204      model->setProblemStatus(0);
1205      static_cast<ClpSimplexOther *> (model)->afterCrunch(*small,whichRow,whichColumn,nBound);
1206    } else if (small->problemStatus()!=3) {
1207      model->setProblemStatus(1);
1208    } else {
1209      if (small->problemStatus()==3) {
1210        // may be problems
1211        small->computeObjectiveValue();
1212        model->setObjectiveValue(small->objectiveValue());
1213        model->setProblemStatus(3);
1214      } else {
1215        model->setProblemStatus(3);
1216      }
1217    }
1218    delete small;
1219  } else {
1220    model->setProblemStatus(1);
1221  }
1222  delete [] whichRow;
1223  delete [] whichColumn;
1224#endif
1225}
1226/*
1227  On input
1228  doAction - 0 just fix in original and return NULL
1229             1 return fixed non-presolved solver
1230             2 as one but use presolve Inside this
1231             3 use presolve and fix ones with large cost
1232             ? do heuristics and set best solution
1233             ? do BAB and just set best solution
1234             10+ then use lastSolution and relax a few
1235             -2 cleanup afterwards if using 2
1236  On output - number fixed
1237*/
1238static OsiClpSolverInterface * 
1239fixVubs(CbcModel & model, int skipZero2,
1240        int & doAction, 
1241        CoinMessageHandler * /*generalMessageHandler*/,
1242        const double * lastSolution, double dextra[6],
1243        int extra[5])
1244{
1245  if (doAction==11&&!lastSolution)
1246    lastSolution = model.bestSolution();
1247  assert (((doAction>=0&&doAction<=3)&&!lastSolution)||(doAction==11&&lastSolution));
1248  double fractionIntFixed = dextra[3];
1249  double fractionFixed = dextra[4];
1250  double fixAbove = dextra[2];
1251  double fixAboveValue = (dextra[5]>0.0) ? dextra[5] : 1.0;
1252  double time1 = CoinCpuTime();
1253  int leaveIntFree = extra[1];
1254  OsiSolverInterface * originalSolver = model.solver();
1255  OsiClpSolverInterface * originalClpSolver = dynamic_cast< OsiClpSolverInterface*> (originalSolver);
1256  ClpSimplex * originalLpSolver = originalClpSolver->getModelPtr();
1257  int * originalColumns=NULL;
1258  OsiClpSolverInterface * clpSolver;
1259  ClpSimplex * lpSolver;
1260  ClpPresolve pinfo;
1261  assert(originalSolver->getObjSense()>0);
1262  if (doAction==2||doAction==3) {
1263    double * saveLB=NULL;
1264    double * saveUB=NULL;
1265    int numberColumns = originalLpSolver->numberColumns();
1266    if (fixAbove>0.0) {
1267      double time1 = CoinCpuTime();
1268      originalClpSolver->initialSolve();
1269      printf("first solve took %g seconds\n",CoinCpuTime()-time1);
1270      double * columnLower = originalLpSolver->columnLower() ;
1271      double * columnUpper = originalLpSolver->columnUpper() ;
1272      const double * solution = originalLpSolver->primalColumnSolution();
1273      saveLB = CoinCopyOfArray(columnLower,numberColumns);
1274      saveUB = CoinCopyOfArray(columnUpper,numberColumns);
1275      const double * objective = originalLpSolver->getObjCoefficients() ;
1276      int iColumn;
1277      int nFix=0;
1278      int nArt=0;
1279      for (iColumn=0;iColumn<numberColumns;iColumn++) {
1280        if (objective[iColumn]>fixAbove) {
1281          if (solution[iColumn]<columnLower[iColumn]+1.0e-8) {
1282            columnUpper[iColumn]=columnLower[iColumn];
1283            nFix++;
1284          } else {
1285            nArt++;
1286          }
1287        } else if (objective[iColumn]<-fixAbove) {
1288          if (solution[iColumn]>columnUpper[iColumn]-1.0e-8) {
1289            columnLower[iColumn]=columnUpper[iColumn];
1290            nFix++;
1291          } else {
1292            nArt++;
1293          }
1294        }
1295      }
1296      printf("%d artificials fixed, %d left as in solution\n",nFix,nArt);
1297      lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
1298      if (!lpSolver||doAction==2) {
1299        // take off fixing in original
1300        memcpy(columnLower,saveLB,numberColumns*sizeof(double));
1301        memcpy(columnUpper,saveUB,numberColumns*sizeof(double));
1302      }
1303      delete [] saveLB;
1304      delete [] saveUB;
1305      if (!lpSolver) {
1306        // try again
1307        pinfo.destroyPresolve();
1308        lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
1309        assert (lpSolver);
1310      }
1311    } else {
1312      lpSolver = pinfo.presolvedModel(*originalLpSolver,1.0e-8,true,10);
1313      assert (lpSolver);
1314    }
1315    clpSolver = new OsiClpSolverInterface(lpSolver,true);
1316    assert(lpSolver == clpSolver->getModelPtr());
1317    numberColumns = lpSolver->numberColumns();
1318    originalColumns = CoinCopyOfArray(pinfo.originalColumns(),numberColumns);
1319    doAction=1;
1320  } else {
1321    OsiSolverInterface * solver = originalSolver->clone();
1322    clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
1323    lpSolver = clpSolver->getModelPtr();
1324  }
1325  // Tighten bounds
1326  lpSolver->tightenPrimalBounds(0.0,11,true);
1327  int numberColumns = clpSolver->getNumCols() ;
1328  double * saveColumnLower = CoinCopyOfArray(lpSolver->columnLower(),numberColumns);
1329  double * saveColumnUpper = CoinCopyOfArray(lpSolver->columnUpper(),numberColumns);
1330  //char generalPrint[200];
1331  const double *objective = lpSolver->getObjCoefficients() ;
1332  double *columnLower = lpSolver->columnLower() ;
1333  double *columnUpper = lpSolver->columnUpper() ;
1334  int numberRows = clpSolver->getNumRows();
1335  int iRow,iColumn;
1336
1337  // Row copy
1338  CoinPackedMatrix matrixByRow(*clpSolver->getMatrixByRow());
1339  const double * elementByRow = matrixByRow.getElements();
1340  const int * column = matrixByRow.getIndices();
1341  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
1342  const int * rowLength = matrixByRow.getVectorLengths();
1343
1344  // Column copy
1345  CoinPackedMatrix  matrixByCol(*clpSolver->getMatrixByCol());
1346  //const double * element = matrixByCol.getElements();
1347  const int * row = matrixByCol.getIndices();
1348  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
1349  const int * columnLength = matrixByCol.getVectorLengths();
1350
1351  const double * rowLower = clpSolver->getRowLower();
1352  const double * rowUpper = clpSolver->getRowUpper();
1353
1354  // Get maximum size of VUB tree
1355  // otherColumn is one fixed to 0 if this one zero
1356  int nEl = matrixByCol.getNumElements();
1357  CoinBigIndex * fixColumn = new CoinBigIndex [numberColumns+1];
1358  int * otherColumn = new int [nEl];
1359  int * fix = new int[numberColumns];
1360  char * mark = new char [numberColumns];
1361  memset(mark,0,numberColumns);
1362  int numberInteger=0;
1363  int numberOther=0;
1364  fixColumn[0]=0;
1365  double large = lpSolver->largeValue(); // treat bounds > this as infinite
1366#ifndef NDEBUG
1367  double large2= 1.0e10*large;
1368#endif
1369  double tolerance = lpSolver->primalTolerance();
1370  int * check = new int[numberRows];
1371  for (iRow=0;iRow<numberRows;iRow++) {
1372    check[iRow]=-2; // don't check
1373    if (rowLower[iRow]<-1.0e6&&rowUpper[iRow]>1.0e6) 
1374      continue;// unlikely
1375    // possible row
1376    int numberPositive=0;
1377    int iPositive=-1;
1378    int numberNegative=0;
1379    int iNegative=-1;
1380    CoinBigIndex rStart = rowStart[iRow];
1381    CoinBigIndex rEnd = rowStart[iRow]+rowLength[iRow];
1382    CoinBigIndex j;
1383    int kColumn;
1384    for (j = rStart; j < rEnd; ++j) {
1385      double value=elementByRow[j];
1386      kColumn = column[j];
1387      if (columnUpper[kColumn]>columnLower[kColumn]) {
1388        if (value > 0.0) { 
1389          numberPositive++;
1390          iPositive=kColumn;
1391        } else {
1392          numberNegative++;
1393          iNegative=kColumn;
1394        }
1395      }
1396    }
1397    if (numberPositive==1&&numberNegative==1)
1398      check[iRow]=-1; // try both
1399    if (numberPositive==1&&rowLower[iRow]>-1.0e20)
1400      check[iRow]=iPositive;
1401    else if (numberNegative==1&&rowUpper[iRow]<1.0e20)
1402      check[iRow]=iNegative;
1403  }
1404  for (iColumn=0;iColumn<numberColumns;iColumn++) {
1405    fix[iColumn]=-1;
1406    if (columnUpper[iColumn] > columnLower[iColumn]+1.0e-8) {
1407      if (clpSolver->isInteger(iColumn))
1408        numberInteger++;
1409      if (columnLower[iColumn]==0.0) {
1410        bool infeasible=false;
1411        fix[iColumn]=0;
1412        // fake upper bound
1413        double saveUpper = columnUpper[iColumn];
1414        columnUpper[iColumn]=0.0;
1415        for (CoinBigIndex i=columnStart[iColumn];
1416             i<columnStart[iColumn]+columnLength[iColumn];i++) {
1417          iRow = row[i];
1418          if (check[iRow]!=-1&&check[iRow]!=iColumn)
1419            continue; // unlikely
1420          // possible row
1421          int infiniteUpper = 0;
1422          int infiniteLower = 0;
1423          double maximumUp = 0.0;
1424          double maximumDown = 0.0;
1425          double newBound;
1426          CoinBigIndex rStart = rowStart[iRow];
1427          CoinBigIndex rEnd = rowStart[iRow]+rowLength[iRow];
1428          CoinBigIndex j;
1429          int kColumn;
1430          // Compute possible lower and upper ranges
1431          for (j = rStart; j < rEnd; ++j) {
1432            double value=elementByRow[j];
1433            kColumn = column[j];
1434            if (value > 0.0) {
1435              if (columnUpper[kColumn] >= large) {
1436                ++infiniteUpper;
1437              } else {
1438                maximumUp += columnUpper[kColumn] * value;
1439              }
1440              if (columnLower[kColumn] <= -large) {
1441                ++infiniteLower;
1442              } else {
1443                maximumDown += columnLower[kColumn] * value;
1444              }
1445            } else if (value<0.0) {
1446              if (columnUpper[kColumn] >= large) {
1447                ++infiniteLower;
1448              } else {
1449                maximumDown += columnUpper[kColumn] * value;
1450              }
1451              if (columnLower[kColumn] <= -large) {
1452                ++infiniteUpper;
1453              } else {
1454                maximumUp += columnLower[kColumn] * value;
1455              }
1456            }
1457          }
1458          // Build in a margin of error
1459          maximumUp += 1.0e-8*fabs(maximumUp);
1460          maximumDown -= 1.0e-8*fabs(maximumDown);
1461          double maxUp = maximumUp+infiniteUpper*1.0e31;
1462          double maxDown = maximumDown-infiniteLower*1.0e31;
1463          if (maxUp <= rowUpper[iRow] + tolerance && 
1464              maxDown >= rowLower[iRow] - tolerance) {
1465            //printf("Redundant row in vubs %d\n",iRow);
1466          } else {
1467            if (maxUp < rowLower[iRow] -100.0*tolerance ||
1468                maxDown > rowUpper[iRow]+100.0*tolerance) {
1469              infeasible=true;
1470              break;
1471            }
1472            double lower = rowLower[iRow];
1473            double upper = rowUpper[iRow];
1474            for (j = rStart; j < rEnd; ++j) {
1475              double value=elementByRow[j];
1476              kColumn = column[j];
1477              double nowLower = columnLower[kColumn];
1478              double nowUpper = columnUpper[kColumn];
1479              if (value > 0.0) {
1480                // positive value
1481                if (lower>-large) {
1482                  if (!infiniteUpper) {
1483                    assert(nowUpper < large2);
1484                    newBound = nowUpper + 
1485                      (lower - maximumUp) / value;
1486                    // relax if original was large
1487                    if (fabs(maximumUp)>1.0e8)
1488                      newBound -= 1.0e-12*fabs(maximumUp);
1489                  } else if (infiniteUpper==1&&nowUpper>large) {
1490                    newBound = (lower -maximumUp) / value;
1491                    // relax if original was large
1492                    if (fabs(maximumUp)>1.0e8)
1493                      newBound -= 1.0e-12*fabs(maximumUp);
1494                  } else {
1495                    newBound = -COIN_DBL_MAX;
1496                  }
1497                  if (newBound > nowLower + 1.0e-12&&newBound>-large) {
1498                    // Tighten the lower bound
1499                    // check infeasible (relaxed)
1500                    if (nowUpper < newBound) { 
1501                      if (nowUpper - newBound < 
1502                          -100.0*tolerance) { 
1503                        infeasible=true;
1504                        break;
1505                      }
1506                    }
1507                  }
1508                } 
1509                if (upper <large) {
1510                  if (!infiniteLower) {
1511                    assert(nowLower >- large2);
1512                    newBound = nowLower + 
1513                      (upper - maximumDown) / value;
1514                    // relax if original was large
1515                    if (fabs(maximumDown)>1.0e8)
1516                    newBound += 1.0e-12*fabs(maximumDown);
1517                  } else if (infiniteLower==1&&nowLower<-large) {
1518                    newBound =   (upper - maximumDown) / value;
1519                    // relax if original was large
1520                    if (fabs(maximumDown)>1.0e8)
1521                      newBound += 1.0e-12*fabs(maximumDown);
1522                  } else {
1523                    newBound = COIN_DBL_MAX;
1524                  }
1525                  if (newBound < nowUpper - 1.0e-12&&newBound<large) {
1526                    // Tighten the upper bound
1527                    // check infeasible (relaxed)
1528                    if (nowLower > newBound) { 
1529                      if (newBound - nowLower < 
1530                          -100.0*tolerance) {
1531                        infeasible=true;
1532                        break;
1533                      } else {
1534                        newBound=nowLower;
1535                      }
1536                    }
1537                    if (!newBound||(clpSolver->isInteger(kColumn)&&newBound<0.999)) {
1538                      // fix to zero
1539                      if (!mark[kColumn]) {
1540                        otherColumn[numberOther++]=kColumn;
1541                        mark[kColumn]=1;
1542                        if (check[iRow]==-1)
1543                          check[iRow]=iColumn;
1544                        else
1545                          assert(check[iRow]==iColumn);
1546                      }
1547                    }
1548                  }
1549                }
1550              } else {
1551                // negative value
1552                if (lower>-large) {
1553                  if (!infiniteUpper) {
1554                    assert(nowLower < large2);
1555                    newBound = nowLower + 
1556                      (lower - maximumUp) / value;
1557                    // relax if original was large
1558                    if (fabs(maximumUp)>1.0e8)
1559                      newBound += 1.0e-12*fabs(maximumUp);
1560                  } else if (infiniteUpper==1&&nowLower<-large) {
1561                    newBound = (lower -maximumUp) / value;
1562                    // relax if original was large
1563                    if (fabs(maximumUp)>1.0e8)
1564                      newBound += 1.0e-12*fabs(maximumUp);
1565                  } else {
1566                    newBound = COIN_DBL_MAX;
1567                  }
1568                  if (newBound < nowUpper - 1.0e-12&&newBound<large) {
1569                    // Tighten the upper bound
1570                    // check infeasible (relaxed)
1571                    if (nowLower > newBound) { 
1572                      if (newBound - nowLower < 
1573                          -100.0*tolerance) {
1574                        infeasible=true;
1575                        break;
1576                      } else { 
1577                        newBound=nowLower;
1578                      }
1579                    }
1580                    if (!newBound||(clpSolver->isInteger(kColumn)&&newBound<0.999)) {
1581                      // fix to zero
1582                      if (!mark[kColumn]) {
1583                        otherColumn[numberOther++]=kColumn;
1584                        mark[kColumn]=1;
1585                        if (check[iRow]==-1)
1586                          check[iRow]=iColumn;
1587                        else
1588                          assert(check[iRow]==iColumn);
1589                      }
1590                    }
1591                  }
1592                }
1593                if (upper <large) {
1594                  if (!infiniteLower) {
1595                    assert(nowUpper < large2);
1596                    newBound = nowUpper + 
1597                      (upper - maximumDown) / value;
1598                    // relax if original was large
1599                    if (fabs(maximumDown)>1.0e8)
1600                      newBound -= 1.0e-12*fabs(maximumDown);
1601                  } else if (infiniteLower==1&&nowUpper>large) {
1602                    newBound =   (upper - maximumDown) / value;
1603                    // relax if original was large
1604                    if (fabs(maximumDown)>1.0e8)
1605                      newBound -= 1.0e-12*fabs(maximumDown);
1606                  } else {
1607                    newBound = -COIN_DBL_MAX;
1608                  }
1609                  if (newBound > nowLower + 1.0e-12&&newBound>-large) {
1610                    // Tighten the lower bound
1611                    // check infeasible (relaxed)
1612                    if (nowUpper < newBound) { 
1613                      if (nowUpper - newBound < 
1614                          -100.0*tolerance) {
1615                        infeasible=true;
1616                        break;
1617                      }
1618                    }
1619                  }
1620                }
1621              }
1622            }
1623          }
1624        }
1625        for (int i=fixColumn[iColumn];i<numberOther;i++)
1626          mark[otherColumn[i]]=0;
1627        // reset bound unless infeasible
1628        if (!infeasible||!clpSolver->isInteger(iColumn))
1629          columnUpper[iColumn]=saveUpper;
1630        else if (clpSolver->isInteger(iColumn))
1631          columnLower[iColumn]=1.0;
1632      }
1633    }
1634    fixColumn[iColumn+1]=numberOther;
1635  }
1636  delete [] check;
1637  delete [] mark;
1638  // Now do reverse way
1639  int * counts = new int [numberColumns];
1640  CoinZeroN(counts,numberColumns);
1641  for (iColumn=0;iColumn<numberColumns;iColumn++) {
1642    for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++)
1643      counts[otherColumn[i]]++;
1644  }
1645  numberOther=0;
1646  CoinBigIndex * fixColumn2 = new CoinBigIndex [numberColumns+1];
1647  int * otherColumn2 = new int [fixColumn[numberColumns]];
1648  fixColumn2[0]=0;
1649  for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1650    numberOther += counts[iColumn];
1651    counts[iColumn]=0;
1652    fixColumn2[iColumn+1]=numberOther;
1653  }
1654  // Create other way
1655  for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1656    for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
1657      int jColumn=otherColumn[i];
1658      CoinBigIndex put = fixColumn2[jColumn] + counts[jColumn];
1659      counts[jColumn]++;
1660      otherColumn2[put]=iColumn;
1661    }
1662  }
1663  // get top layer i.e. those which are not fixed by any other
1664  int kLayer=0;
1665  while (true) {
1666    int numberLayered=0;
1667    for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1668      if (fix[iColumn]==kLayer) {
1669        for (int i=fixColumn2[iColumn];i<fixColumn2[iColumn+1];i++) {
1670          int jColumn=otherColumn2[i];
1671          if (fix[jColumn]==kLayer) {
1672            fix[iColumn]=kLayer+100;
1673          }
1674        }
1675      }
1676      if (fix[iColumn]==kLayer) {
1677        numberLayered++;
1678      }
1679    }
1680    if (numberLayered) {
1681      kLayer+=100;
1682    } else {
1683      break;
1684    }
1685  }
1686  for (int iPass=0;iPass<2;iPass++) {
1687    for (int jLayer=0;jLayer<kLayer;jLayer++) {
1688      int check[]={-1,0,1,2,3,4,5,10,50,100,500,1000,5000,10000,COIN_INT_MAX};
1689      int nCheck = static_cast<int> (sizeof(check)/sizeof(int));
1690      int countsI[20];
1691      int countsC[20];
1692      assert (nCheck<=20);
1693      memset(countsI,0,nCheck*sizeof(int));
1694      memset(countsC,0,nCheck*sizeof(int));
1695      check[nCheck-1]=numberColumns;
1696      int numberLayered=0;
1697      int numberInteger=0;
1698      for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1699        if (fix[iColumn]==jLayer) {
1700          numberLayered++;
1701          int nFix = fixColumn[iColumn+1]-fixColumn[iColumn];
1702          if (iPass) {
1703            // just integers
1704            nFix=0;
1705            for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
1706              int jColumn=otherColumn[i];
1707              if (clpSolver->isInteger(jColumn))
1708                nFix++;
1709            }
1710          }
1711          int iFix;
1712          for (iFix=0;iFix<nCheck;iFix++) {
1713            if (nFix<=check[iFix])
1714              break;
1715          }
1716          assert (iFix<nCheck);
1717          if (clpSolver->isInteger(iColumn)) {
1718            numberInteger++;
1719            countsI[iFix]++;
1720          } else {
1721            countsC[iFix]++;
1722          }
1723        }
1724      }
1725      if (numberLayered) {
1726        printf("%d (%d integer) at priority %d\n",numberLayered,numberInteger,1+(jLayer/100));
1727        char buffer[50];
1728        for (int i=1;i<nCheck;i++) {
1729          if (countsI[i]||countsC[i]) {
1730            if (i==1)
1731              sprintf(buffer," ==    zero            ");
1732            else if (i<nCheck-1)
1733              sprintf(buffer,"> %6d and <= %6d ",check[i-1],check[i]);
1734            else
1735              sprintf(buffer,"> %6d                ",check[i-1]);
1736            printf("%s %8d integers and %8d continuous\n",buffer,countsI[i],countsC[i]);
1737          }
1738        }
1739      }
1740    }
1741  }
1742  delete [] counts;
1743  // Now do fixing
1744  {
1745    // switch off presolve and up weight
1746    ClpSolve solveOptions;
1747    //solveOptions.setPresolveType(ClpSolve::presolveOff,0);
1748    solveOptions.setSolveType(ClpSolve::usePrimalorSprint);
1749    //solveOptions.setSpecialOption(1,3,30); // sprint
1750    int numberColumns = lpSolver->numberColumns();
1751    int iColumn;
1752    bool allSlack=true; 
1753    for (iColumn=0;iColumn<numberColumns;iColumn++) {
1754      if (lpSolver->getColumnStatus(iColumn)==ClpSimplex::basic) {
1755        allSlack=false;
1756        break;
1757      } 
1758    }
1759    if (allSlack)
1760      solveOptions.setSpecialOption(1,2,50); // idiot
1761    lpSolver->setInfeasibilityCost(1.0e11);
1762    lpSolver->defaultFactorizationFrequency();
1763    if (doAction!=11)
1764      lpSolver->initialSolve(solveOptions);
1765    double * columnLower = lpSolver->columnLower();
1766    double * columnUpper = lpSolver->columnUpper();
1767    double * fullSolution = lpSolver->primalColumnSolution();
1768    const double * dj = lpSolver->dualColumnSolution();
1769    int iPass=0;
1770#define MAXPROB 2
1771    ClpSimplex models[MAXPROB];
1772    int pass[MAXPROB];
1773    int kPass=-1;
1774    int kLayer=0;
1775    int skipZero=0;
1776    if (skipZero2==-1)
1777      skipZero2=40; //-1;
1778    /* 0 fixed to 0 by choice
1779       1 lb of 1 by choice
1780       2 fixed to 0 by another
1781       3 as 2 but this go
1782       -1 free
1783    */
1784    char * state = new char [numberColumns];
1785    for (iColumn=0;iColumn<numberColumns;iColumn++) 
1786      state[iColumn]=-1;
1787    while (true) {
1788      double largest=-0.1;
1789      double smallest=1.1;
1790      int iLargest=-1;
1791      int iSmallest=-1;
1792      int atZero=0;
1793      int atOne=0;
1794      int toZero=0;
1795      int toOne=0;
1796      int numberFree=0;
1797      int numberGreater=0;
1798      columnLower = lpSolver->columnLower();
1799      columnUpper = lpSolver->columnUpper();
1800      fullSolution = lpSolver->primalColumnSolution();
1801      if (doAction==11) {
1802        {
1803          double * columnLower = lpSolver->columnLower();
1804          double * columnUpper = lpSolver->columnUpper();
1805          //      lpSolver->dual();
1806          memcpy(columnLower,saveColumnLower,numberColumns*sizeof(double));
1807          memcpy(columnUpper,saveColumnUpper,numberColumns*sizeof(double));
1808          //      lpSolver->dual();
1809          int iColumn;
1810          for (iColumn=0;iColumn<numberColumns;iColumn++) {
1811            if (columnUpper[iColumn] > columnLower[iColumn]+1.0e-8) {
1812              if (clpSolver->isInteger(iColumn)) {
1813                double value = lastSolution[iColumn];
1814                int iValue = static_cast<int> (value+0.5);
1815                assert (fabs(value-static_cast<double> (iValue))<1.0e-3);
1816                assert (iValue>=columnLower[iColumn]&&
1817                        iValue<=columnUpper[iColumn]);
1818                columnLower[iColumn]=iValue;
1819                columnUpper[iColumn]=iValue;
1820              }
1821            }
1822          }
1823          lpSolver->initialSolve(solveOptions);
1824          memcpy(columnLower,saveColumnLower,numberColumns*sizeof(double));
1825          memcpy(columnUpper,saveColumnUpper,numberColumns*sizeof(double));
1826        }
1827        for (iColumn=0;iColumn<numberColumns;iColumn++) {
1828          if (columnUpper[iColumn] > columnLower[iColumn]+1.0e-8) {
1829            if (clpSolver->isInteger(iColumn)) {
1830              double value = lastSolution[iColumn];
1831              int iValue = static_cast<int> (value+0.5);
1832              assert (fabs(value-static_cast<double> (iValue))<1.0e-3);
1833              assert (iValue>=columnLower[iColumn]&&
1834                      iValue<=columnUpper[iColumn]);
1835              if (!fix[iColumn]) {
1836                if (iValue==0) {
1837                  state[iColumn]=0;
1838                  assert (!columnLower[iColumn]);
1839                  columnUpper[iColumn]=0.0;
1840                } else if (iValue==1) {
1841                  state[iColumn]=1;
1842                  columnLower[iColumn]=1.0;
1843                } else {
1844                  // leave fixed
1845                  columnLower[iColumn]=iValue;
1846                  columnUpper[iColumn]=iValue;
1847                }
1848              } else if (iValue==0) {
1849                state[iColumn]=10;
1850                columnUpper[iColumn]=0.0;
1851              } else {
1852                // leave fixed
1853                columnLower[iColumn]=iValue;
1854                columnUpper[iColumn]=iValue;
1855              }
1856            }
1857          }
1858        }
1859        int jLayer=0;
1860        int nFixed=-1;
1861        int nTotalFixed=0;
1862        while (nFixed) {
1863          nFixed=0;
1864          for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1865            if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) {
1866              for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
1867                int jColumn=otherColumn[i];
1868                if (columnUpper[jColumn]) {
1869                  bool canFix=true;
1870                  for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) {
1871                    int kColumn=otherColumn2[k];
1872                    if (state[kColumn]==1) {
1873                      canFix=false;
1874                      break;
1875                    }
1876                  }
1877                  if (canFix) {
1878                    columnUpper[jColumn]=0.0;
1879                    nFixed++;
1880                  }
1881                }
1882              }
1883            }
1884          }
1885          nTotalFixed += nFixed;
1886          jLayer += 100;
1887        }
1888        printf("This fixes %d variables in lower priorities\n",nTotalFixed);
1889        break;
1890      }
1891      for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1892        if (!clpSolver->isInteger(iColumn)||fix[iColumn]>kLayer)
1893          continue;
1894        // skip if fixes nothing
1895        if (fixColumn[iColumn+1]-fixColumn[iColumn]<=skipZero2)
1896          continue;
1897        double value = fullSolution[iColumn];
1898        if (value>1.00001) {
1899          numberGreater++;
1900          continue;
1901        }
1902        double lower = columnLower[iColumn];
1903        double upper = columnUpper[iColumn];
1904        if (lower==upper) {
1905          if (lower)
1906            atOne++;
1907          else
1908            atZero++;
1909          continue;
1910        }
1911        if (value<1.0e-7) {
1912          toZero++;
1913          columnUpper[iColumn]=0.0;
1914          state[iColumn]=10;
1915          continue;
1916        }
1917        if (value>1.0-1.0e-7) {
1918          toOne++;
1919          columnLower[iColumn]=1.0;
1920          state[iColumn]=1;
1921          continue;
1922        }
1923        numberFree++;
1924        // skip if fixes nothing
1925        if (fixColumn[iColumn+1]-fixColumn[iColumn]<=skipZero)
1926          continue;
1927        if (value<smallest) {
1928          smallest=value;
1929          iSmallest=iColumn;
1930        }
1931        if (value>largest) {
1932          largest=value;
1933          iLargest=iColumn;
1934        }
1935      }
1936      if (toZero||toOne)
1937        printf("%d at 0 fixed and %d at one fixed\n",toZero,toOne);
1938      printf("%d variables free, %d fixed to 0, %d to 1 - smallest %g, largest %g\n",
1939             numberFree,atZero,atOne,smallest,largest);
1940      if (numberGreater&&!iPass)
1941        printf("%d variables have value > 1.0\n",numberGreater);
1942      //skipZero2=0; // leave 0 fixing
1943      int jLayer=0;
1944      int nFixed=-1;
1945      int nTotalFixed=0;
1946      while (nFixed) {
1947        nFixed=0;
1948        for ( iColumn=0;iColumn<numberColumns;iColumn++) {
1949          if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) {
1950            for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
1951              int jColumn=otherColumn[i];
1952              if (columnUpper[jColumn]) {
1953                bool canFix=true;
1954                for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) {
1955                  int kColumn=otherColumn2[k];
1956                  if (state[kColumn]==1) {
1957                    canFix=false;
1958                    break;
1959                  }
1960                }
1961                if (canFix) {
1962                  columnUpper[jColumn]=0.0;
1963                  nFixed++;
1964                }
1965              }
1966            }
1967          }
1968        }
1969        nTotalFixed += nFixed;
1970        jLayer += 100;
1971      }
1972      printf("This fixes %d variables in lower priorities\n",nTotalFixed);
1973      if (iLargest<0||numberFree<=leaveIntFree)
1974        break;
1975      double movement;
1976      int way;
1977      if (smallest<=1.0-largest&&smallest<0.2&&largest<fixAboveValue) {
1978        columnUpper[iSmallest]=0.0;
1979        state[iSmallest]=0;
1980        movement=smallest;
1981        way=-1;
1982      } else {
1983        columnLower[iLargest]=1.0;
1984        state[iLargest]=1;
1985        movement=1.0-largest;
1986        way=1;
1987      }
1988      double saveObj = lpSolver->objectiveValue();
1989      iPass++;
1990      kPass = iPass%MAXPROB;
1991      models[kPass]=*lpSolver;
1992      if (way==-1) {
1993        // fix others
1994        for (int i=fixColumn[iSmallest];i<fixColumn[iSmallest+1];i++) {
1995          int jColumn=otherColumn[i];
1996          if (state[jColumn]==-1) {
1997            columnUpper[jColumn]=0.0;
1998            state[jColumn]=3;
1999          }
2000        }
2001      }
2002      pass[kPass]=iPass;
2003      double maxCostUp = COIN_DBL_MAX;
2004      objective = lpSolver->getObjCoefficients() ;
2005      if (way==-1)
2006        maxCostUp= (1.0-movement)*objective[iSmallest];
2007      lpSolver->setDualObjectiveLimit(saveObj+maxCostUp);
2008      crunchIt(lpSolver);
2009      double moveObj = lpSolver->objectiveValue()-saveObj;
2010      printf("movement %s was %g costing %g\n",
2011             (way==-1) ? "down" : "up",movement,moveObj);
2012      if (way==-1&&(moveObj>=(1.0-movement)*objective[iSmallest]||lpSolver->status())) {
2013        // go up
2014        columnLower = models[kPass].columnLower();
2015        columnUpper = models[kPass].columnUpper();
2016        columnLower[iSmallest]=1.0;
2017        columnUpper[iSmallest]=saveColumnUpper[iSmallest];
2018        *lpSolver=models[kPass];
2019        columnLower = lpSolver->columnLower();
2020        columnUpper = lpSolver->columnUpper();
2021        fullSolution = lpSolver->primalColumnSolution();
2022        dj = lpSolver->dualColumnSolution();
2023        columnLower[iSmallest]=1.0;
2024        columnUpper[iSmallest]=saveColumnUpper[iSmallest];
2025        state[iSmallest]=1;
2026        // unfix others
2027        for (int i=fixColumn[iSmallest];i<fixColumn[iSmallest+1];i++) {
2028          int jColumn=otherColumn[i];
2029          if (state[jColumn]==3) {
2030            columnUpper[jColumn]=saveColumnUpper[jColumn];
2031            state[jColumn]=-1;
2032          }
2033        }
2034        crunchIt(lpSolver);
2035      }
2036      models[kPass]=*lpSolver;
2037    }
2038    lpSolver->dual();
2039    printf("Fixing took %g seconds\n",CoinCpuTime()-time1);
2040    columnLower = lpSolver->columnLower();
2041    columnUpper = lpSolver->columnUpper();
2042    fullSolution = lpSolver->primalColumnSolution();
2043    dj = lpSolver->dualColumnSolution();
2044    int * sort = new int[numberColumns];
2045    double * dsort = new double[numberColumns];
2046    int chunk=20;
2047    int iRelax=0;
2048    //double fractionFixed=6.0/8.0;
2049    // relax while lots fixed
2050    while (true) {
2051      if (skipZero2>10&&doAction<10)
2052        break;
2053      iRelax++;
2054      int n=0;
2055      double sum0=0.0;
2056      double sum00=0.0;
2057      double sum1=0.0;
2058      for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2059        if (!clpSolver->isInteger(iColumn)||fix[iColumn]>kLayer)
2060          continue;
2061        // skip if fixes nothing
2062        if (fixColumn[iColumn+1]-fixColumn[iColumn]==0&&doAction<10)
2063          continue;
2064        double djValue = dj[iColumn];
2065        if (state[iColumn]==1) {
2066          assert (columnLower[iColumn]);
2067          assert (fullSolution[iColumn]>0.1);
2068          if (djValue>0.0) {
2069            //printf("YY dj of %d at %g is %g\n",iColumn,value,djValue);
2070            sum1 += djValue;
2071            sort[n]=iColumn;
2072            dsort[n++]=-djValue;
2073          } else {
2074            //printf("dj of %d at %g is %g\n",iColumn,value,djValue);
2075          }
2076        } else if (state[iColumn]==0||state[iColumn]==10) {
2077          assert (fullSolution[iColumn]<0.1);
2078          assert (!columnUpper[iColumn]);
2079          double otherValue=0.0;
2080          int nn=0;
2081          for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
2082            int jColumn=otherColumn[i];
2083            if (columnUpper[jColumn]==0.0) {
2084              if (dj[jColumn]<-1.0e-5) {
2085                nn++;
2086                otherValue += dj[jColumn]; // really need to look at rest
2087              }
2088            }
2089          }
2090          if (djValue<-1.0e-2||otherValue<-1.0e-2) {
2091            //printf("XX dj of %d at %g is %g - %d out of %d contribute %g\n",iColumn,value,djValue,
2092            // nn,fixColumn[iColumn+1]-fixColumn[iColumn],otherValue);
2093            if (djValue<1.0e-8) {
2094              sum0 -= djValue;
2095              sum00 -= otherValue;
2096              sort[n]=iColumn;
2097              if (djValue<-1.0e-2) 
2098                dsort[n++]=djValue+otherValue;
2099              else
2100                dsort[n++]=djValue+0.001*otherValue;
2101            }
2102          } else {
2103            //printf("dj of %d at %g is %g - no contribution from %d\n",iColumn,value,djValue,
2104            //   fixColumn[iColumn+1]-fixColumn[iColumn]);
2105          }
2106        }
2107      }
2108      CoinSort_2(dsort,dsort+n,sort);
2109      double * originalColumnLower = saveColumnLower;
2110      double * originalColumnUpper = saveColumnUpper;
2111      double * lo = CoinCopyOfArray(columnLower,numberColumns);
2112      double * up = CoinCopyOfArray(columnUpper,numberColumns);
2113      for (int k=0;k<CoinMin(chunk,n);k++) {
2114        iColumn = sort[k];
2115        state[iColumn]=-2;
2116      }
2117      memcpy(columnLower,originalColumnLower,numberColumns*sizeof(double));
2118      memcpy(columnUpper,originalColumnUpper,numberColumns*sizeof(double));
2119      int nFixed=0;
2120      int nFixed0=0;
2121      int nFixed1=0;
2122      for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2123        if (state[iColumn]==0||state[iColumn]==10) {
2124          columnUpper[iColumn]=0.0;
2125          assert (lo[iColumn]==0.0);
2126          nFixed++;
2127          nFixed0++;
2128          for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
2129            int jColumn=otherColumn[i];
2130            if (columnUpper[jColumn]) {
2131              bool canFix=true;
2132              for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) {
2133                int kColumn=otherColumn2[k];
2134                if (state[kColumn]==1||state[kColumn]==-2) {
2135                  canFix=false;
2136                  break;
2137                }
2138              }
2139              if (canFix) {
2140                columnUpper[jColumn]=0.0;
2141                assert (lo[jColumn]==0.0);
2142                nFixed++;
2143              }
2144            }
2145          }
2146        } else if (state[iColumn]==1) {
2147          columnLower[iColumn]=1.0;
2148          nFixed1++;
2149        }
2150      }
2151      printf("%d fixed %d orig 0 %d 1\n",nFixed,nFixed0,nFixed1);
2152      int jLayer=0;
2153      nFixed=-1;
2154      int nTotalFixed=0;
2155      while (nFixed) {
2156        nFixed=0;
2157        for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2158          if (columnUpper[iColumn]==0.0&&fix[iColumn]==jLayer) {
2159            for (int i=fixColumn[iColumn];i<fixColumn[iColumn+1];i++) {
2160              int jColumn=otherColumn[i];
2161              if (columnUpper[jColumn]) {
2162                bool canFix=true;
2163                for (int k=fixColumn2[jColumn];k<fixColumn2[jColumn+1];k++) {
2164                  int kColumn=otherColumn2[k];
2165                  if (state[kColumn]==1||state[kColumn]==-2) {
2166                    canFix=false;
2167                    break;
2168                  }
2169                }
2170                if (canFix) {
2171                  columnUpper[jColumn]=0.0;
2172                  assert (lo[jColumn]==0.0);
2173                  nFixed++;
2174                }
2175              }
2176            }
2177          }
2178        }
2179        nTotalFixed += nFixed;
2180        jLayer += 100;
2181      }
2182      nFixed=0;
2183      int nFixedI=0;
2184      for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2185        if (columnLower[iColumn]==columnUpper[iColumn]) {
2186          if (clpSolver->isInteger(iColumn))
2187            nFixedI++;
2188          nFixed++;
2189        }
2190      }
2191      printf("This fixes %d variables in lower priorities - total %d (%d integer) - all target %d, int target %d\n",
2192             nTotalFixed,nFixed,nFixedI,static_cast<int>(fractionFixed*numberColumns),static_cast<int> (fractionIntFixed*numberInteger));
2193      int nBad=0;
2194      int nRelax=0;
2195      for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2196        if (lo[iColumn]<columnLower[iColumn]||
2197            up[iColumn]>columnUpper[iColumn]) {
2198          printf("bad %d old %g %g, new %g %g\n",iColumn,lo[iColumn],up[iColumn],
2199                 columnLower[iColumn],columnUpper[iColumn]);
2200          nBad++;
2201        }
2202        if (lo[iColumn]>columnLower[iColumn]||
2203            up[iColumn]<columnUpper[iColumn]) {
2204          nRelax++;
2205        }
2206      }
2207      printf("%d relaxed\n",nRelax);
2208      if (iRelax>20&&nRelax==chunk)
2209        nRelax=0;
2210      if (iRelax>50)
2211        nRelax=0;
2212      assert (!nBad);
2213      delete [] lo;
2214      delete [] up;
2215      lpSolver->primal(1);
2216      if (nFixed<fractionFixed*numberColumns||nFixedI<fractionIntFixed*numberInteger||!nRelax)
2217        break;
2218    }
2219    delete [] state;
2220    delete [] sort;
2221    delete [] dsort;
2222  }
2223  delete [] fix;
2224  delete [] fixColumn;
2225  delete [] otherColumn;
2226  delete [] otherColumn2;
2227  delete [] fixColumn2;
2228  // See if was presolved
2229  if (originalColumns) {
2230    columnLower = lpSolver->columnLower();
2231    columnUpper = lpSolver->columnUpper();
2232    for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2233      saveColumnLower[iColumn] = columnLower[iColumn];
2234      saveColumnUpper[iColumn] = columnUpper[iColumn];
2235    }
2236    pinfo.postsolve(true);
2237    columnLower = originalLpSolver->columnLower();
2238    columnUpper = originalLpSolver->columnUpper();
2239    double * newColumnLower = lpSolver->columnLower();
2240    double * newColumnUpper = lpSolver->columnUpper();
2241    for (iColumn=0;iColumn<numberColumns;iColumn++) {
2242      int jColumn = originalColumns[iColumn];
2243      columnLower[jColumn] = CoinMax(columnLower[jColumn],newColumnLower[iColumn]);
2244      columnUpper[jColumn] = CoinMin(columnUpper[jColumn],newColumnUpper[iColumn]);
2245    }
2246    numberColumns = originalLpSolver->numberColumns();
2247    delete [] originalColumns;
2248  }
2249  delete [] saveColumnLower;
2250  delete [] saveColumnUpper;
2251  if (!originalColumns) {
2252    // Basis
2253    memcpy(originalLpSolver->statusArray(),lpSolver->statusArray(),numberRows+numberColumns);
2254    memcpy(originalLpSolver->primalColumnSolution(),lpSolver->primalColumnSolution(),numberColumns*sizeof(double));
2255    memcpy(originalLpSolver->primalRowSolution(),lpSolver->primalRowSolution(),numberRows*sizeof(double));
2256    // Fix in solver
2257    columnLower = lpSolver->columnLower();
2258    columnUpper = lpSolver->columnUpper();
2259  }
2260  double * originalColumnLower = originalLpSolver->columnLower();
2261  double * originalColumnUpper = originalLpSolver->columnUpper();
2262  // number fixed
2263  doAction=0;
2264  for ( iColumn=0;iColumn<numberColumns;iColumn++) {
2265    originalColumnLower[iColumn] = columnLower[iColumn];
2266    originalColumnUpper[iColumn] = columnUpper[iColumn];
2267    if (columnLower[iColumn]==columnUpper[iColumn])
2268      doAction++;
2269  }
2270  printf("%d fixed by vub preprocessing\n",doAction);
2271  if (originalColumns) {
2272    originalLpSolver->initialSolve();
2273  }
2274  delete clpSolver;
2275  return NULL;
2276}
2277#endif
2278#ifdef COIN_HAS_LINK
2279/*  Returns OsiSolverInterface (User should delete)
2280    On entry numberKnapsack is maximum number of Total entries
2281*/
2282static OsiSolverInterface * 
2283expandKnapsack(CoinModel & model, int * whichColumn, int * knapsackStart, 
2284               int * knapsackRow, int &numberKnapsack,
2285               CglStored & stored, int logLevel,
2286               int fixedPriority, int SOSPriority,CoinModel & tightenedModel)
2287{
2288  int maxTotal = numberKnapsack;
2289  // load from coin model
2290  OsiSolverLink *si = new OsiSolverLink();
2291  OsiSolverInterface * finalModel=NULL;
2292  si->setDefaultMeshSize(0.001);
2293  // need some relative granularity
2294  si->setDefaultBound(100.0);
2295  si->setDefaultMeshSize(0.01);
2296  si->setDefaultBound(100000.0);
2297  si->setIntegerPriority(1000);
2298  si->setBiLinearPriority(10000);
2299  si->load(model,true,logLevel);
2300  // get priorities
2301  const int * priorities=model.priorities();
2302  int numberColumns = model.numberColumns();
2303  if (priorities) {
2304    OsiObject ** objects = si->objects();
2305    int numberObjects = si->numberObjects();
2306    for (int iObj = 0;iObj<numberObjects;iObj++) {
2307      int iColumn = objects[iObj]->columnNumber();
2308      if (iColumn>=0&&iColumn<numberColumns) {
2309#ifndef NDEBUG
2310        OsiSimpleInteger * obj =
2311          dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
2312#endif
2313        assert (obj);
2314        int iPriority = priorities[iColumn];
2315        if (iPriority>0)
2316          objects[iObj]->setPriority(iPriority);
2317      }
2318    }
2319    if (fixedPriority>0) {
2320      si->setFixedPriority(fixedPriority);
2321    }
2322    if (SOSPriority<0)
2323      SOSPriority=100000;
2324  } 
2325  CoinModel coinModel=*si->coinModel();
2326  assert(coinModel.numberRows()>0);
2327  tightenedModel = coinModel;
2328  int numberRows = coinModel.numberRows();
2329  // Mark variables
2330  int * whichKnapsack = new int [numberColumns];
2331  int iRow,iColumn;
2332  for (iColumn=0;iColumn<numberColumns;iColumn++) 
2333    whichKnapsack[iColumn]=-1;
2334  int kRow;
2335  bool badModel=false;
2336  // analyze
2337  if (logLevel>1) {
2338    for (iRow=0;iRow<numberRows;iRow++) {
2339      /* Just obvious one at first
2340         positive non unit coefficients
2341         all integer
2342         positive rowUpper
2343         for now - linear (but further down in code may use nonlinear)
2344         column bounds should be tight
2345      */
2346      //double lower = coinModel.getRowLower(iRow);
2347      double upper = coinModel.getRowUpper(iRow);
2348      if (upper<1.0e10) {
2349        CoinModelLink triple=coinModel.firstInRow(iRow);
2350        bool possible=true;
2351        int n=0;
2352        int n1=0;
2353        while (triple.column()>=0) {
2354          int iColumn = triple.column();
2355          const char *  el = coinModel.getElementAsString(iRow,iColumn);
2356          if (!strcmp("Numeric",el)) {
2357            if (coinModel.columnLower(iColumn)==coinModel.columnUpper(iColumn)) {
2358              triple=coinModel.next(triple);
2359              continue; // fixed
2360            }
2361            double value=coinModel.getElement(iRow,iColumn);
2362            if (value<0.0) {
2363              possible=false;
2364            } else {
2365              n++;
2366              if (value==1.0)
2367                n1++;
2368              if (coinModel.columnLower(iColumn)<0.0)
2369                possible=false;
2370              if (!coinModel.isInteger(iColumn))
2371                possible=false;
2372              if (whichKnapsack[iColumn]>=0)
2373                possible=false;
2374            }
2375          } else {
2376            possible=false; // non linear
2377          } 
2378          triple=coinModel.next(triple);
2379        }
2380        if (n-n1>1&&possible) {
2381          double lower = coinModel.getRowLower(iRow);
2382          double upper = coinModel.getRowUpper(iRow);
2383          CoinModelLink triple=coinModel.firstInRow(iRow);
2384          while (triple.column()>=0) {
2385            int iColumn = triple.column();
2386            lower -= coinModel.columnLower(iColumn)*triple.value();
2387            upper -= coinModel.columnLower(iColumn)*triple.value();
2388            triple=coinModel.next(triple);
2389          }
2390          printf("%d is possible %g <=",iRow,lower);
2391          // print
2392          triple=coinModel.firstInRow(iRow);
2393          while (triple.column()>=0) {
2394            int iColumn = triple.column();
2395            if (coinModel.columnLower(iColumn)!=coinModel.columnUpper(iColumn))
2396              printf(" (%d,el %g up %g)",iColumn,triple.value(),
2397                     coinModel.columnUpper(iColumn)-coinModel.columnLower(iColumn));
2398            triple=coinModel.next(triple);
2399          }
2400          printf(" <= %g\n",upper);
2401        }
2402      }
2403    }
2404  }
2405  numberKnapsack=0;
2406  for (kRow=0;kRow<numberRows;kRow++) {
2407    iRow=kRow;
2408    /* Just obvious one at first
2409       positive non unit coefficients
2410       all integer
2411       positive rowUpper
2412       for now - linear (but further down in code may use nonlinear)
2413       column bounds should be tight
2414    */
2415    //double lower = coinModel.getRowLower(iRow);
2416    double upper = coinModel.getRowUpper(iRow);
2417    if (upper<1.0e10) {
2418      CoinModelLink triple=coinModel.firstInRow(iRow);
2419      bool possible=true;
2420      int n=0;
2421      int n1=0;
2422      while (triple.column()>=0) {
2423        int iColumn = triple.column();
2424        const char *  el = coinModel.getElementAsString(iRow,iColumn);
2425        if (!strcmp("Numeric",el)) {
2426          if (coinModel.columnLower(iColumn)==coinModel.columnUpper(iColumn)) {
2427            triple=coinModel.next(triple);
2428            continue; // fixed
2429          }
2430          double value=coinModel.getElement(iRow,iColumn);
2431          if (value<0.0) {
2432            possible=false;
2433          } else {
2434            n++;
2435            if (value==1.0)
2436              n1++;
2437            if (coinModel.columnLower(iColumn)<0.0)
2438              possible=false;
2439            if (!coinModel.isInteger(iColumn))
2440              possible=false;
2441            if (whichKnapsack[iColumn]>=0)
2442              possible=false;
2443          }
2444        } else {
2445          possible=false; // non linear
2446        } 
2447        triple=coinModel.next(triple);
2448      }
2449      if (n-n1>1&&possible) {
2450        // try
2451        CoinModelLink triple=coinModel.firstInRow(iRow);
2452        while (triple.column()>=0) {
2453          int iColumn = triple.column();
2454          if (coinModel.columnLower(iColumn)!=coinModel.columnUpper(iColumn))
2455            whichKnapsack[iColumn]=numberKnapsack;
2456          triple=coinModel.next(triple);
2457        }
2458        knapsackRow[numberKnapsack++]=iRow;
2459      }
2460    }
2461  }
2462  if (logLevel>0)
2463    printf("%d out of %d candidate rows are possible\n",numberKnapsack,numberRows);
2464  // Check whether we can get rid of nonlinearities
2465  /* mark rows
2466     -2 in knapsack and other variables
2467     -1 not involved
2468     n only in knapsack n
2469  */
2470  int * markRow = new int [numberRows];
2471  for (iRow=0;iRow<numberRows;iRow++) 
2472    markRow[iRow]=-1;
2473  int canDo=1; // OK and linear
2474  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2475    CoinModelLink triple=coinModel.firstInColumn(iColumn);
2476    int iKnapsack = whichKnapsack[iColumn];
2477    bool linear=true;
2478    // See if quadratic objective
2479    const char * expr = coinModel.getColumnObjectiveAsString(iColumn);
2480    if (strcmp(expr,"Numeric")) {
2481      linear=false;
2482    }
2483    while (triple.row()>=0) {
2484      int iRow = triple.row();
2485      if (iKnapsack>=0) {
2486        if (markRow[iRow]==-1) {
2487          markRow[iRow]=iKnapsack;
2488        } else if (markRow[iRow]!=iKnapsack) {
2489          markRow[iRow]=-2;
2490        }
2491      }
2492      const char * expr = coinModel.getElementAsString(iRow,iColumn);
2493      if (strcmp(expr,"Numeric")) {
2494        linear=false;
2495      }
2496      triple=coinModel.next(triple);
2497    }
2498    if (!linear) {
2499      if (whichKnapsack[iColumn]<0) {
2500        canDo=0;
2501        break;
2502      } else {
2503        canDo=2;
2504      }
2505    }
2506  }
2507  int * markKnapsack = NULL;
2508  double * coefficient = NULL;
2509  double * linear = NULL;
2510  int * whichRow = NULL;
2511  int * lookupRow = NULL;
2512  badModel=(canDo==0);
2513  if (numberKnapsack&&canDo) {
2514    /* double check - OK if
2515       no nonlinear
2516       nonlinear only on columns in knapsack
2517       nonlinear only on columns in knapsack * ONE other - same for all in knapsack
2518       AND that is only row connected to knapsack
2519       (theoretically could split knapsack if two other and small numbers)
2520       also ONE could be ONE expression - not just a variable
2521    */
2522    int iKnapsack;
2523    markKnapsack = new int [numberKnapsack];
2524    coefficient = new double [numberKnapsack];
2525    linear = new double [numberColumns];
2526    for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) 
2527      markKnapsack[iKnapsack]=-1;
2528    if (canDo==2) {
2529      for (iRow=-1;iRow<numberRows;iRow++) {
2530        int numberOdd;
2531        CoinPackedMatrix * row = coinModel.quadraticRow(iRow,linear,numberOdd);
2532        if (row) {
2533          // see if valid
2534          const double * element = row->getElements();
2535          const int * column = row->getIndices();
2536          const CoinBigIndex * columnStart = row->getVectorStarts();
2537          const int * columnLength = row->getVectorLengths();
2538          int numberLook = row->getNumCols();
2539          for (int i=0;i<numberLook;i++) {
2540            int iKnapsack=whichKnapsack[i];
2541            if (iKnapsack<0) {
2542              // might be able to swap - but for now can't have knapsack in
2543              for (int j=columnStart[i];j<columnStart[i]+columnLength[i];j++) {
2544                int iColumn = column[j];
2545                if (whichKnapsack[iColumn]>=0) {
2546                  canDo=0; // no good
2547                  badModel=true;
2548                  break;
2549                }
2550              }
2551            } else {
2552              // OK if in same knapsack - or maybe just one
2553              int marked=markKnapsack[iKnapsack];
2554              for (int j=columnStart[i];j<columnStart[i]+columnLength[i];j++) {
2555                int iColumn = column[j];
2556                if (whichKnapsack[iColumn]!=iKnapsack&&whichKnapsack[iColumn]>=0) {
2557                  canDo=0; // no good
2558                  badModel=true;
2559                  break;
2560                } else if (marked==-1) {
2561                  markKnapsack[iKnapsack]=iColumn;
2562                  marked=iColumn;
2563                  coefficient[iKnapsack]=element[j];
2564                  coinModel.associateElement(coinModel.columnName(iColumn),1.0);
2565                } else if (marked!=iColumn) {
2566                  badModel=true;
2567                  canDo=0; // no good
2568                  break;
2569                } else {
2570                  // could manage with different coefficients - but for now ...
2571                  assert(coefficient[iKnapsack]==element[j]);
2572                }
2573              }
2574            }
2575          }
2576          delete row;
2577        }
2578      }
2579    }
2580    if (canDo) {
2581      // for any rows which are cuts
2582      whichRow = new int [numberRows];
2583      lookupRow = new int [numberRows];
2584      bool someNonlinear=false;
2585      double maxCoefficient=1.0;
2586      for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
2587        if (markKnapsack[iKnapsack]>=0) {
2588          someNonlinear=true;
2589          int iColumn = markKnapsack[iKnapsack];
2590          maxCoefficient = CoinMax(maxCoefficient,fabs(coefficient[iKnapsack]*coinModel.columnUpper(iColumn)));
2591        }
2592      }
2593      if (someNonlinear) {
2594        // associate all columns to stop possible error messages
2595        for (iColumn=0;iColumn<numberColumns;iColumn++) {
2596          coinModel.associateElement(coinModel.columnName(iColumn),1.0);
2597        }
2598      }
2599      ClpSimplex tempModel;
2600      tempModel.loadProblem(coinModel);
2601      // Create final model - first without knapsacks
2602      int nCol=0;
2603      int nRow=0;
2604      for (iRow=0;iRow<numberRows;iRow++) {
2605        if (markRow[iRow]<0) {
2606          lookupRow[iRow]=nRow;
2607          whichRow[nRow++]=iRow;
2608        } else {
2609          lookupRow[iRow]=-1;
2610        }
2611      }
2612      for (iColumn=0;iColumn<numberColumns;iColumn++) {
2613        if (whichKnapsack[iColumn]<0)
2614          whichColumn[nCol++]=iColumn;
2615      }
2616      ClpSimplex finalModelX(&tempModel,nRow,whichRow,nCol,whichColumn,false,false,false);
2617      OsiClpSolverInterface finalModelY(&finalModelX,true);
2618      finalModel = finalModelY.clone();
2619      finalModelY.releaseClp();
2620      // Put back priorities
2621      const int * priorities=model.priorities();
2622      if (priorities) {
2623        finalModel->findIntegers(false);
2624        OsiObject ** objects = finalModel->objects();
2625        int numberObjects = finalModel->numberObjects();
2626        for (int iObj = 0;iObj<numberObjects;iObj++) {
2627          int iColumn = objects[iObj]->columnNumber();
2628          if (iColumn>=0&&iColumn<nCol) {
2629#ifndef NDEBUG
2630            OsiSimpleInteger * obj =
2631              dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
2632#endif
2633            assert (obj);
2634            int iPriority = priorities[whichColumn[iColumn]];
2635            if (iPriority>0)
2636              objects[iObj]->setPriority(iPriority);
2637          }
2638        }
2639      }
2640      for (iRow=0;iRow<numberRows;iRow++) {
2641        whichRow[iRow]=iRow;
2642      }
2643      int numberOther=finalModel->getNumCols();
2644      int nLargest=0;
2645      int nelLargest=0;
2646      int nTotal=0;
2647      for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
2648        iRow = knapsackRow[iKnapsack];
2649        int nCreate = maxTotal;
2650        int nelCreate=coinModel.expandKnapsack(iRow,nCreate,NULL,NULL,NULL,NULL);
2651        if (nelCreate<0)
2652          badModel=true;
2653        nTotal+=nCreate;
2654        nLargest = CoinMax(nLargest,nCreate);
2655        nelLargest = CoinMax(nelLargest,nelCreate);
2656      }
2657      if (nTotal>maxTotal) 
2658        badModel=true;
2659      if (!badModel) {
2660        // Now arrays for building
2661        nelLargest = CoinMax(nelLargest,nLargest)+1;
2662        double * buildObj = new double [nLargest];
2663        double * buildElement = new double [nelLargest];
2664        int * buildStart = new int[nLargest+1];
2665        int * buildRow = new int[nelLargest];
2666        // alow for integers in knapsacks
2667        OsiObject ** object = new OsiObject * [numberKnapsack+nTotal];
2668        int nSOS=0;
2669        int nObj=numberKnapsack;
2670        for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
2671          knapsackStart[iKnapsack]=finalModel->getNumCols();
2672          iRow = knapsackRow[iKnapsack];
2673          int nCreate = 10000;
2674          coinModel.expandKnapsack(iRow,nCreate,buildObj,buildStart,buildRow,buildElement);
2675          // Redo row numbers
2676          for (iColumn=0;iColumn<nCreate;iColumn++) {
2677            for (int j=buildStart[iColumn];j<buildStart[iColumn+1];j++) {
2678              int jRow=buildRow[j];
2679              jRow=lookupRow[jRow];
2680              assert (jRow>=0&&jRow<nRow);
2681              buildRow[j]=jRow;
2682            }
2683          }
2684          finalModel->addCols(nCreate,buildStart,buildRow,buildElement,NULL,NULL,buildObj);
2685          int numberFinal=finalModel->getNumCols();
2686          for (iColumn=numberOther;iColumn<numberFinal;iColumn++) {
2687            if (markKnapsack[iKnapsack]<0) {
2688              finalModel->setColUpper(iColumn,maxCoefficient);
2689              finalModel->setInteger(iColumn);
2690            } else {
2691              finalModel->setColUpper(iColumn,maxCoefficient+1.0);
2692              finalModel->setInteger(iColumn);
2693            }
2694            OsiSimpleInteger * sosObject = new OsiSimpleInteger(finalModel,iColumn);
2695            sosObject->setPriority(1000000);
2696            object[nObj++]=sosObject;
2697            buildRow[iColumn-numberOther]=iColumn;
2698            buildElement[iColumn-numberOther]=1.0;
2699          }
2700          if (markKnapsack[iKnapsack]<0) {
2701            // convexity row
2702            finalModel->addRow(numberFinal-numberOther,buildRow,buildElement,1.0,1.0);
2703          } else {
2704            int iColumn = markKnapsack[iKnapsack];
2705            int n=numberFinal-numberOther;
2706            buildRow[n]=iColumn;
2707            buildElement[n++]=-fabs(coefficient[iKnapsack]);
2708            // convexity row (sort of)
2709            finalModel->addRow(n,buildRow,buildElement,0.0,0.0);
2710            OsiSOS * sosObject = new OsiSOS(finalModel,n-1,buildRow,NULL,1);
2711            sosObject->setPriority(iKnapsack+SOSPriority);
2712            // Say not integral even if is (switch off heuristics)
2713            sosObject->setIntegerValued(false);
2714            object[nSOS++]=sosObject;
2715          }
2716          numberOther=numberFinal;
2717        }
2718        finalModel->addObjects(nObj,object);
2719        for (iKnapsack=0;iKnapsack<nObj;iKnapsack++) 
2720          delete object[iKnapsack];
2721        delete [] object;
2722        // Can we move any rows to cuts
2723        const int * cutMarker = coinModel.cutMarker();
2724        if (cutMarker&&0) {
2725          printf("AMPL CUTS OFF until global cuts fixed\n");
2726          cutMarker=NULL;
2727        }
2728        if (cutMarker) {
2729          // Row copy
2730          const CoinPackedMatrix * matrixByRow = finalModel->getMatrixByRow();
2731          const double * elementByRow = matrixByRow->getElements();
2732          const int * column = matrixByRow->getIndices();
2733          const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
2734          const int * rowLength = matrixByRow->getVectorLengths();
2735         
2736          const double * rowLower = finalModel->getRowLower();
2737          const double * rowUpper = finalModel->getRowUpper();
2738          int nDelete=0;
2739          for (iRow=0;iRow<numberRows;iRow++) {
2740            if (cutMarker[iRow]&&lookupRow[iRow]>=0) {
2741              int jRow=lookupRow[iRow];
2742              whichRow[nDelete++]=jRow;
2743              int start = rowStart[jRow];
2744              stored.addCut(rowLower[jRow],rowUpper[jRow],
2745                            rowLength[jRow],column+start,elementByRow+start);
2746            }
2747          }
2748          finalModel->deleteRows(nDelete,whichRow);
2749        }
2750        knapsackStart[numberKnapsack]=finalModel->getNumCols();
2751        delete [] buildObj;
2752        delete [] buildElement;
2753        delete [] buildStart;
2754        delete [] buildRow;
2755        finalModel->writeMps("full");
2756      }
2757    }
2758  }
2759  delete [] whichKnapsack;
2760  delete [] markRow;
2761  delete [] markKnapsack;
2762  delete [] coefficient;
2763  delete [] linear;
2764  delete [] whichRow;
2765  delete [] lookupRow;
2766  delete si;
2767  si=NULL;
2768  if (!badModel&&finalModel) {
2769    finalModel->setDblParam(OsiObjOffset,coinModel.objectiveOffset());
2770    return finalModel;
2771  } else {
2772    delete finalModel;
2773    printf("can't make knapsacks - did you set fixedPriority (extra1)\n");
2774    return NULL;
2775  }
2776}
2777#if NEW_STYLE_SOLVER
2778// Fills in original solution (coinModel length)
2779static void 
2780afterKnapsack(const CoinModel & coinModel2, const int * whichColumn, const int * knapsackStart, 
2781                   const int * knapsackRow, int numberKnapsack,
2782                   const double * knapsackSolution, double * solution, int logLevel)
2783{
2784  CoinModel coinModel = coinModel2;
2785  int numberColumns = coinModel.numberColumns();
2786  int iColumn;
2787  // associate all columns to stop possible error messages
2788  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2789    coinModel.associateElement(coinModel.columnName(iColumn),1.0);
2790  }
2791  CoinZeroN(solution,numberColumns);
2792  int nCol=knapsackStart[0];
2793  for (iColumn=0;iColumn<nCol;iColumn++) {
2794    int jColumn = whichColumn[iColumn];
2795    solution[jColumn]=knapsackSolution[iColumn];
2796  }
2797  int * buildRow = new int [numberColumns]; // wild overkill
2798  double * buildElement = new double [numberColumns];
2799  int iKnapsack;
2800  for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
2801    int k=-1;
2802    double value=0.0;
2803    for (iColumn=knapsackStart[iKnapsack];iColumn<knapsackStart[iKnapsack+1];iColumn++) {
2804      if (knapsackSolution[iColumn]>1.0e-5) {
2805        if (k>=0) {
2806          printf("Two nonzero values for knapsack %d at (%d,%g) and (%d,%g)\n",iKnapsack,
2807                 k,knapsackSolution[k],iColumn,knapsackSolution[iColumn]);
2808          abort();
2809        }
2810        k=iColumn;
2811        value=floor(knapsackSolution[iColumn]+0.5);
2812        assert (fabs(value-knapsackSolution[iColumn])<1.0e-5);
2813      }
2814    }
2815    if (k>=0) {
2816      int iRow = knapsackRow[iKnapsack];
2817      int nCreate = 10000;
2818      int nel=coinModel.expandKnapsack(iRow,nCreate,NULL,NULL,buildRow,buildElement,k-knapsackStart[iKnapsack]);
2819      assert (nel);
2820      if (logLevel>0) 
2821        printf("expanded column %d in knapsack %d has %d nonzero entries:\n",
2822               k-knapsackStart[iKnapsack],iKnapsack,nel);
2823      for (int i=0;i<nel;i++) {
2824        int jColumn = buildRow[i];
2825        double value = buildElement[i];
2826        if (logLevel>0) 
2827          printf("%d - original %d has value %g\n",i,jColumn,value);
2828        solution[jColumn]=value;
2829      }
2830    }
2831  }
2832  delete [] buildRow;
2833  delete [] buildElement;
2834#if 0
2835  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2836    if (solution[iColumn]>1.0e-5&&coinModel.isInteger(iColumn))
2837      printf("%d %g\n",iColumn,solution[iColumn]);
2838  }
2839#endif
2840}
2841#endif
2842#endif
2843#if 0
2844static int outDupRow(OsiSolverInterface * solver)
2845{
2846  CglDuplicateRow dupCuts(solver);
2847  CglTreeInfo info;
2848  info.level = 0;
2849  info.pass = 0;
2850  int numberRows = solver->getNumRows();
2851  info.formulation_rows = numberRows;
2852  info.inTree = false;
2853  info.strengthenRow= NULL;
2854  info.pass = 0;
2855  OsiCuts cs;
2856  dupCuts.generateCuts(*solver,cs,info);
2857  const int * duplicate = dupCuts.duplicate();
2858  // Get rid of duplicate rows
2859  int * which = new int[numberRows];
2860  int numberDrop=0;
2861  for (int iRow=0;iRow<numberRows;iRow++) {
2862    if (duplicate[iRow]==-2||duplicate[iRow]>=0)
2863      which[numberDrop++]=iRow;
2864  }
2865  if (numberDrop) {
2866    solver->deleteRows(numberDrop,which);
2867  }
2868  delete [] which;
2869  // see if we have any column cuts
2870  int numberColumnCuts = cs.sizeColCuts() ;
2871  const double * columnLower = solver->getColLower();
2872  const double * columnUpper = solver->getColUpper();
2873  for (int k = 0;k<numberColumnCuts;k++) {
2874    OsiColCut * thisCut = cs.colCutPtr(k) ;
2875    const CoinPackedVector & lbs = thisCut->lbs() ;
2876    const CoinPackedVector & ubs = thisCut->ubs() ;
2877    int j ;
2878    int n ;
2879    const int * which ;
2880    const double * values ;
2881    n = lbs.getNumElements() ;
2882    which = lbs.getIndices() ;
2883    values = lbs.getElements() ;
2884    for (j = 0;j<n;j++) {
2885      int iColumn = which[j] ;
2886      if (values[j]>columnLower[iColumn])
2887        solver->setColLower(iColumn,values[j]) ;
2888    }
2889    n = ubs.getNumElements() ;
2890    which = ubs.getIndices() ;
2891    values = ubs.getElements() ;
2892    for (j = 0;j<n;j++) {
2893      int iColumn = which[j] ;
2894      if (values[j]<columnUpper[iColumn])
2895        solver->setColUpper(iColumn,values[j]) ;
2896    }
2897  }
2898  return numberDrop;
2899}
2900#endif
2901#ifdef COIN_DEVELOP
2902void checkSOS(CbcModel * babModel, const OsiSolverInterface * solver)
2903#else
2904void checkSOS(CbcModel * /*babModel*/, const OsiSolverInterface * /*solver*/)
2905#endif
2906{
2907#ifdef COIN_DEVELOP
2908  if (!babModel->ownObjects())
2909    return;
2910#if COIN_DEVELOP>2
2911  //const double *objective = solver->getObjCoefficients() ;
2912  const double *columnLower = solver->getColLower() ;
2913  const double * columnUpper = solver->getColUpper() ;
2914  const double * solution = solver->getColSolution();
2915  //int numberColumns = solver->getNumCols() ;
2916  //int numberRows = solver->getNumRows();
2917  //double direction = solver->getObjSense();
2918  //int iRow,iColumn;
2919#endif
2920
2921  // Row copy
2922  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
2923  //const double * elementByRow = matrixByRow.getElements();
2924  //const int * column = matrixByRow.getIndices();
2925  //const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
2926  const int * rowLength = matrixByRow.getVectorLengths();
2927
2928  // Column copy
2929  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
2930  const double * element = matrixByCol.getElements();
2931  const int * row = matrixByCol.getIndices();
2932  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
2933  const int * columnLength = matrixByCol.getVectorLengths();
2934
2935  const double * rowLower = solver->getRowLower();
2936  const double * rowUpper = solver->getRowUpper();
2937  OsiObject ** objects = babModel->objects();
2938  int numberObjects = babModel->numberObjects();
2939  for (int iObj = 0;iObj<numberObjects;iObj++) {
2940    CbcSOS * objSOS =
2941      dynamic_cast <CbcSOS *>(objects[iObj]) ;
2942    if (objSOS) {
2943      int n=objSOS->numberMembers();
2944      const int * which = objSOS->members();
2945#if COIN_DEVELOP>2
2946      const double * weight = objSOS->weights();
2947#endif
2948      int type = objSOS->sosType();
2949      // convexity row?
2950      int iColumn;
2951      iColumn=which[0];
2952      int j;
2953      int convex=-1;
2954      for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
2955        int iRow = row[j];
2956        double value = element[j];
2957        if (rowLower[iRow]==1.0&&rowUpper[iRow]==1.0&&
2958            value==1.0) {
2959          // possible
2960          if (rowLength[iRow]==n) {
2961            if (convex==-1)
2962              convex=iRow;
2963            else
2964              convex=-2;
2965          }
2966        }
2967      }
2968      printf ("set %d of type %d has %d members - possible convexity row %d\n",
2969              iObj,type,n,convex);
2970      for (int i=0;i<n;i++) {
2971        iColumn = which[i];
2972        int convex2=-1;
2973        for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
2974          int iRow = row[j];
2975          if (iRow==convex) {
2976            double value = element[j];
2977            if (value==1.0) {
2978              convex2=iRow;
2979            }
2980          }
2981        }
2982        if (convex2<0&&convex>=0) {
2983          printf("odd convexity row\n");
2984          convex=-2;
2985        }
2986#if COIN_DEVELOP>2
2987        printf("col %d has weight %g and value %g, bounds %g %g\n",
2988               iColumn,weight[i],solution[iColumn],columnLower[iColumn],
2989               columnUpper[iColumn]);
2990#endif
2991      }
2992    }
2993  }
2994#endif
2995}
2996#if NEW_STYLE_SOLVER==0
2997int callCbc1(const char * input2, CbcModel & model, int callBack(CbcModel * currentSolver, int whereFrom))
2998{
2999  char * input = strdup(input2);
3000  int length = strlen(input);
3001  bool blank = input[0]=='0';
3002  int n=blank ? 0 : 1;
3003  for (int i=0;i<length;i++) {
3004    if (blank) {
3005      // look for next non blank
3006      if (input[i]==' ') {
3007        continue;
3008      } else {
3009        n++;
3010        blank=false;
3011      }
3012    } else {
3013      // look for next blank
3014      if (input[i]!=' ') {
3015        continue;
3016      } else {
3017        blank=true;
3018      }
3019    }
3020  }
3021  char ** argv = new char * [n+2];
3022  argv[0]=strdup("cbc");
3023  int i=0;
3024  while(input[i]==' ')
3025    i++;
3026  for (int j=0;j<n;j++) {
3027    int saveI=i;
3028    for (;i<length;i++) {
3029      // look for next blank
3030      if (input[i]!=' ') {
3031        continue;
3032      } else {
3033        break;
3034      }
3035    }
3036    input[i]='\0';
3037    argv[j+1]=strdup(input+saveI);
3038    while(input[i]==' ')
3039      i++;
3040  }
3041  argv[n+1]=strdup("-quit");
3042  free(input);
3043  totalTime=0.0;
3044  currentBranchModel = NULL;
3045  CbcOrClpRead_mode=1;
3046  CbcOrClpReadCommand=stdin;
3047  noPrinting=false;
3048  int returnCode = CbcMain1(n+2,const_cast<const char **>(argv),model,callBack);
3049  for (int k=0;k<n+2;k++)
3050    free(argv[k]);
3051  delete [] argv;
3052  return returnCode;
3053}
3054int callCbc1(const std::string input2, CbcModel & babSolver)
3055{
3056  char * input3 = strdup(input2.c_str());
3057  int returnCode=callCbc1(input3,babSolver);
3058  free(input3);
3059  return returnCode;
3060}
3061int callCbc(const char * input2, CbcModel & babSolver)
3062{
3063  CbcMain0(babSolver);
3064  return callCbc1(input2,babSolver);
3065}
3066int callCbc(const std::string input2, CbcModel & babSolver)
3067{
3068  char * input3 = strdup(input2.c_str());
3069  CbcMain0(babSolver);
3070  int returnCode=callCbc1(input3,babSolver);
3071  free(input3);
3072  return returnCode;
3073}
3074int callCbc(const char * input2, OsiClpSolverInterface& solver1) 
3075{
3076  CbcModel model(solver1);
3077  return callCbc(input2,model);
3078}
3079int callCbc(const char * input2)
3080{
3081  {
3082    OsiClpSolverInterface solver1;
3083    return callCbc(input2,solver1);
3084  }
3085}
3086int callCbc(const std::string input2, OsiClpSolverInterface& solver1) 
3087{
3088  char * input3 = strdup(input2.c_str());
3089  int returnCode=callCbc(input3,solver1);
3090  free(input3);
3091  return returnCode;
3092}
3093int callCbc(const std::string input2) 
3094{
3095  char * input3 = strdup(input2.c_str());
3096  OsiClpSolverInterface solver1;
3097  int returnCode=callCbc(input3,solver1);
3098  free(input3);
3099  return returnCode;
3100}
3101static int dummyCallBack(CbcModel * /*model*/, int /*whereFrom*/)
3102{
3103  return 0;
3104}
3105int CbcMain1 (int argc, const char *argv[],
3106              CbcModel  & model)
3107{
3108  return CbcMain1(argc,argv,model,dummyCallBack);
3109}
3110int callCbc1(const std::string input2, CbcModel & babSolver, int callBack(CbcModel * currentSolver, int whereFrom))
3111{
3112  char * input3 = strdup(input2.c_str());
3113  int returnCode=callCbc1(input3,babSolver,callBack);
3114  free(input3);
3115  return returnCode;
3116}
3117int callCbc1(const char * input2, CbcModel & model)
3118{
3119  return callCbc1(input2,model,dummyCallBack);
3120}
3121int CbcMain (int argc, const char *argv[],
3122             CbcModel  & model)
3123{
3124  CbcMain0(model);
3125  return CbcMain1(argc,argv,model);
3126}
3127#define CBCMAXPARAMETERS 200
3128static CbcOrClpParam parameters[CBCMAXPARAMETERS];
3129static int numberParameters=0 ;
3130void CbcMain0 (CbcModel  & model)
3131{
3132#ifndef CBC_OTHER_SOLVER
3133  OsiClpSolverInterface * originalSolver = dynamic_cast<OsiClpSolverInterface *> (model.solver());
3134#elif CBC_OTHER_SOLVER==1
3135  OsiCpxSolverInterface * originalSolver = dynamic_cast<OsiCpxSolverInterface *> (model.solver());
3136  // Dummy solvers
3137  OsiClpSolverInterface dummySolver;
3138  ClpSimplex * lpSolver = dummySolver.getModelPtr();
3139  OsiCpxSolverInterface * clpSolver = originalSolver;
3140#endif
3141  assert (originalSolver);
3142  CoinMessageHandler * generalMessageHandler = originalSolver->messageHandler();
3143  generalMessageHandler->setPrefix(true);
3144#ifndef CBC_OTHER_SOLVER
3145  OsiSolverInterface * solver = model.solver();
3146  OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
3147  ClpSimplex * lpSolver = clpSolver->getModelPtr();
3148  lpSolver->setPerturbation(50);
3149  lpSolver->messageHandler()->setPrefix(false);
3150#endif
3151  establishParams(numberParameters,parameters) ;
3152  const char dirsep =  CoinFindDirSeparator();
3153  std::string directory;
3154  std::string dirSample;
3155  std::string dirNetlib;
3156  std::string dirMiplib;
3157  if (dirsep == '/') {
3158    directory = "./";
3159    dirSample = "../../Data/Sample/";
3160    dirNetlib = "../../Data/Netlib/";
3161    dirMiplib = "../../Data/miplib3/";
3162  } else {
3163    directory = ".\\";
3164    dirSample = "..\\..\\Data\\Sample\\";
3165    dirNetlib = "..\\..\\Data\\Netlib\\";
3166    dirMiplib = "..\\..\\Data\\miplib3\\";
3167  }
3168  std::string defaultDirectory = directory;
3169  std::string importFile ="";
3170  std::string exportFile ="default.mps";
3171  std::string importBasisFile ="";
3172  std::string importPriorityFile ="";
3173  std::string debugFile="";
3174  std::string printMask="";
3175  std::string exportBasisFile ="default.bas";
3176  std::string saveFile ="default.prob";
3177  std::string restoreFile ="default.prob";
3178  std::string solutionFile ="stdout";
3179  std::string solutionSaveFile ="solution.file";
3180  int doIdiot=-1;
3181  int outputFormat=2;
3182  int substitution=3;
3183  int dualize=3;
3184  int preSolve=5;
3185  int doSprint=-1;
3186  int testOsiParameters=-1;
3187  parameters[whichParam(BASISIN,numberParameters,parameters)].setStringValue(importBasisFile);
3188  parameters[whichParam(PRIORITYIN,numberParameters,parameters)].setStringValue(importPriorityFile);
3189  parameters[whichParam(BASISOUT,numberParameters,parameters)].setStringValue(exportBasisFile);
3190  parameters[whichParam(DEBUG,numberParameters,parameters)].setStringValue(debugFile);
3191  parameters[whichParam(PRINTMASK,numberParameters,parameters)].setStringValue(printMask);
3192  parameters[whichParam(DIRECTORY,numberParameters,parameters)].setStringValue(directory);
3193  parameters[whichParam(DIRSAMPLE,numberParameters,parameters)].setStringValue(dirSample);
3194  parameters[whichParam(DIRNETLIB,numberParameters,parameters)].setStringValue(dirNetlib);
3195  parameters[whichParam(DIRMIPLIB,numberParameters,parameters)].setStringValue(dirMiplib);
3196  parameters[whichParam(DUALBOUND,numberParameters,parameters)].setDoubleValue(lpSolver->dualBound());
3197  parameters[whichParam(DUALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->dualTolerance());
3198  parameters[whichParam(EXPORT,numberParameters,parameters)].setStringValue(exportFile);
3199  parameters[whichParam(IDIOT,numberParameters,parameters)].setIntValue(doIdiot);
3200  parameters[whichParam(IMPORT,numberParameters,parameters)].setStringValue(importFile);
3201  parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].setDoubleValue(1.0e-8);
3202  int slog = whichParam(SOLVERLOGLEVEL,numberParameters,parameters);
3203  int log = whichParam(LOGLEVEL,numberParameters,parameters);
3204  parameters[slog].setIntValue(1);
3205  clpSolver->messageHandler()->setLogLevel(1) ;
3206  model.messageHandler()->setLogLevel(1);
3207  lpSolver->setLogLevel(1);
3208  parameters[log].setIntValue(1);
3209  parameters[whichParam(MAXFACTOR,numberParameters,parameters)].setIntValue(lpSolver->factorizationFrequency());
3210  parameters[whichParam(MAXITERATION,numberParameters,parameters)].setIntValue(lpSolver->maximumIterations());
3211  parameters[whichParam(OUTPUTFORMAT,numberParameters,parameters)].setIntValue(outputFormat);
3212  parameters[whichParam(PRESOLVEPASS,numberParameters,parameters)].setIntValue(preSolve);
3213  parameters[whichParam(PERTVALUE,numberParameters,parameters)].setIntValue(lpSolver->perturbation());
3214  parameters[whichParam(PRIMALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->primalTolerance());
3215  parameters[whichParam(PRIMALWEIGHT,numberParameters,parameters)].setDoubleValue(lpSolver->infeasibilityCost());
3216  parameters[whichParam(RESTORE,numberParameters,parameters)].setStringValue(restoreFile);
3217  parameters[whichParam(SAVE,numberParameters,parameters)].setStringValue(saveFile);
3218  //parameters[whichParam(TIMELIMIT,numberParameters,parameters)].setDoubleValue(1.0e8);
3219  parameters[whichParam(TIMELIMIT_BAB,numberParameters,parameters)].setDoubleValue(1.0e8);
3220  parameters[whichParam(SOLUTION,numberParameters,parameters)].setStringValue(solutionFile);
3221  parameters[whichParam(SAVESOL,numberParameters,parameters)].setStringValue(solutionSaveFile);
3222  parameters[whichParam(SPRINT,numberParameters,parameters)].setIntValue(doSprint);
3223  parameters[whichParam(SUBSTITUTION,numberParameters,parameters)].setIntValue(substitution);
3224  parameters[whichParam(DUALIZE,numberParameters,parameters)].setIntValue(dualize);
3225  model.setNumberBeforeTrust(10);
3226  parameters[whichParam(NUMBERBEFORE,numberParameters,parameters)].setIntValue(5);
3227  parameters[whichParam(MAXNODES,numberParameters,parameters)].setIntValue(model.getMaximumNodes());
3228  model.setNumberStrong(5);
3229  parameters[whichParam(STRONGBRANCHING,numberParameters,parameters)].setIntValue(model.numberStrong());
3230  parameters[whichParam(INFEASIBILITYWEIGHT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcInfeasibilityWeight));
3231  parameters[whichParam(INTEGERTOLERANCE,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcIntegerTolerance));
3232  parameters[whichParam(INCREMENT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcCutoffIncrement));
3233  parameters[whichParam(TESTOSI,numberParameters,parameters)].setIntValue(testOsiParameters);
3234  parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].setIntValue(1003);
3235  initialPumpTune=1003;
3236#ifdef CBC_THREAD
3237  parameters[whichParam(THREADS,numberParameters,parameters)].setIntValue(0);
3238#endif
3239  // Set up likely cut generators and defaults
3240  parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("sos");
3241  parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(1057);
3242  parameters[whichParam(CUTPASSINTREE,numberParameters,parameters)].setIntValue(1);
3243  parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].setIntValue(-1);
3244  parameters[whichParam(MAXHOTITS,numberParameters,parameters)].setIntValue(100);
3245  parameters[whichParam(CUTSSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
3246  parameters[whichParam(HEURISTICSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
3247  parameters[whichParam(NODESTRATEGY,numberParameters,parameters)].setCurrentOption("fewest");
3248  parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
3249  parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
3250  parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
3251#ifdef ZERO_HALF_CUTS
3252  parameters[whichParam(ZEROHALFCUTS,numberParameters,parameters)].setCurrentOption("off");
3253#endif
3254  parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption("off");
3255  parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption("ifmove");
3256  parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
3257  parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
3258  parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption("root");
3259  parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption("off");
3260  parameters[whichParam(RESIDCUTS,numberParameters,parameters)].setCurrentOption("off");
3261  parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption("on");
3262  parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption("on");
3263  parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption("on");
3264  parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
3265  parameters[whichParam(CROSSOVER2,numberParameters,parameters)].setCurrentOption("off");
3266  parameters[whichParam(PIVOTANDFIX,numberParameters,parameters)].setCurrentOption("off");
3267  parameters[whichParam(RANDROUND,numberParameters,parameters)].setCurrentOption("off");
3268  parameters[whichParam(NAIVE,numberParameters,parameters)].setCurrentOption("off");
3269  parameters[whichParam(RINS,numberParameters,parameters)].setCurrentOption("off");
3270  parameters[whichParam(DINS,numberParameters,parameters)].setCurrentOption("off");
3271  parameters[whichParam(RENS,numberParameters,parameters)].setCurrentOption("off");
3272  parameters[whichParam(LOCALTREE,numberParameters,parameters)].setCurrentOption("off");
3273  parameters[whichParam(COSTSTRATEGY,numberParameters,parameters)].setCurrentOption("off");
3274}
3275#endif
3276/* 1 - add heuristics to model
3277   2 - do heuristics (and set cutoff and best solution)
3278   3 - for miplib test so skip some
3279*/
3280#if NEW_STYLE_SOLVER==0
3281static int doHeuristics(CbcModel * model,int type)
3282#else
3283int 
3284  CbcSolver::doHeuristics(CbcModel * model,int type)
3285#endif
3286{
3287#if NEW_STYLE_SOLVER==0
3288  CbcOrClpParam * parameters_ = parameters;
3289  int numberParameters_ = numberParameters;
3290  bool noPrinting_ = noPrinting;
3291#endif
3292  char generalPrint[10000];
3293  CoinMessages generalMessages = model->messages();
3294  CoinMessageHandler * generalMessageHandler = model->messageHandler();
3295  //generalMessageHandler->setPrefix(false);
3296  bool anyToDo=false;
3297  int logLevel = parameters_[whichParam(LOGLEVEL,numberParameters_,parameters_)].intValue();
3298  int useFpump = parameters_[whichParam(FPUMP,numberParameters_,parameters_)].currentOptionAsInteger();
3299  int useRounding = parameters_[whichParam(ROUNDING,numberParameters_,parameters_)].currentOptionAsInteger();
3300  int useGreedy = parameters_[whichParam(GREEDY,numberParameters_,parameters_)].currentOptionAsInteger();
3301  int useCombine = parameters_[whichParam(COMBINE,numberParameters_,parameters_)].currentOptionAsInteger();
3302  int useCrossover = parameters_[whichParam(CROSSOVER2,numberParameters_,parameters_)].currentOptionAsInteger();
3303  int usePivot = parameters_[whichParam(PIVOTANDFIX,numberParameters_,parameters_)].currentOptionAsInteger();
3304  int useRand = parameters_[whichParam(RANDROUND,numberParameters_,parameters_)].currentOptionAsInteger();
3305  int useRINS = parameters_[whichParam(RINS,numberParameters_,parameters_)].currentOptionAsInteger();
3306  int useRENS = parameters_[whichParam(RENS,numberParameters_,parameters_)].currentOptionAsInteger();
3307  int useDINS = parameters_[whichParam(DINS,numberParameters_,parameters_)].currentOptionAsInteger();
3308  int useDIVING2 = parameters_[whichParam(DIVINGS,numberParameters_,parameters_)].currentOptionAsInteger();
3309  int useNaive = parameters_[whichParam(NAIVE,numberParameters_,parameters_)].currentOptionAsInteger();
3310  int kType = (type<10) ? type : 1;
3311  assert (kType==1||kType==2);
3312  // FPump done first as it only works if no solution
3313  if (useFpump>=kType&&useFpump<=kType+1) {
3314    anyToDo=true;
3315    CbcHeuristicFPump heuristic4(*model);
3316    double dextra3 = parameters_[whichParam(SMALLBAB,numberParameters_,parameters_)].doubleValue();
3317    heuristic4.setFractionSmall(dextra3);
3318    double dextra1 = parameters_[whichParam(ARTIFICIALCOST,numberParameters_,parameters_)].doubleValue();
3319    if (dextra1)
3320      heuristic4.setArtificialCost(dextra1);
3321    heuristic4.setMaximumPasses(parameters_[whichParam(FPUMPITS,numberParameters_,parameters_)].intValue());
3322    if (parameters_[whichParam(FPUMPITS,numberParameters_,parameters_)].intValue()==21)
3323      heuristic4.setIterationRatio(1.0);
3324    int pumpTune=parameters_[whichParam(FPUMPTUNE,numberParameters_,parameters_)].intValue();
3325    int pumpTune2=parameters_[whichParam(FPUMPTUNE2,numberParameters_,parameters_)].intValue();
3326    if (pumpTune>0) {
3327      bool printStuff = (pumpTune!=initialPumpTune||logLevel>1)
3328        &&!noPrinting_;
3329      if (printStuff) {
3330        generalMessageHandler->message(CBC_GENERAL,generalMessages)
3331          << "Options for feasibility pump - "
3332          <<CoinMessageEol;
3333      }
3334      /*
3335        >=10000000 for using obj
3336        >=1000000 use as accumulate switch
3337        >=1000 use index+1 as number of large loops
3338        >=100 use dextra1 as cutoff
3339        %100 == 10,20 etc for experimentation
3340        1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
3341        4 and static continuous, 5 as 3 but no internal integers
3342        6 as 3 but all slack basis!
3343      */
3344      double value = model->solver()->getObjSense()*model->solver()->getObjValue();
3345      int w = pumpTune/10;
3346      int i = w % 10;
3347      w /= 10;
3348      int c = w % 10;
3349      w /= 10;
3350      int r = w;
3351      int accumulate = r/1000;
3352      r -= 1000*accumulate;
3353      if (accumulate>=10) {
3354        int which = accumulate/10;
3355        accumulate -= 10*which;
3356        which--;
3357        // weights and factors
3358        double weight[]={0.01,0.01,0.1,0.1,0.5,0.5,1.0,1.0,5.0,5.0};
3359        double factor[] = {0.1,0.5,0.1,0.5,0.1,0.5,0.1,0.5,0.1,0.5};
3360        heuristic4.setInitialWeight(weight[which]);
3361        heuristic4.setWeightFactor(factor[which]);
3362        if (printStuff) {
3363          sprintf(generalPrint,"Initial weight for objective %g, decay factor %g",
3364                  weight[which],factor[which]);
3365          generalMessageHandler->message(CBC_GENERAL,generalMessages)
3366            << generalPrint
3367            <<CoinMessageEol;
3368        }
3369         
3370      }
3371      // fake cutoff
3372      if (c) {
3373        double cutoff;
3374        model->solver()->getDblParam(OsiDualObjectiveLimit,cutoff);
3375        cutoff = CoinMin(cutoff,value + 0.05*fabs(value)*c);
3376        double fakeCutoff = parameters_[whichParam(FAKECUTOFF,numberParameters_,parameters_)].doubleValue();
3377        if (fakeCutoff)
3378          cutoff=fakeCutoff;
3379        heuristic4.setFakeCutoff(cutoff);
3380        if (printStuff) {
3381          sprintf(generalPrint,"Fake cutoff of %g",cutoff);
3382          generalMessageHandler->message(CBC_GENERAL,generalMessages)
3383            << generalPrint
3384            <<CoinMessageEol;
3385        }
3386      }
3387      int offRandomEtc=0;
3388      if (pumpTune2) {
3389        if ((pumpTune2/100)!=0) {
3390          offRandomEtc=100*(pumpTune2/100);
3391          if (printStuff) {
3392            generalMessageHandler->message(CBC_GENERAL,generalMessages)
3393              << "Not using randomized objective"
3394              <<CoinMessageEol;
3395          }
3396        }
3397        int maxAllowed=pumpTune2%100;
3398        if (maxAllowed) {
3399          offRandomEtc += 1000*maxAllowed;
3400          if (printStuff) {
3401            sprintf(generalPrint,"Fixing if same for %d passes",
3402                    maxAllowed);
3403            generalMessageHandler->message(CBC_GENERAL,generalMessages)
3404              << generalPrint
3405              <<CoinMessageEol;
3406          }
3407        }
3408      }
3409      if (r) {
3410        // also set increment
3411        //double increment = (0.01*i+0.005)*(fabs(value)+1.0e-12);
3412        double increment = 0.0;
3413        double fakeIncrement = parameters_[whichParam(FAKEINCREMENT,numberParameters_,parameters_)].doubleValue();
3414        if (fakeIncrement)
3415          increment = fakeIncrement;
3416        heuristic4.setAbsoluteIncrement(increment);
3417        heuristic4.setAccumulate(accumulate);
3418        heuristic4.setMaximumRetries(r+1);
3419        if (printStuff) {
3420          if (increment) {
3421            sprintf(generalPrint,"Increment of %g",increment);
3422            generalMessageHandler->message(CBC_GENERAL,generalMessages)
3423              << generalPrint
3424              <<CoinMessageEol;
3425          }
3426          if (accumulate) {
3427            sprintf(generalPrint,"Accumulate of %d",accumulate);
3428            generalMessageHandler->message(CBC_GENERAL,generalMessages)
3429              << generalPrint
3430              <<CoinMessageEol;
3431          }
3432          sprintf(generalPrint,"%d retries",r+1);
3433          generalMessageHandler->message(CBC_GENERAL,generalMessages)
3434            << generalPrint
3435            <<CoinMessageEol;
3436        }
3437      }
3438      if (i+offRandomEtc) {
3439        heuristic4.setFeasibilityPumpOptions(i*10+offRandomEtc);
3440        if (printStuff) {
3441          sprintf(generalPrint,"Feasibility pump options of %d",
3442                  i*10+offRandomEtc);
3443          generalMessageHandler->message(CBC_GENERAL,generalMessages)
3444            << generalPrint
3445            <<CoinMessageEol;
3446        }
3447      } 
3448      pumpTune = pumpTune%100;
3449      if (pumpTune==6)
3450        pumpTune =13;
3451      heuristic4.setWhen((pumpTune%10)+10);
3452      if (printStuff) {
3453        sprintf(generalPrint,"Tuning (fixing) %d",pumpTune%10);
3454        generalMessageHandler->message(CBC_GENERAL,generalMessages)
3455          << generalPrint
3456          <<CoinMessageEol;
3457      }
3458    }
3459    heuristic4.setHeuristicName("feasibility pump");
3460    //#define ROLF
3461#ifdef ROLF   
3462    CbcHeuristicFPump pump(*model);
3463    pump.setMaximumTime(60);
3464    pump.setMaximumPasses(100);
3465    pump.setMaximumRetries(1);
3466    pump.setFixOnReducedCosts(0);
3467    pump.setHeuristicName("Feasibility pump");
3468    pump.setFractionSmall(1.0);
3469    pump.setWhen(13);
3470    model->addHeuristic(&pump);
3471#else
3472    model->addHeuristic(&heuristic4);
3473#endif
3474  }
3475  if (useRounding>=type&&useRounding>=kType&&useRounding<=kType+1) {
3476    CbcRounding heuristic1(*model);
3477    heuristic1.setHeuristicName("rounding");
3478    model->addHeuristic(&heuristic1) ;
3479    anyToDo=true;
3480  }
3481  if (useGreedy>=type&&useGreedy>=kType&&useGreedy<=kType+1) {
3482    CbcHeuristicGreedyCover heuristic3(*model);
3483    heuristic3.setHeuristicName("greedy cover");
3484    CbcHeuristicGreedyEquality heuristic3a(*model);
3485    heuristic3a.setHeuristicName("greedy equality");
3486    model->addHeuristic(&heuristic3);
3487    model->addHeuristic(&heuristic3a);
3488    anyToDo=true;
3489  }
3490  if (useRENS>=kType&&useRENS<=kType+1) {
3491    CbcHeuristicRENS heuristic6(*model);
3492    heuristic6.setHeuristicName("RENS");
3493    heuristic6.setFractionSmall(0.4);
3494    heuristic6.setFeasibilityPumpOptions(1008003);
3495    int nodes []={-2,50,50,50,200,1000,10000};
3496    heuristic6.setNumberNodes(nodes[useRENS]);
3497    model->addHeuristic(&heuristic6) ;
3498    anyToDo=true;
3499  }
3500  if (useNaive>=kType&&useNaive<=kType+1) {
3501    CbcHeuristicNaive heuristic5b(*model);
3502    heuristic5b.setHeuristicName("Naive");
3503    heuristic5b.setFractionSmall(0.4);
3504    heuristic5b.setNumberNodes(50);
3505    model->addHeuristic(&heuristic5b) ;
3506    anyToDo=true;
3507  }
3508  int useDIVING=0;
3509  {
3510    int useD;
3511    useD = parameters_[whichParam(DIVINGV,numberParameters_,parameters_)].currentOptionAsInteger();
3512    useDIVING |= 1*((useD>=kType) ? 1 : 0);
3513    useD = parameters_[whichParam(DIVINGG,numberParameters_,parameters_)].currentOptionAsInteger();
3514    useDIVING |= 2*((useD>=kType) ? 1 : 0);
3515    useD = parameters_[whichParam(DIVINGF,numberParameters_,parameters_)].currentOptionAsInteger();
3516    useDIVING |= 4*((useD>=kType) ? 1 : 0);
3517    useD = parameters_[whichParam(DIVINGC,numberParameters_,parameters_)].currentOptionAsInteger();
3518    useDIVING |= 8*((useD>=kType) ? 1 : 0);
3519    useD = parameters_[whichParam(DIVINGL,numberParameters_,parameters_)].currentOptionAsInteger();
3520    useDIVING |= 16*((useD>=kType) ? 1 : 0);
3521    useD = parameters_[whichParam(DIVINGP,numberParameters_,parameters_)].currentOptionAsInteger();
3522    useDIVING |= 32*((useD>=kType) ? 1 : 0);
3523  }
3524  if (useDIVING2>=kType&&useDIVING2<=kType+1) {
3525    int diveOptions=parameters_[whichParam(DIVEOPT,numberParameters_,parameters_)].intValue();
3526    if (diveOptions<0||diveOptions>10)
3527      diveOptions=2;
3528    CbcHeuristicJustOne heuristicJustOne(*model);
3529    heuristicJustOne.setHeuristicName("DiveAny");
3530    heuristicJustOne.setWhen(diveOptions);
3531    // add in others
3532    CbcHeuristicDiveCoefficient heuristicDC(*model);
3533    heuristicDC.setHeuristicName("DiveCoefficient");
3534    heuristicJustOne.addHeuristic(&heuristicDC,1.0) ;
3535    CbcHeuristicDiveFractional heuristicDF(*model);
3536    heuristicDF.setHeuristicName("DiveFractional");
3537    heuristicJustOne.addHeuristic(&heuristicDF,1.0) ;
3538    CbcHeuristicDiveGuided heuristicDG(*model);
3539    heuristicDG.setHeuristicName("DiveGuided");
3540    heuristicJustOne.addHeuristic(&heuristicDG,1.0) ;
3541    CbcHeuristicDiveLineSearch heuristicDL(*model);
3542    heuristicDL.setHeuristicName("DiveLineSearch");
3543    heuristicJustOne.addHeuristic(&heuristicDL,1.0) ;
3544    CbcHeuristicDivePseudoCost heuristicDP(*model);
3545    heuristicDP.setHeuristicName("DivePseudoCost");
3546    heuristicJustOne.addHeuristic(&heuristicDP,1.0) ;
3547    CbcHeuristicDiveVectorLength heuristicDV(*model);
3548    heuristicDV.setHeuristicName("DiveVectorLength");
3549    heuristicJustOne.addHeuristic(&heuristicDV,1.0) ;
3550    // Now normalize probabilities
3551    heuristicJustOne.normalizeProbabilities();
3552    model->addHeuristic(&heuristicJustOne) ;
3553  }
3554 
3555  if (useDIVING>0) {
3556    int diveOptions=parameters_[whichParam(DIVEOPT,numberParameters_,parameters_)].intValue();
3557    if (diveOptions<0||diveOptions>10)
3558      diveOptions=2;
3559    if ((useDIVING&1)!=0) {
3560      CbcHeuristicDiveVectorLength heuristicDV(*model);
3561      heuristicDV.setHeuristicName("DiveVectorLength");
3562      heuristicDV.setWhen(diveOptions);
3563      model->addHeuristic(&heuristicDV) ;
3564    }
3565    if ((useDIVING&2)!=0) {
3566      CbcHeuristicDiveGuided heuristicDG(*model);
3567      heuristicDG.setHeuristicName("DiveGuided");
3568      heuristicDG.setWhen(diveOptions);
3569      model->addHeuristic(&heuristicDG) ;
3570    }
3571    if ((useDIVING&4)!=0) {
3572      CbcHeuristicDiveFractional heuristicDF(*model);
3573      heuristicDF.setHeuristicName("DiveFractional");
3574      heuristicDF.setWhen(diveOptions);
3575      model->addHeuristic(&heuristicDF) ;
3576    }
3577    if ((useDIVING&8)!=0) {
3578      CbcHeuristicDiveCoefficient heuristicDC(*model);
3579      heuristicDC.setHeuristicName("DiveCoefficient");
3580      heuristicDC.setWhen(diveOptions);
3581      model->addHeuristic(&heuristicDC) ;
3582    }
3583    if ((useDIVING&16)!=0) {
3584      CbcHeuristicDiveLineSearch heuristicDL(*model);
3585      heuristicDL.setHeuristicName("DiveLineSearch");
3586      heuristicDL.setWhen(diveOptions);
3587      model->addHeuristic(&heuristicDL) ;
3588    }
3589    if ((useDIVING&32)!=0) {
3590      CbcHeuristicDivePseudoCost heuristicDP(*model);
3591      heuristicDP.setHeuristicName("DivePseudoCost");
3592      heuristicDP.setWhen(diveOptions);
3593      model->addHeuristic(&heuristicDP) ;
3594    }
3595    anyToDo=true;
3596  }
3597  if (usePivot>=type&&usePivot<=kType+1) {
3598    CbcHeuristicPivotAndFix heuristic(*model);
3599    heuristic.setHeuristicName("pivot and fix");
3600    heuristic.setFractionSmall(10.0); // normally 0.5
3601    model->addHeuristic(&heuristic);
3602    anyToDo=true;
3603  }
3604  if (useRand>=type&&useRand<=kType+1) {
3605    CbcHeuristicRandRound heuristic(*model);
3606    heuristic.setHeuristicName("randomized rounding");
3607    heuristic.setFractionSmall(10.0); // normally 0.5
3608    model->addHeuristic(&heuristic);
3609    anyToDo=true;
3610  }
3611  if (useDINS>=kType&&useDINS<=kType+1) {
3612    CbcHeuristicDINS heuristic5a(*model);
3613    heuristic5a.setHeuristicName("DINS");
3614    heuristic5a.setFractionSmall(0.6);
3615    if (useDINS<4)
3616      heuristic5a.setDecayFactor(5.0);
3617    else
3618      heuristic5a.setDecayFactor(1.5);
3619    heuristic5a.setNumberNodes(1000);
3620    model->addHeuristic(&heuristic5a) ;
3621    anyToDo=true;
3622  }
3623  if (useRINS>=kType&&useRINS<=kType+1) {
3624    CbcHeuristicRINS heuristic5(*model);
3625    heuristic5.setHeuristicName("RINS");
3626    if (useRINS<4) {
3627      heuristic5.setFractionSmall(0.5);
3628      heuristic5.setDecayFactor(5.0);
3629    } else {
3630      heuristic5.setFractionSmall(0.6);
3631      heuristic5.setDecayFactor(1.5);
3632    }
3633    model->addHeuristic(&heuristic5) ;
3634    anyToDo=true;
3635  }
3636  if (useCombine>=type&&useCombine>=kType&&useCombine<=kType+1) {
3637    CbcHeuristicLocal heuristic2(*model);
3638    heuristic2.setHeuristicName("combine solutions");
3639    heuristic2.setFractionSmall(0.5);
3640    heuristic2.setSearchType(1);
3641    model->addHeuristic(&heuristic2);
3642    anyToDo=true;
3643  }
3644  if (useCrossover>=kType&&useCrossover<=kType+1) {
3645    CbcHeuristicCrossover heuristic2a(*model);
3646    heuristic2a.setHeuristicName("crossover");
3647    heuristic2a.setFractionSmall(0.3);
3648    // just fix at lower
3649    heuristic2a.setWhen(11);
3650    model->addHeuristic(&heuristic2a);
3651    model->setMaximumSavedSolutions(5);
3652    anyToDo=true;
3653  }
3654  int heurSwitches=parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()%100;
3655  if (heurSwitches) {
3656    for (int iHeur=0;iHeur<model->numberHeuristics();iHeur++) {
3657      CbcHeuristic * heuristic = model->heuristic(iHeur);
3658      heuristic->setSwitches(heurSwitches);
3659    }
3660  }
3661  if (type==2&&anyToDo) {
3662    // Do heuristics
3663#if 1
3664    // clean copy
3665    CbcModel model2(*model);
3666    // But get rid of heuristics in model
3667    model->doHeuristicsAtRoot(2);
3668    if (logLevel<=1)
3669      model2.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
3670    OsiBabSolver defaultC;
3671    //solver_->setAuxiliaryInfo(&defaultC);
3672    model2.passInSolverCharacteristics(&defaultC);
3673    // Save bounds
3674    int numberColumns = model2.solver()->getNumCols();
3675    model2.createContinuousSolver();
3676    bool cleanModel = !model2.numberIntegers()&&!model2.numberObjects();
3677    model2.findIntegers(false);
3678    int heurOptions=(parameters_[whichParam(HOPTIONS,numberParameters_,parameters_)].intValue()/100)%100;
3679    if (heurOptions==0||heurOptions==2) {
3680      model2.doHeuristicsAtRoot(1);
3681    } else if (heurOptions==1||heurOptions==3) {
3682      model2.setMaximumNodes(-1);
3683      CbcStrategyDefault strategy(0,5,5);
3684      strategy.setupPreProcessing(1,0);
3685      model2.setStrategy(strategy);
3686      model2.branchAndBound();
3687    }
3688    if (cleanModel)
3689      model2.zapIntegerInformation(false);
3690    if (model2.bestSolution()) {
3691      double value = model2.getMinimizationObjValue();
3692      model->setCutoff(value);
3693      model->setBestSolution(model2.bestSolution(),numberColumns,value);
3694      model->setSolutionCount(1);
3695      model->setNumberHeuristicSolutions(1);
3696    }
3697#else
3698    if (logLevel<=1)
3699      model->solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
3700    OsiBabSolver defaultC;
3701    //solver_->setAuxiliaryInfo(&defaultC);
3702    model->passInSolverCharacteristics(&defaultC);
3703    // Save bounds
3704    int numberColumns = model->solver()->getNumCols();
3705    model->createContinuousSolver();
3706    bool cleanModel = !model->numberIntegers()&&!model->numberObjects();
3707    model->findIntegers(false);
3708    model->doHeuristicsAtRoot(1);
3709    if (cleanModel)
3710      model->zapIntegerInformation(false);
3711#endif
3712    return 0;
3713  } else {
3714    return 0;
3715  }
3716} 
3717
3718int CbcClpUnitTest (const CbcModel & saveModel,
3719                    std::string& dirMiplib, int testSwitch,
3720                    double * stuff);
3721#if NEW_STYLE_SOLVER
3722/* This takes a list of commands, does "stuff" and returns
3723   returnMode -
3724   0 model and solver untouched - babModel updated
3725   1 model updated - just with solution basis etc
3726   2 model updated i.e. as babModel (babModel NULL) (only use without preprocessing)
3727*/
3728int 
3729CbcSolver::solve(const char * input2, int returnMode)
3730{
3731  char * input = strdup(input2);
3732  int length = strlen(input);
3733  bool blank = input[0]=='0';
3734  int n=blank ? 0 : 1;
3735  for (int i=0;i<length;i++) {
3736    if (blank) {
3737      // look for next non blank
3738      if (input[i]==' ') {
3739        continue;
3740      } else {
3741        n++;
3742        blank=false;
3743      }
3744    } else {
3745      // look for next blank
3746      if (input[i]!=' ') {
3747        continue;
3748      } else {
3749        blank=true;
3750      }
3751    }
3752  }
3753  char ** argv = new char * [n+2];
3754  argv[0]=strdup("cbc");
3755  int i=0;
3756  while(input[i]==' ')
3757    i++;
3758  for (int j=0;j<n;j++) {
3759    int saveI=i;
3760    for (;i<length;i++) {
3761      // look for next blank
3762      if (input[i]!=' ') {
3763        continue;
3764      } else {
3765        break;
3766      }
3767    }
3768    char save = input[i];
3769    input[i]='\0';
3770    argv[j+1]=strdup(input+saveI);
3771    input[i]=save;
3772    while(input[i]==' ')
3773      i++;
3774  }
3775  argv[n+1]=strdup("-quit");
3776  free(input);
3777  int returnCode = solve(n+2,const_cast<const char **>(argv),returnMode);
3778  for (int k=0;k<n+2;k++)
3779    free(argv[k]);
3780  delete [] argv;
3781  return returnCode;
3782}
3783#endif
3784/* Meaning of whereFrom:
3785   1 after initial solve by dualsimplex etc
3786   2 after preprocessing
3787   3 just before branchAndBound (so user can override)
3788   4 just after branchAndBound (before postprocessing)
3789   5 after postprocessing
3790   6 after a user called heuristic phase
3791*/
3792
3793#if NEW_STYLE_SOLVER==0
3794int CbcMain1 (int argc, const char *argv[],
3795              CbcModel  & model,
3796              int callBack(CbcModel * currentSolver, int whereFrom))
3797#else
3798/* This takes a list of commands, does "stuff" and returns
3799   returnMode -
3800   0 model and solver untouched - babModel updated
3801   1 model updated - just with solution basis etc
3802   2 model updated i.e. as babModel (babModel NULL)
3803*/
3804int 
3805  CbcSolver::solve (int argc, const char *argv[], int returnMode)
3806#endif
3807{
3808#if NEW_STYLE_SOLVER==0
3809  CbcOrClpParam * parameters_ = parameters;
3810  int numberParameters_ = numberParameters;
3811  CbcModel & model_ = model;
3812  CbcModel * babModel_ = NULL;
3813  int returnMode=1;
3814  CbcOrClpRead_mode=1;
3815  int statusUserFunction_[1];
3816  int numberUserFunctions_=1; // to allow for ampl
3817#else
3818  delete babModel_;
3819  babModel_ = NULL;
3820  CbcOrClpRead_mode=1;
3821  delete [] statusUserFunction_;
3822  statusUserFunction_ = new int [numberUserFunctions_];
3823  int iUser;
3824#endif
3825  // Statistics
3826  double statistics_seconds=0.0, statistics_obj=0.0;
3827  double statistics_sys_seconds=0.0, statistics_elapsed_seconds=0.0;
3828  CoinWallclockTime();
3829  double statistics_continuous=0.0, statistics_tighter=0.0;
3830  double statistics_cut_time=0.0;
3831  int statistics_nodes=0, statistics_iterations=0;
3832  int statistics_nrows=0, statistics_ncols=0;
3833  int statistics_nprocessedrows=0, statistics_nprocessedcols=0;
3834  std::string statistics_result;
3835  int * statistics_number_cuts=NULL;
3836  const char ** statistics_name_generators=NULL;
3837  int statistics_number_generators=0;
3838  memset(statusUserFunction_,0,numberUserFunctions_*sizeof(int));
3839  /* Note
3840     This is meant as a stand-alone executable to do as much of coin as possible.
3841     It should only have one solver known to it.
3842  */
3843  CoinMessageHandler * generalMessageHandler = model_.messageHandler();
3844  generalMessageHandler->setPrefix(false);
3845#ifndef CBC_OTHER_SOLVER
3846  OsiClpSolverInterface * originalSolver = dynamic_cast<OsiClpSolverInterface *> (model_.solver());
3847  assert (originalSolver);
3848  // Move handler across if not default
3849  if (!originalSolver->defaultHandler()&&originalSolver->getModelPtr()->defaultHandler())
3850    originalSolver->getModelPtr()->passInMessageHandler(originalSolver->messageHandler());
3851  CoinMessages generalMessages = originalSolver->getModelPtr()->messages();
3852  char generalPrint[10000];
3853  if (originalSolver->getModelPtr()->logLevel()==0)
3854    noPrinting=true;
3855#elif CBC_OTHER_SOLVER==1
3856  OsiCpxSolverInterface * originalSolver = dynamic_cast<OsiCpxSolverInterface *> (model_.solver());
3857  assert (originalSolver);
3858  OsiClpSolverInterface dummySolver;
3859  OsiCpxSolverInterface * clpSolver = originalSolver;
3860  CoinMessages generalMessages = dummySolver.getModelPtr()->messages();
3861  char generalPrint[10000];
3862  noPrinting=true;
3863#endif
3864#if NEW_STYLE_SOLVER==0
3865  bool noPrinting_=noPrinting;
3866#endif
3867  // Say not in integer
3868  int integerStatus=-1;
3869  // Say no resolve after cuts
3870  model_.setResolveAfterTakeOffCuts(false);
3871  // see if log in list
3872  for (int i=1;i<argc;i++) {
3873    if (!strncmp(argv[i],"log",3)) {
3874      const char * equals = strchr(argv[i],'=');
3875      if (equals&&atoi(equals+1)!=0) 
3876        noPrinting_=false;
3877      else
3878        noPrinting_=true;
3879      break;
3880    } else if (!strncmp(argv[i],"-log",4)&&i<argc-1) {
3881      if (atoi(argv[i+1])!=0) 
3882        noPrinting_=false;
3883      else
3884        noPrinting_=true;
3885      break;
3886    }
3887  }
3888  double time0;
3889  {
3890    double time1 = CoinCpuTime(),time2;
3891    time0=time1;
3892    bool goodModel=(originalSolver->getNumCols()) ? true : false;
3893
3894    // register signal handler
3895    //CoinSighandler_t saveSignal=signal(SIGINT,signal_handler);
3896    signal(SIGINT,signal_handler);
3897    // Set up all non-standard stuff
3898    int cutPass=-1234567;
3899    int cutPassInTree=-1234567;
3900    int tunePreProcess=0;
3901    int testOsiParameters=-1;
3902    // 0 normal, 1 from ampl or MIQP etc (2 allows cuts)
3903    int complicatedInteger=0;
3904    OsiSolverInterface * solver = model_.solver();
3905    if (noPrinting_) 
3906      setCbcOrClpPrinting(false);
3907#ifndef CBC_OTHER_SOLVER
3908    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
3909    ClpSimplex * lpSolver = clpSolver->getModelPtr();
3910    if (noPrinting_) {
3911      lpSolver->setLogLevel(0);
3912    }
3913#else
3914    ClpSimplex * lpSolver = NULL;
3915#endif
3916    // For priorities etc
3917    int * priorities=NULL;
3918    int * branchDirection=NULL;
3919    double * pseudoDown=NULL;
3920    double * pseudoUp=NULL;
3921    double * solutionIn = NULL;
3922    int * prioritiesIn = NULL;
3923    int numberSOS = 0;
3924    int * sosStart = NULL;
3925    int * sosIndices = NULL;
3926    char * sosType = NULL;
3927    double * sosReference = NULL;
3928    int * cut=NULL;
3929    int * sosPriority=NULL;
3930    CglStored storedAmpl;
3931    CoinModel * coinModel = NULL;
3932    CoinModel saveCoinModel;
3933    CoinModel saveTightenedModel;
3934    int * whichColumn = NULL;
3935    int * knapsackStart=NULL;
3936    int * knapsackRow=NULL;
3937    int numberKnapsack=0;
3938#if NEW_STYLE_SOLVER
3939    int numberInputs=0;
3940    readMode_=CbcOrClpRead_mode;
3941    for (iUser=0;iUser<numberUserFunctions_;iUser++) {
3942      int status = userFunction_[iUser]->importData(this,argc,const_cast<char **>(argv));
3943      if (status>=0) {
3944        if (!status) {
3945          numberInputs++;
3946          statusUserFunction_[iUser]=1;
3947          goodModel=true;
3948          solver = model_.solver();
3949#ifndef CBC_OTHER_SOLVER
3950          clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
3951          lpSolver = clpSolver->getModelPtr();
3952#endif
3953        } else {
3954          printf("Bad input from user function %s\n",userFunction_[iUser]->name().c_str());
3955          abort();
3956        }
3957      }
3958    }
3959    if (numberInputs>1) {
3960      printf("Two or more user inputs!\n");
3961      abort();
3962    }
3963    if (originalCoinModel_)
3964      complicatedInteger=1;
3965    testOsiParameters = intValue(TESTOSI);
3966    if (noPrinting_) {
3967      model_.messageHandler()->setLogLevel(0);
3968      setCbcOrClpPrinting(false);
3969    }
3970    CbcOrClpRead_mode=readMode_;
3971#endif
3972#ifdef COIN_HAS_ASL
3973    ampl_info info;
3974    {
3975    memset(&info,0,sizeof(info));
3976    if (argc>2&&!strcmp(argv[2],"-AMPL")) {
3977      statusUserFunction_[0]=1;
3978      // see if log in list
3979      noPrinting_=true;
3980      for (int i=1;i<argc;i++) {
3981        if (!strncmp(argv[i],"log",3)) {
3982          const char * equals = strchr(argv[i],'=');
3983          if (equals&&atoi(equals+1)>0) {
3984            noPrinting_=false;
3985            info.logLevel=atoi(equals+1);
3986            int log = whichParam(LOGLEVEL,numberParameters_,parameters_);
3987            parameters_[log].setIntValue(info.logLevel);
3988            // mark so won't be overWritten
3989            info.numberRows=-1234567;
3990            break;
3991          }
3992        }
3993      }
3994
3995      union { void * voidModel; CoinModel * model; } coinModelStart;
3996      coinModelStart.model=NULL;
3997      int returnCode = readAmpl(&info,argc, const_cast<char **>(argv),& coinModelStart.voidModel);
3998      coinModel=coinModelStart.model;
3999      if (returnCode)
4000        return returnCode;
4001      CbcOrClpRead_mode=2; // so will start with parameters
4002      // see if log in list (including environment)
4003      for (int i=1;i<info.numberArguments;i++) {
4004        if (!strcmp(info.arguments[i],"log")) {
4005          if (i<info.numberArguments-1&&atoi(info.arguments[i+1])>0)
4006            noPrinting_=false;
4007          break;
4008        }
4009      }
4010      if (noPrinting_) {
4011        model_.messageHandler()->setLogLevel(0);
4012        setCbcOrClpPrinting(false);
4013      }
4014      if (!noPrinting_)
4015        printf("%d rows, %d columns and %d elements\n",
4016               info.numberRows,info.numberColumns,info.numberElements);
4017#ifdef COIN_HAS_LINK
4018      if (!coinModel) {
4019#endif
4020      solver->loadProblem(info.numberColumns,info.numberRows,info.starts,
4021                          info.rows,info.elements,
4022                          info.columnLower,info.columnUpper,info.objective,
4023                          info.rowLower,info.rowUpper);
4024      // take off cuts if ampl wants that
4025      if (info.cut&&0) {
4026        printf("AMPL CUTS OFF until global cuts fixed\n");
4027        info.cut=NULL;
4028      }
4029      if (info.cut) {
4030        int numberRows = info.numberRows;
4031        int * whichRow = new int [numberRows];
4032        // Row copy
4033        const CoinPackedMatrix * matrixByRow = solver->getMatrixByRow();
4034        const double * elementByRow = matrixByRow->getElements();
4035        const int * column = matrixByRow->getIndices();
4036        const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
4037        const int * rowLength = matrixByRow->getVectorLengths();
4038       
4039        const double * rowLower = solver->getRowLower();
4040        const double * rowUpper = solver->getRowUpper();
4041        int nDelete=0;
4042        for (int iRow=0;iRow<numberRows;iRow++) {
4043          if (info.cut[iRow]) {
4044            whichRow[nDelete++]=iRow;
4045            int start = rowStart[iRow];
4046            storedAmpl.addCut(rowLower[iRow],rowUpper[iRow],
4047                          rowLength[iRow],column+start,elementByRow+start);
4048          }
4049        }
4050        solver->deleteRows(nDelete,whichRow);
4051        delete [] whichRow;
4052      }
4053#ifdef COIN_HAS_LINK
4054      } else {
4055#ifndef CBC_OTHER_SOLVER
4056        // save
4057        saveCoinModel = *coinModel;
4058        // load from coin model
4059        OsiSolverLink solver1;
4060        OsiSolverInterface * solver2 = solver1.clone();
4061        model_.assignSolver(solver2,false);
4062        OsiSolverLink * si =
4063          dynamic_cast<OsiSolverLink *>(model_.solver()) ;
4064        assert (si != NULL);
4065        si->setDefaultMeshSize(0.001);
4066        // need some relative granularity
4067        si->setDefaultBound(100.0);
4068        double dextra3 = parameters_[whichParam(DEXTRA3,numberParameters_,parameters_)].doubleValue();
4069        if (dextra3)
4070          si->setDefaultMeshSize(dextra3);
4071        si->setDefaultBound(100000.0);
4072        si->setIntegerPriority(1000);
4073        si->setBiLinearPriority(10000);
4074        CoinModel * model2 = reinterpret_cast<CoinModel *> (coinModel);
4075        int logLevel = parameters_[whichParam(LOGLEVEL,numberParameters_,parameters_)].intValue();
4076        si->load(*model2,true,logLevel);
4077        // redo
4078        solver = model_.solver();
4079        clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
4080        lpSolver = clpSolver->getModelPtr();
4081        clpSolver->messageHandler()->setLogLevel(0) ;
4082        testOsiParameters=0;
4083        parameters_[whichParam(TESTOSI,numberParameters_,parameters_)].setIntValue(0);
4084        complicatedInteger=1;
4085        if (info.cut) {
4086          printf("Sorry - can't do cuts with LOS as ruins delicate row order\n");
4087          abort();
4088          int numberRows = info.numberRows;
4089          int * whichRow = new int [numberRows];
4090          // Row copy
4091          const CoinPackedMatrix * matrixByRow = solver->getMatrixByRow();
4092          const double * elementByRow = matrixByRow->getElements();
4093          const int * column = matrixByRow->getIndices();
4094          const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
4095          const int * rowLength = matrixByRow->getVectorLengths();
4096         
4097          const double * rowLower = solver->getRowLower();
4098          const double * rowUpper = solver->getRowUpper();
4099          int nDelete=0;
4100          for (int iRow=0;iRow<numberRows;iRow++) {
4101            if (info.cut[iRow]) {
4102              whichRow[nDelete++]=iRow;
4103              int start = rowStart[iRow];
4104              storedAmpl.addCut(rowLower[iRow],rowUpper[iRow],
4105                                rowLength[iRow],column+start,elementByRow+start);
4106            }
4107          }
4108          solver->deleteRows(nDelete,whichRow);
4109          // and special matrix
4110          si->cleanMatrix()->deleteRows(nDelete,whichRow);
4111          delete [] whichRow;
4112        }
4113#endif
4114      }
4115#endif
4116      // If we had a solution use it
4117      if (info.primalSolution) {
4118        solver->setColSolution(info.primalSolution);
4119      }
4120      // status
4121      if (info.rowStatus) {
4122        unsigned char * statusArray = lpSolver->statusArray();
4123        int i;
4124        for (i=0;i<info.numberColumns;i++)
4125          statusArray[i]= info.columnStatus[i];
4126        statusArray+=info.numberColumns;
4127        for (i=0;i<info.numberRows;i++)
4128          statusArray[i]= info.rowStatus[i];
4129        CoinWarmStartBasis * basis = lpSolver->getBasis();
4130        solver->setWarmStart(basis);
4131        delete basis;
4132      }
4133      freeArrays1(&info);
4134      // modify objective if necessary
4135      solver->setObjSense(info.direction);
4136      solver->setDblParam(OsiObjOffset,info.offset);
4137      if (info.offset) {
4138        sprintf(generalPrint,"Ampl objective offset is %g",
4139                info.offset);
4140        generalMessageHandler->message(CLP_GENERAL,generalMessages)
4141          << generalPrint
4142          <<CoinMessageEol;
4143      }
4144      // Set integer variables (unless nonlinear when set)
4145      if (!info.nonLinear) {
4146        for (int i=info.numberColumns-info.numberIntegers;
4147             i<info.numberColumns;i++)
4148          solver->setInteger(i);
4149      }
4150      goodModel=true;
4151      // change argc etc
4152      argc = info.numberArguments;
4153      argv = const_cast<const char **>(info.arguments);
4154    }
4155    }
4156#endif   
4157    // default action on import
4158    int allowImportErrors=0;
4159    int keepImportNames=1;
4160    int doIdiot=-1;
4161    int outputFormat=2;
4162    int slpValue=-1;
4163    int cppValue=-1;
4164    int printOptions=0;
4165    int printMode=0;
4166    int presolveOptions=0;
4167    int substitution=3;
4168    int dualize=3;
4169    int doCrash=0;
4170    int doVector=0;
4171    int doSprint=-1;
4172    int doScaling=4;
4173    // set reasonable defaults
4174    int preSolve=5;
4175    int preProcess=4;
4176    bool useStrategy=false;
4177    bool preSolveFile=false;
4178    bool strongChanged=false;
4179    bool pumpChanged=false;
4180   
4181    double djFix=1.0e100;
4182    double tightenFactor=0.0;
4183    const char dirsep =  CoinFindDirSeparator();
4184    std::string directory;
4185    std::string dirSample;
4186    std::string dirNetlib;
4187    std::string dirMiplib;
4188    if (dirsep == '/') {
4189      directory = "./";
4190      dirSample = "../../Data/Sample/";
4191      dirNetlib = "../../Data/Netlib/";
4192      dirMiplib = "../../Data/miplib3/";
4193    } else {
4194      directory = ".\\";
4195      dirSample = "..\\..\\Data\\Sample\\";
4196      dirNetlib = "..\\..\\Data\\Netlib\\";
4197      dirMiplib = "..\\..\\Data\\miplib3\\";
4198    }
4199    std::string defaultDirectory = directory;
4200    std::string importFile ="";
4201    std::string exportFile ="default.mps";
4202    std::string importBasisFile ="";
4203    std::string importPriorityFile ="";
4204    std::string debugFile="";
4205    std::string printMask="";
4206    double * debugValues = NULL;
4207    int numberDebugValues = -1;
4208    int basisHasValues=0;
4209    std::string exportBasisFile ="default.bas";
4210    std::string saveFile ="default.prob";
4211    std::string restoreFile ="default.prob";
4212    std::string solutionFile ="stdout";
4213    std::string solutionSaveFile ="solution.file";
4214    int slog = whichParam(SOLVERLOGLEVEL,numberParameters_,parameters_);
4215    int log = whichParam(LOGLEVEL,numberParameters_,parameters_);
4216#ifndef CBC_OTHER_SOLVER
4217    double normalIncrement=model_.getCutoffIncrement();;
4218#endif
4219    if (testOsiParameters>=0) {
4220      // trying nonlinear - switch off some stuff
4221      preProcess=0;
4222    }
4223    // Set up likely cut generators and defaults
4224    int nodeStrategy=0;
4225    bool dominatedCuts=false;
4226    int doSOS=1;
4227    int verbose=0;
4228    CglGomory gomoryGen;
4229    // try larger limit
4230    gomoryGen.setLimitAtRoot(1000);
4231    gomoryGen.setLimit(50);
4232    // set default action (0=off,1=on,2=root)
4233    int gomoryAction=3;
4234
4235    CglProbing probingGen;
4236    probingGen.setUsingObjective(1);
4237    probingGen.setMaxPass(1);
4238    probingGen.setMaxPassRoot(1);
4239    // Number of unsatisfied variables to look at
4240    probingGen.setMaxProbe(10);
4241    probingGen.setMaxProbeRoot(50);
4242    // How far to follow the consequences
4243    probingGen.setMaxLook(10);
4244    probingGen.setMaxLookRoot(50);
4245    probingGen.setMaxLookRoot(10);
4246    // Only look at rows with fewer than this number of elements
4247    probingGen.setMaxElements(200);
4248    probingGen.setMaxElementsRoot(300);
4249    probingGen.setRowCuts(3);
4250    // set default action (0=off,1=on,2=root)
4251    int probingAction=1;
4252
4253    CglKnapsackCover knapsackGen;
4254    //knapsackGen.switchOnExpensive();
4255    //knapsackGen.setMaxInKnapsack(100);
4256    // set default action (0=off,1=on,2=root)
4257    int knapsackAction=3;
4258
4259    CglRedSplit redsplitGen;
4260    //redsplitGen.setLimit(100);
4261    // set default action (0=off,1=on,2=root)
4262    // Off as seems to give some bad cuts
4263    int redsplitAction=0;
4264
4265    CglFakeClique cliqueGen(NULL,false);
4266    //CglClique cliqueGen(false,true);
4267    cliqueGen.setStarCliqueReport(false);
4268    cliqueGen.setRowCliqueReport(false);
4269    cliqueGen.setMinViolation(0.1);
4270    // set default action (0=off,1=on,2=root)
4271    int cliqueAction=3;
4272
4273    CglMixedIntegerRounding2 mixedGen;
4274    // set default action (0=off,1=on,2=root)
4275    int mixedAction=3;
4276
4277    CglFlowCover flowGen;
4278    // set default action (0=off,1=on,2=root)
4279    int flowAction=3;
4280
4281    CglTwomir twomirGen;
4282    twomirGen.setMaxElements(250);
4283    // set default action (0=off,1=on,2=root)
4284    int twomirAction=2;
4285#ifndef DEBUG_MALLOC
4286    CglLandP landpGen;
4287#endif
4288    // set default action (0=off,1=on,2=root)
4289    int landpAction=0;
4290    CglResidualCapacity residualCapacityGen;
4291    // set default action (0=off,1=on,2=root)
4292    int residualCapacityAction=0;
4293
4294#ifdef ZERO_HALF_CUTS
4295    CglZeroHalf zerohalfGen;
4296    //zerohalfGen.switchOnExpensive();
4297    // set default action (0=off,1=on,2=root)
4298    int zerohalfAction=0;
4299#endif
4300
4301    // Stored cuts
4302    //bool storedCuts = false;
4303
4304    int useCosts=0;
4305    // don't use input solution
4306    int useSolution=-1;
4307   
4308    // total number of commands read
4309    int numberGoodCommands=0;
4310    // Set false if user does anything advanced
4311    bool defaultSettings=true;
4312
4313    // Hidden stuff for barrier
4314    int choleskyType = 0;
4315    int gamma=0;
4316    int scaleBarrier=0;
4317    int doKKT=0;
4318    int crossover=2; // do crossover unless quadratic
4319    // For names
4320    int lengthName = 0;
4321    std::vector<std::string> rowNames;
4322    std::vector<std::string> columnNames;
4323    // Default strategy stuff
4324    {
4325      // try changing tolerance at root
4326#define MORE_CUTS
4327#ifdef MORE_CUTS
4328      gomoryGen.setAwayAtRoot(0.005);
4329      twomirGen.setAwayAtRoot(0.005);
4330      twomirGen.setAway(0.01);
4331#else
4332      gomoryGen.setAwayAtRoot(0.01);
4333      twomirGen.setAwayAtRoot(0.01);
4334      twomirGen.setAway(0.01);
4335#endif
4336      int iParam;
4337      iParam = whichParam(DIVEOPT,numberParameters_,parameters_);
4338      parameters_[iParam].setIntValue(3);
4339      iParam = whichParam(FPUMPITS,numberParameters_,parameters_);
4340      parameters_[iParam].setIntValue(30);
4341      iParam = whichParam(FPUMPTUNE,numberParameters_,parameters_);
4342      parameters_[iParam].setIntValue(1005043);
4343      initialPumpTune=1005043;
4344      iParam = whichParam(PROCESSTUNE,numberParameters_,parameters_);
4345      parameters_[iParam].setIntValue(6);
4346      tunePreProcess=6;
4347      iParam = whichParam(DIVINGC,numberParameters_,parameters_);
4348      parameters_[iParam].setCurrentOption("on");
4349      iParam = whichParam(RINS,numberParameters_,parameters_);
4350      parameters_[iParam].setCurrentOption("on");
4351      iParam = whichParam(PROBINGCUTS,numberParameters_,parameters_);
4352      parameters_[iParam].setCurrentOption("on");
4353      probingAction=1;
4354    }
4355    std::string field;
4356    if (!noPrinting_) {
4357      sprintf(generalPrint,"Coin Cbc and Clp Solver version %s, build %s",
4358              CBCVERSION,__DATE__);
4359      generalMessageHandler->message(CLP_GENERAL,generalMessages)
4360        << generalPrint
4361        <<CoinMessageEol;
4362      // Print command line
4363      if (argc>1) {
4364        bool foundStrategy=false;
4365        sprintf(generalPrint,"command line - ");
4366        for (int i=0;i<argc;i++) {
4367          if (!argv[i])
4368            break;
4369          if (strstr(argv[i],"strat"))
4370            foundStrategy=true;
4371          sprintf(generalPrint+strlen(generalPrint),"%s ",argv[i]);
4372        }
4373        if (!foundStrategy)
4374          sprintf(generalPrint+strlen(generalPrint),"(default strategy 1)");
4375        generalMessageHandler->message(CLP_GENERAL,generalMessages)
4376          << generalPrint
4377          <<CoinMessageEol;
4378      }
4379    }
4380    while (1) {
4381      // next command
4382      field=CoinReadGetCommand(argc,argv);
4383      // Reset time
4384      time1 = CoinCpuTime();
4385      // adjust field if has odd trailing characters
4386      char temp [200];
4387      strcpy(temp,field.c_str());
4388      int length = strlen(temp);
4389      for (int k=length-1;k>=0;k--) {
4390        if (temp[k]<' ')
4391          length--;
4392        else
4393          break;
4394      }
4395      temp[length]='\0';
4396      field=temp;
4397      // exit if null or similar
4398      if (!field.length()) {
4399        if (numberGoodCommands==1&&goodModel) {
4400          // we just had file name - do branch and bound
4401          field="branch";
4402        } else if (!numberGoodCommands) {
4403          // let's give the sucker a hint
4404          std::cout
4405            <<"CoinSolver takes input from arguments ( - switches to stdin)"
4406            <<std::endl
4407            <<"Enter ? for list of commands or help"<<std::endl;
4408          field="-";
4409        } else {
4410          break;
4411        }
4412      }
4413     
4414      // see if ? at end
4415      int numberQuery=0;
4416      if (field!="?"&&field!="???") {
4417        int length = field.length();
4418        int i;
4419        for (i=length-1;i>0;i--) {
4420          if (field[i]=='?') 
4421            numberQuery++;
4422          else
4423            break;
4424        }
4425        field=field.substr(0,length-numberQuery);
4426      }
4427      // find out if valid command
4428      int iParam;
4429      int numberMatches=0;
4430      int firstMatch=-1;
4431      for ( iParam=0; iParam<numberParameters_; iParam++ ) {
4432        int match = parameters_[iParam].matches(field);
4433        if (match==1) {
4434          numberMatches = 1;
4435          firstMatch=iParam;
4436          break;
4437        } else {
4438          if (match&&firstMatch<0)
4439            firstMatch=iParam;
4440          numberMatches += match>>1;
4441        }
4442      }
4443      if (iParam<numberParameters_&&!numberQuery) {
4444        // found
4445        CbcOrClpParam found = parameters_[iParam];
4446        CbcOrClpParameterType type = found.type();
4447        int valid;
4448        numberGoodCommands++;
4449        if (type==BAB&&goodModel) {
4450          // check if any integers
4451#ifndef CBC_OTHER_SOLVER
4452#ifdef COIN_HAS_ASL
4453          if (info.numberSos&&doSOS&&statusUserFunction_[0]) {
4454            // SOS
4455            numberSOS = info.numberSos;
4456          }
4457#endif
4458          lpSolver = clpSolver->getModelPtr();
4459          if (!lpSolver->integerInformation()&&!numberSOS&&
4460              !clpSolver->numberSOS()&&!model_.numberObjects()&&!clpSolver->numberObjects())
4461            type=DUALSIMPLEX;
4462#endif
4463        }
4464        if (type==GENERALQUERY) {
4465          bool evenHidden=false;
4466          int printLevel = 
4467            parameters_[whichParam(ALLCOMMANDS,
4468                                   numberParameters_,parameters_)].currentOptionAsInteger();
4469          int convertP[]={2,1,0};
4470          printLevel=convertP[printLevel];
4471          if ((verbose&8)!=0) {
4472            // even hidden
4473            evenHidden = true;
4474            verbose &= ~8;
4475          }
4476#ifdef COIN_HAS_ASL
4477          if (verbose<4&&statusUserFunction_[0])
4478            verbose +=4;
4479#endif
4480          if (verbose<4) {
4481            std::cout<<"In argument list keywords have leading - "
4482              ", -stdin or just - switches to stdin"<<std::endl;
4483            std::cout<<"One command per line (and no -)"<<std::endl;
4484            std::cout<<"abcd? gives list of possibilities, if only one + explanation"<<std::endl;
4485            std::cout<<"abcd?? adds explanation, if only one fuller help"<<std::endl;
4486            std::cout<<"abcd without value (where expected) gives current value"<<std::endl;
4487            std::cout<<"abcd value sets value"<<std::endl;
4488            std::cout<<"Commands are:"<<std::endl;
4489          } else {
4490            std::cout<<"Cbc options are set within AMPL with commands like:"<<std::endl<<std::endl;
4491            std::cout<<"         option cbc_options \"cuts=root log=2 feas=on slog=1\""<<std::endl<<std::endl;
4492            std::cout<<"only maximize, dual, primal, help and quit are recognized without ="<<std::endl;
4493          }
4494          int maxAcross=10;
4495          if ((verbose%4)!=0)
4496            maxAcross=1;
4497          int limits[]={1,51,101,151,201,251,301,351,401};
4498          std::vector<std::string> types;
4499          types.push_back("Double parameters:");
4500          types.push_back("Branch and Cut double parameters:");
4501          types.push_back("Integer parameters:");
4502          types.push_back("Branch and Cut integer parameters:");
4503          types.push_back("Keyword parameters:");
4504          types.push_back("Branch and Cut keyword parameters:");
4505          types.push_back("Actions or string parameters:");
4506          types.push_back("Branch and Cut actions:");
4507          int iType;
4508          for (iType=0;iType<8;iType++) {
4509            int across=0;
4510            int lengthLine=0;
4511            if ((verbose%4)!=0)
4512              std::cout<<std::endl;
4513            std::cout<<types[iType]<<std::endl;
4514            if ((verbose&2)!=0)
4515              std::cout<<std::endl;
4516            for ( iParam=0; iParam<numberParameters_; iParam++ ) {
4517              int type = parameters_[iParam].type();
4518              //printf("%d type %d limits %d %d display %d\n",iParam,
4519              //     type,limits[iType],limits[iType+1],parameters_[iParam].displayThis());
4520              if ((parameters_[iParam].displayThis()>=printLevel||evenHidden)&&
4521                  type>=limits[iType]
4522                  &&type<limits[iType+1]) {
4523                // but skip if not useful for ampl (and in ampl mode)
4524                if (verbose>=4&&(parameters_[iParam].whereUsed()&4)==0)
4525                  continue;
4526                if (!across) {
4527                  if ((verbose&2)!=0) 
4528                    std::cout<<"Command ";
4529                }
4530                int length = parameters_[iParam].lengthMatchName()+1;
4531                if (lengthLine+length>80) {
4532                  std::cout<<std::endl;
4533                  across=0;
4534                  lengthLine=0;
4535                }
4536                std::cout<<" "<<parameters_[iParam].matchName();
4537                lengthLine += length;
4538                across++;
4539                if (across==maxAcross) {
4540                  across=0;
4541                  if ((verbose%4)!=0) {
4542                    // put out description as well
4543                    if ((verbose&1)!=0) 
4544                      std::cout<<" "<<parameters_[iParam].shortHelp();
4545                    std::cout<<std::endl;
4546                    if ((verbose&2)!=0) {
4547                      std::cout<<"---- description"<<std::endl;
4548                      parameters_[iParam].printLongHelp();
4549                      std::cout<<"----"<<std::endl<<std::endl;
4550                    }
4551                  } else {
4552                    std::cout<<std::endl;
4553                  }
4554                }
4555              }
4556            }
4557            if (across)
4558              std::cout<<std::endl;
4559          }
4560        } else if (type==FULLGENERALQUERY) {
4561          std::cout<<"Full list of commands is:"<<std::endl;
4562          int maxAcross=5;
4563          int limits[]={1,51,101,151,201,251,301,351,401};
4564          std::vector<std::string> types;
4565          types.push_back("Double parameters:");
4566          types.push_back("Branch and Cut double parameters:");
4567          types.push_back("Integer parameters:");
4568          types.push_back("Branch and Cut integer parameters:");
4569          types.push_back("Keyword parameters:");
4570          types.push_back("Branch and Cut keyword parameters:");
4571          types.push_back("Actions or string parameters:");
4572          types.push_back("Branch and Cut actions:");
4573          int iType;
4574          for (iType=0;iType<8;iType++) {
4575            int across=0;
4576            std::cout<<types[iType]<<"  ";
4577            for ( iParam=0; iParam<numberParameters_; iParam++ ) {
4578              int type = parameters_[iParam].type();
4579              if (type>=limits[iType]
4580                  &&type<limits[iType+1]) {
4581                if (!across)
4582                  std::cout<<"  ";
4583                std::cout<<parameters_[iParam].matchName()<<"  ";
4584                across++;
4585                if (across==maxAcross) {
4586                  std::cout<<std::endl;
4587                  across=0;
4588                }
4589              }
4590            }
4591            if (across)
4592              std::cout<<std::endl;
4593          }
4594        } else if (type<101) {
4595          // get next field as double
4596          double value = CoinReadGetDoubleField(argc,argv,&valid);
4597          if (!valid) {
4598            if (type<51) {
4599              int returnCode;
4600              const char * message = 
4601                parameters_[iParam].setDoubleParameterWithMessage(lpSolver,value,returnCode);
4602              if (!noPrinting_&&strlen(message)) {
4603                generalMessageHandler->message(CLP_GENERAL,generalMessages)
4604                  << message
4605                  <<CoinMessageEol;
4606              }
4607            } else if (type<81) {
4608              int returnCode;
4609              const char * message = 
4610                parameters_[iParam].setDoubleParameterWithMessage(model_,value,returnCode);
4611              if (!noPrinting_&&strlen(message)) {
4612                generalMessageHandler->message(CLP_GENERAL,generalMessages)
4613                  << message
4614                  <<CoinMessageEol;
4615              }
4616            } else {
4617              int returnCode;
4618              const char * message = 
4619                parameters_[iParam].setDoubleParameterWithMessage(lpSolver,value,returnCode);
4620              if (!noPrinting_&&strlen(message)) {
4621                generalMessageHandler->message(CLP_GENERAL,generalMessages)
4622                  << message
4623                  <<CoinMessageEol;
4624              }
4625              switch(type) {
4626              case DJFIX:
4627                djFix=value;
4628#ifndef CBC_OTHER_SOLVER
4629                if (goodModel&&djFix<1.0e20) {
4630                  // do some fixing
4631                  clpSolver = dynamic_cast< OsiClpSolverInterface*> (model_.solver());
4632                  clpSolver->initialSolve();
4633                  lpSolver = clpSolver->getModelPtr();
4634                  int numberColumns = lpSolver->numberColumns();
4635                  int i;
4636                  const char * type = lpSolver->integerInformation();
4637                  double * lower = lpSolver->columnLower();
4638                  double * upper = lpSolver->columnUpper();
4639                  double * solution = lpSolver->primalColumnSolution();
4640                  double * dj = lpSolver->dualColumnSolution();
4641                  int numberFixed=0;
4642                  double dextra4 = parameters_[whichParam(DEXTRA4,numberParameters_,parameters_)].doubleValue();
4643                  if (dextra4)
4644                    printf("Multiple for continuous dj fixing is %g\n",dextra4);
4645                  for (i=0;i<numberColumns;i++) {
4646                    double djValue = dj[i];
4647                    if (!type[i])
4648                      djValue *= dextra4;
4649                    if (type[i]||dextra4) {
4650                      double value = solution[i];
4651                      if (value<lower[i]+1.0e-5&&djValue>djFix) {
4652                        solution[i]=lower[i];
4653                        upper[i]=lower[i];
4654                        numberFixed++;
4655                      } else if (value>upper[i]-1.0e-5&&djValue<-djFix) {
4656                        solution[i]=upper[i];
4657                        lower[i]=upper[i];
4658                        numberFixed++;
4659                      }
4660                    }
4661                  }
4662                  sprintf(generalPrint,"%d columns fixed\n",numberFixed);
4663                  generalMessageHandler->message(CLP_GENERAL,generalMessages)
4664                    << generalPrint
4665                    <<CoinMessageEol;
4666                }
4667#endif
4668                break;
4669              case TIGHTENFACTOR:
4670                tightenFactor=value;
4671                if(!complicatedInteger)
4672                  defaultSettings=false; // user knows what she is doing
4673                break;
4674              default:
4675                break;
4676              }
4677            }
4678          } else if (valid==1) {
4679            std::cout<<" is illegal for double parameter "<<parameters_[iParam].name()<<" value remains "<<
4680              parameters_[iParam].doubleValue()<<std::endl;
4681          } else {
4682            std::cout<<parameters_[iParam].name()<<" has value "<<
4683              parameters_[iParam].doubleValue()<<std::endl;
4684          }
4685        } else if (type<201) {
4686          // get next field as int
4687          int value = CoinReadGetIntField(argc,argv,&valid);
4688          if (!valid) {
4689            if (type<151) {
4690              if (parameters_[iParam].type()==PRESOLVEPASS)
4691                preSolve = value;
4692              else if (parameters_[iParam].type()==IDIOT)
4693                doIdiot = value;
4694              else if (parameters_[iParam].type()==SPRINT)
4695                doSprint = value;
4696              else if (parameters_[iParam].type()==OUTPUTFORMAT)
4697                outputFormat = value;
4698              else if (parameters_[iParam].type()==SLPVALUE)
4699                slpValue = value;
4700              else if (parameters_[iParam].type()==CPP)
4701                cppValue = value;
4702              else if (parameters_[iParam].type()==PRESOLVEOPTIONS)
4703                presolveOptions = value;
4704              else if (parameters_[iParam].type()==PRINTOPTIONS)
4705                printOptions = value;
4706              else if (parameters_[iParam].type()==SUBSTITUTION)
4707                substitution = value;
4708              else if (parameters_[iParam].type()==DUALIZE)
4709                dualize = value;
4710              else if (parameters_[iParam].type()==PROCESSTUNE)
4711                tunePreProcess = value;
4712              else if (parameters_[iParam].type()==USESOLUTION)
4713                useSolution = value;
4714              else if (parameters_[iParam].type()==VERBOSE)
4715                verbose = value;
4716              int returnCode;
4717              const char * message = 
4718                parameters_[iParam].setIntParameterWithMessage(lpSolver,value,returnCode);
4719              if (!noPrinting_&&strlen(message)) {
4720                generalMessageHandler->message(CLP_GENERAL,generalMessages)
4721                  << message
4722                  <<CoinMessageEol;
4723              }
4724            } else {
4725              if (parameters_[iParam].type()==CUTPASS)
4726                cutPass = value;
4727              else if (parameters_[iParam].type()==CUTPASSINTREE)
4728                cutPassInTree = value;
4729              else if (parameters_[iParam].type()==STRONGBRANCHING||
4730                       parameters_[iParam].type()==NUMBERBEFORE)
4731                strongChanged=true;
4732              else if (parameters_[iParam].type()==FPUMPTUNE||
4733                       parameters_[iParam].type()==FPUMPITS)
4734                pumpChanged=true;
4735              else if (parameters_[iParam].type()==EXPERIMENT) {
4736                if (value>=1) {
4737                  if (!noPrinting_) {
4738                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
4739                      <<"switching on dense factorization if small, and maybe fast fathoming"
4740                      <<CoinMessageEol;
4741                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
4742                      <<"Gomory cuts using tolerance of 0.01 at root"
4743                      <<CoinMessageEol;
4744                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
4745                      <<"Possible restart after 100 nodes if can fix many"
4746                      <<CoinMessageEol;
4747                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
4748                      <<"extra options - -diving C -diveopt 3 -rins on -tune 6 -probing on -passf 30!"
4749                      <<CoinMessageEol;
4750                  }
4751                  // try changing tolerance at root
4752                  //gomoryGen.setAwayAtRoot(0.01);
4753                  int iParam;
4754                  iParam = whichParam(DIVEOPT,numberParameters_,parameters_);
4755                  parameters_[iParam].setIntValue(3);
4756                  iParam = whichParam(FPUMPITS,numberParameters_,parameters_);
4757                  parameters_[iParam].setIntValue(30);
4758                  iParam = whichParam(FPUMPTUNE,numberParameters_,parameters_);
4759                  parameters_[iParam].setIntValue(1005043);
4760                  initialPumpTune=1005043;
4761                  iParam = whichParam(PROCESSTUNE,numberParameters_,parameters_);
4762                  parameters_[iParam].setIntValue(6);
4763                  tunePreProcess=6;
4764                  iParam = whichParam(DIVINGC,numberParameters_,parameters_);
4765                  parameters_[iParam].setCurrentOption("on");
4766                  iParam = whichParam(RINS,numberParameters_,parameters_);
4767                  parameters_[iParam].setCurrentOption("on");
4768                  iParam = whichParam(PROBINGCUTS,numberParameters_,parameters_);
4769                  parameters_[iParam].setCurrentOption("on");
4770                  probingAction=1;
4771                }
4772              } else if (parameters_[iParam].type()==STRATEGY) {
4773                if (value==0) {
4774                  gomoryGen.setAwayAtRoot(0.05);
4775                  int iParam;
4776                  iParam = whichParam(DIVEOPT,numberParameters_,parameters_);
4777                  parameters_[iParam].setIntValue(-1);
4778                  iParam = whichParam(FPUMPITS,numberParameters_,parameters_);
4779                  parameters_[iParam].setIntValue(20);
4780                  iParam = whichParam(FPUMPTUNE,numberParameters_,parameters_);
4781                  parameters_[iParam].setIntValue(1003);
4782                  initialPumpTune=1003;
4783                  iParam = whichParam(PROCESSTUNE,numberParameters_,parameters_);
4784                  parameters_[iParam].setIntValue(-1);
4785                  tunePreProcess=0;
4786                  iParam = whichParam(DIVINGC,numberParameters_,parameters_);
4787                  parameters_[iParam].setCurrentOption("off");
4788                  iParam = whichParam(RINS,numberParameters_,parameters_);
4789                  parameters_[iParam].setCurrentOption("off");
4790                  iParam = whichParam(PROBINGCUTS,numberParameters_,parameters_);
4791                  parameters_[iParam].setCurrentOption("on");
4792                  probingAction=1;
4793                }
4794              }
4795              int returnCode;
4796              const char * message = 
4797                parameters_[iParam].setIntParameterWithMessage(model_,value,returnCode);
4798              if (!noPrinting_&&strlen(message)) {
4799                generalMessageHandler->message(CLP_GENERAL,generalMessages)
4800                  << message
4801                  <<CoinMessageEol;
4802              }
4803            }
4804          } else if (valid==1) {
4805            std::cout<<" is illegal for integer parameter "<<parameters_[iParam].name()<<" value remains "<<
4806              parameters_[iParam].intValue()<<std::endl;
4807          } else {
4808            std::cout<<parameters_[iParam].name()<<" has value "<<
4809              parameters_[iParam].intValue()<<std::endl;
4810          }
4811        } else if (type<301) {
4812          // one of several strings
4813          std::string value = CoinReadGetString(argc,argv);
4814          int action = parameters_[iParam].parameterOption(value);
4815          if (action<0) {
4816            if (value!="EOL") {
4817              // no match
4818              parameters_[iParam].printOptions();
4819            } else {
4820              // print current value
4821              std::cout<<parameters_[iParam].name()<<" has value "<<
4822                parameters_[iParam].currentOption()<<std::endl;
4823            }
4824          } else {
4825            const char * message = 
4826              parameters_[iParam].setCurrentOptionWithMessage(action);
4827            if (!noPrinting_&&strlen(message)) {
4828              generalMessageHandler->message(CLP_GENERAL,generalMessages)
4829                << message
4830                <<CoinMessageEol;
4831            }
4832            // for now hard wired
4833            switch (type) {
4834            case DIRECTION:
4835              if (action==0)
4836                lpSolver->setOptimizationDirection(1);
4837              else if (action==1)
4838                lpSolver->setOptimizationDirection(-1);
4839              else
4840                lpSolver->setOptimizationDirection(0);
4841              break;
4842            case DUALPIVOT:
4843              if (action==0) {
4844                ClpDualRowSteepest steep(3);
4845                lpSolver->setDualRowPivotAlgorithm(steep);
4846              } else if (action==1) {
4847                ClpDualRowDantzig dantzig;
4848                //ClpDualRowSteepest dantzig(5);
4849                lpSolver->setDualRowPivotAlgorithm(dantzig);
4850              } else if (action==2) {
4851                // partial steep
4852                ClpDualRowSteepest steep(2);
4853                lpSolver->setDualRowPivotAlgorithm(steep);
4854              } else {
4855                ClpDualRowSteepest steep;
4856                lpSolver->setDualRowPivotAlgorithm(steep);
4857              }
4858              break;
4859            case PRIMALPIVOT:
4860              if (action==0) {
4861                ClpPrimalColumnSteepest steep(3);
4862                lpSolver->setPrimalColumnPivotAlgorithm(steep);
4863              } else if (action==1) {
4864                ClpPrimalColumnSteepest steep(0);
4865                lpSolver->setPrimalColumnPivotAlgorithm(steep);
4866              } else if (action==2) {
4867                ClpPrimalColumnDantzig dantzig;
4868                lpSolver->setPrimalColumnPivotAlgorithm(dantzig);
4869              } else if (action==3) {
4870                ClpPrimalColumnSteepest steep(4);
4871                lpSolver->setPrimalColumnPivotAlgorithm(steep);
4872              } else if (action==4) {
4873                ClpPrimalColumnSteepest steep(1);
4874                lpSolver->setPrimalColumnPivotAlgorithm(steep);
4875              } else if (action==5) {
4876                ClpPrimalColumnSteepest steep(2);
4877                lpSolver->setPrimalColumnPivotAlgorithm(steep);
4878              } else if (action==6) {
4879                ClpPrimalColumnSteepest steep(10);
4880                lpSolver->setPrimalColumnPivotAlgorithm(steep);
4881              }
4882              break;
4883            case SCALING:
4884              lpSolver->scaling(action);
4885              solver->setHintParam(OsiDoScale,action!=0,OsiHintTry);
4886              doScaling = action;
4887              break;
4888            case AUTOSCALE:
4889              lpSolver->setAutomaticScaling(action!=0);
4890              break;
4891            case SPARSEFACTOR:
4892              lpSolver->setSparseFactorization((1-action)!=0);
4893              break;
4894            case BIASLU:
4895              lpSolver->factorization()->setBiasLU(action);
4896              break;
4897            case PERTURBATION:
4898              if (action==0)
4899                lpSolver->setPerturbation(50);
4900              else
4901                lpSolver->setPerturbation(100);
4902              break;
4903            case ERRORSALLOWED:
4904              allowImportErrors = action;
4905              break;
4906            case INTPRINT:
4907              printMode=action;
4908              break;
4909              //case ALGORITHM:
4910              //algorithm  = action;
4911              //defaultSettings=false; // user knows what she is doing
4912              //abort();
4913              //break;
4914            case KEEPNAMES:
4915              keepImportNames = 1-action;
4916              break;
4917            case PRESOLVE:
4918              if (action==0)
4919                preSolve = 5;
4920              else if (action==1)
4921                preSolve=0;
4922              else if (action==2)
4923                preSolve=10;
4924              else
4925                preSolveFile=true;
4926              break;
4927            case PFI:
4928              lpSolver->factorization()->setForrestTomlin(action==0);
4929              break;
4930            case FACTORIZATION:
4931              lpSolver->factorization()->forceOtherFactorization(action);
4932              break;
4933            case CRASH:
4934              doCrash=action;
4935              break;
4936            case VECTOR:
4937              doVector=action;
4938              break;
4939            case MESSAGES:
4940              lpSolver->messageHandler()->setPrefix(action!=0);
4941              break;
4942            case CHOLESKY:
4943              choleskyType = action;
4944              break;
4945            case GAMMA:
4946              gamma=action;
4947              break;
4948            case BARRIERSCALE:
4949              scaleBarrier=action;
4950              break;
4951            case KKT:
4952              doKKT=action;
4953              break;
4954            case CROSSOVER:
4955              crossover=action;
4956              break;
4957            case SOS:
4958              doSOS=action;
4959              break;
4960            case GOMORYCUTS:
4961              defaultSettings=false; // user knows what she is doing
4962              gomoryAction = action;
4963              break;
4964            case PROBINGCUTS:
4965              defaultSettings=false; // user knows what she is doing
4966              probingAction = action;
4967              break;
4968            case KNAPSACKCUTS:
4969              defaultSettings=false; // user knows what she is doing
4970              knapsackAction = action;
4971              break;
4972            case REDSPLITCUTS:
4973              defaultSettings=false; // user knows what she is doing
4974              redsplitAction = action;
4975              break;
4976            case CLIQUECUTS:
4977              defaultSettings=false; // user knows what she is doing
4978              cliqueAction = action;
4979              break;
4980            case FLOWCUTS:
4981              defaultSettings=false; // user knows what she is doing
4982              flowAction = action;
4983              break;
4984            case MIXEDCUTS:
4985              defaultSettings=false; // user knows what she is doing
4986              mixedAction = action;
4987              break;
4988            case TWOMIRCUTS:
4989              defaultSettings=false; // user knows what she is doing
4990              twomirAction = action;
4991              break;
4992            case LANDPCUTS:
4993              defaultSettings=false; // user knows what she is doing
4994              landpAction = action;
4995              break;
4996            case RESIDCUTS:
4997              defaultSettings=false; // user knows what she is doing
4998              residualCapacityAction = action;
4999              break;
5000#ifdef ZERO_HALF_CUTS
5001            case ZEROHALFCUTS:
5002              defaultSettings=false; // user knows what she is doing
5003              zerohalfAction = action;
5004              break;
5005#endif
5006            case ROUNDING:
5007              defaultSettings=false; // user knows what she is doing
5008              break;
5009            case FPUMP:
5010              defaultSettings=false; // user knows what she is doing
5011              break;
5012            case RINS:
5013              break;
5014            case DINS:
5015              break;
5016            case RENS:
5017              break;
5018            case CUTSSTRATEGY:
5019              gomoryAction = action;
5020              probingAction = action;
5021              knapsackAction = action;
5022#ifdef ZERO_HALF_CUTS
5023              zerohalfAction = action;
5024#endif
5025              cliqueAction = action;
5026              flowAction = action;
5027              mixedAction = action;
5028              twomirAction = action;
5029              //landpAction = action;
5030              parameters_[whichParam(GOMORYCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5031              parameters_[whichParam(PROBINGCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5032              parameters_[whichParam(KNAPSACKCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5033              parameters_[whichParam(CLIQUECUTS,numberParameters_,parameters_)].setCurrentOption(action);
5034              parameters_[whichParam(FLOWCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5035              parameters_[whichParam(MIXEDCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5036              parameters_[whichParam(TWOMIRCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5037#ifdef ZERO_HALF_CUTS
5038              parameters_[whichParam(ZEROHALFCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5039#endif
5040              if (!action) {
5041                redsplitAction = action;
5042                parameters_[whichParam(REDSPLITCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5043                landpAction = action;
5044                parameters_[whichParam(LANDPCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5045                residualCapacityAction = action;
5046                parameters_[whichParam(RESIDCUTS,numberParameters_,parameters_)].setCurrentOption(action);
5047              }
5048              break;
5049            case HEURISTICSTRATEGY:
5050              parameters_[whichParam(ROUNDING,numberParameters_,parameters_)].setCurrentOption(action);
5051              parameters_[whichParam(GREEDY,numberParameters_,parameters_)].setCurrentOption(action);
5052              parameters_[whichParam(COMBINE,numberParameters_,parameters_)].setCurrentOption(action);
5053              //parameters_[whichParam(LOCALTREE,numberParameters_,parameters_)].setCurrentOption(action);
5054              parameters_[whichParam(FPUMP,numberParameters_,parameters_)].setCurrentOption(action);
5055              parameters_[whichParam(DIVINGC,numberParameters_,parameters_)].setCurrentOption(action);
5056              parameters_[whichParam(RINS,numberParameters_,parameters_)].setCurrentOption(action);
5057              break;
5058            case GREEDY:
5059            case DIVINGS:
5060            case DIVINGC:
5061            case DIVINGF:
5062            case DIVINGG:
5063            case DIVINGL:
5064            case DIVINGP:
5065            case DIVINGV:
5066            case COMBINE:
5067            case PIVOTANDFIX:
5068            case RANDROUND:
5069            case LOCALTREE:
5070            case NAIVE:
5071            case CPX:
5072              defaultSettings=false; // user knows what she is doing
5073              break;
5074            case COSTSTRATEGY:
5075              useCosts=action;
5076              break;
5077            case NODESTRATEGY:
5078              nodeStrategy=action;
5079              break;
5080            case PREPROCESS:
5081              preProcess = action;
5082              break;
5083            default:
5084              //abort();
5085              break;
5086            }
5087          }
5088        } else {
5089          // action
5090          if (type==EXIT) {
5091#ifdef COIN_HAS_ASL
5092            if(statusUserFunction_[0]) {
5093              if (info.numberIntegers||info.numberBinary) {
5094                // integer
5095              } else {
5096                // linear
5097              }
5098              writeAmpl(&info);
5099              freeArrays2(&info);
5100              freeArgs(&info);
5101            }
5102#endif
5103            break; // stop all
5104          }
5105          switch (type) {
5106          case DUALSIMPLEX:
5107          case PRIMALSIMPLEX:
5108          case SOLVECONTINUOUS:
5109          case BARRIER:
5110            if (goodModel) {
5111              // Say not in integer
5112              integerStatus=-1;
5113              double objScale = 
5114                parameters_[whichParam(OBJSCALE2,numberParameters_,parameters_)].doubleValue();
5115              if (objScale!=1.0) {
5116                int iColumn;
5117                int numberColumns=lpSolver->numberColumns();
5118                double * dualColumnSolution = 
5119                  lpSolver->dualColumnSolution();
5120                ClpObjective * obj = lpSolver->objectiveAsObject();
5121                assert(dynamic_cast<ClpLinearObjective *> (obj));
5122                double offset;
5123                double * objective = obj->gradient(NULL,NULL,offset,true);
5124                for (iColumn=0;iColumn<numberColumns;iColumn++) {
5125                  dualColumnSolution[iColumn] *= objScale;
5126                  objective[iColumn] *= objScale;;
5127                }
5128                int iRow;
5129                int numberRows=lpSolver->numberRows();
5130                double * dualRowSolution = 
5131                  lpSolver->dualRowSolution();
5132                for (iRow=0;iRow<numberRows;iRow++) 
5133                  dualRowSolution[iRow] *= objScale;
5134                lpSolver->setObjectiveOffset(objScale*lpSolver->objectiveOffset());
5135              }
5136              ClpSolve::SolveType method;
5137              ClpSolve::PresolveType presolveType;
5138              ClpSimplex * model2 = lpSolver;
5139              if (dualize) {
5140                bool tryIt=true;
5141                double fractionColumn=1.0;
5142                double fractionRow=1.0;
5143                if (dualize==3) {
5144                  dualize=1;
5145                  int numberColumns=lpSolver->numberColumns();
5146                  int numberRows=lpSolver->numberRows();
5147                  if (numberRows<50000||5*numberColumns>numberRows) {
5148                    tryIt=false;
5149                  } else {
5150                    fractionColumn=0.1;
5151                    fractionRow=0.1;
5152                  }
5153                }
5154                if (tryIt) {
5155                  model2 = static_cast<ClpSimplexOther *> (model2)->dualOfModel(fractionRow,fractionColumn);
5156                  if (model2) {
5157                    sprintf(generalPrint,"Dual of model has %d rows and %d columns",
5158                            model2->numberRows(),model2->numberColumns());
5159                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
5160                      << generalPrint
5161                      <<CoinMessageEol;
5162                    model2->setOptimizationDirection(1.0);
5163                  } else {
5164                    model2 = lpSolver;
5165                    dualize=0;
5166                  }
5167                } else {
5168                  dualize=0;
5169                }
5170              }
5171              if (noPrinting_)
5172                lpSolver->setLogLevel(0);
5173              ClpSolve solveOptions;
5174              solveOptions.setPresolveActions(presolveOptions);
5175              solveOptions.setSubstitution(substitution);
5176              if (preSolve!=5&&preSolve) {
5177                presolveType=ClpSolve::presolveNumber;
5178                if (preSolve<0) {
5179                  preSolve = - preSolve;
5180                  if (preSolve<=100) {
5181                    presolveType=ClpSolve::presolveNumber;
5182                    sprintf(generalPrint,"Doing %d presolve passes - picking up non-costed slacks",
5183                           preSolve);
5184                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
5185                      << generalPrint
5186                      <<CoinMessageEol;
5187                    solveOptions.setDoSingletonColumn(true);
5188                  } else {
5189                    preSolve -=100;
5190                    presolveType=ClpSolve::presolveNumberCost;
5191                    sprintf(generalPrint,"Doing %d presolve passes - picking up costed slacks",
5192                           preSolve);
5193                    generalMessageHandler->message(CLP_GENERAL,generalMessages)
5194                      << generalPrint
5195                      <<CoinMessageEol;
5196                  }
5197                } 
5198              } else if (preSolve) {
5199                presolveType=ClpSolve::presolveOn;
5200              } else {
5201                presolveType=ClpSolve::presolveOff;
5202              }
5203              solveOptions.setPresolveType(presolveType,preSolve);
5204              if (type==DUALSIMPLEX||type==SOLVECONTINUOUS) {
5205                method=ClpSolve::useDual;
5206              } else if (type==PRIMALSIMPLEX) {
5207                method=ClpSolve::usePrimalorSprint;
5208              } else {
5209                method = ClpSolve::useBarrier;
5210                if (crossover==1) {
5211                  method=ClpSolve::useBarrierNoCross;
5212                } else if (crossover==2) {
5213                  ClpObjective * obj = lpSolver->objectiveAsObject();
5214                  if (obj->type()>1) {
5215                    method=ClpSolve::useBarrierNoCross;
5216                    presolveType=ClpSolve::presolveOff;
5217                    solveOptions.setPresolveType(presolveType,preSolve);
5218                  } 
5219                }
5220              }
5221              solveOptions.setSolveType(method);
5222              if(preSolveFile)
5223                presolveOptions |= 0x40000000;
5224              solveOptions.setSpecialOption(4,presolveOptions);
5225              solveOptions.setSpecialOption(5,printOptions);
5226              if (doVector) {
5227                ClpMatrixBase * matrix = lpSolver->clpMatrix();
5228                if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
5229                  ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
5230                  clpMatrix->makeSpecialColumnCopy();
5231                }
5232              }
5233              if (method==ClpSolve::useDual) {
5234                // dual
5235                if (doCrash)
5236                  solveOptions.setSpecialOption(0,1,doCrash); // crash
5237                else if (doIdiot)
5238                  solveOptions.setSpecialOption(0,2,doIdiot); // possible idiot
5239              } else if (method==ClpSolve::usePrimalorSprint) {
5240                // primal
5241                // if slp turn everything off
5242                if (slpValue>0) {
5243                  doCrash=false;
5244                  doSprint=0;
5245                  doIdiot=-1;
5246                  solveOptions.setSpecialOption(1,10,slpValue); // slp
5247                  method=ClpSolve::usePrimal;
5248                }
5249                if (doCrash) {
5250                  solveOptions.setSpecialOption(1,1,doCrash); // crash
5251                } else if (doSprint>0) {
5252                  // sprint overrides idiot
5253                  solveOptions.setSpecialOption(1,3,doSprint); // sprint
5254                } else if (doIdiot>0) {
5255                  solveOptions.setSpecialOption(1,2,doIdiot); // idiot
5256                } else if (slpValue<=0) {
5257                  if (doIdiot==0) {
5258                    if (doSprint==0)
5259                      solveOptions.setSpecialOption(1,4); // all slack
5260                    else
5261                      solveOptions.setSpecialOption(1,9); // all slack or sprint
5262                  } else {
5263                    if (doSprint==0)
5264                      solveOptions.setSpecialOption(1,8); // all slack or idiot
5265                    else
5266                      solveOptions.setSpecialOption(1,7); // initiative
5267                  }
5268                }
5269                if (basisHasValues==-1)
5270                  solveOptions.setSpecialOption(1,11); // switch off values
5271              } else if (method==ClpSolve::useBarrier||method==ClpSolve::useBarrierNoCross) {
5272                int barrierOptions = choleskyType;
5273                if (scaleBarrier)
5274                  barrierOptions |= 8;
5275                if (doKKT)
5276                  barrierOptions |= 16;
5277                if (gamma)
5278                  barrierOptions |= 32*gamma;
5279                if (crossover==3) 
5280                  barrierOptions |= 256; // try presolve in crossover
5281                solveOptions.setSpecialOption(4,barrierOptions);
5282              }
5283              model2->setMaximumSeconds(model_.getMaximumSeconds());
5284#ifdef COIN_HAS_LINK
5285              OsiSolverInterface * coinSolver = model_.solver();
5286              OsiSolverLink * linkSolver = dynamic_cast< OsiSolverLink*> (coinSolver);
5287              if (!linkSolver) {
5288                model2->initialSolve(solveOptions);
5289              } else {
5290                // special solver
5291                int testOsiOptions = parameters_[whichParam(TESTOSI,numberParameters_,parameters_)].intValue();
5292                double * solution = NULL;
5293                if (testOsiOptions<10) {
5294                  solution = linkSolver->nonlinearSLP(slpValue>0 ? slpValue : 20 ,1.0e-5);
5295                } else if (testOsiOptions>=10) {
5296                  CoinModel coinModel = *linkSolver->coinModel();
5297                  ClpSimplex * tempModel = approximateSolution(coinModel,slpValue>0 ? slpValue : 50 ,1.0e-5,0);
5298                  assert (tempModel);
5299                  solution = CoinCopyOfArray(tempModel->primalColumnSolution(),coinModel.numberColumns());
5300                  model2->setObjectiveValue(tempModel->objectiveValue());
5301                  model2->setProblemStatus(tempModel->problemStatus());
5302                  model2->setSecondaryStatus(tempModel->secondaryStatus());
5303                  delete tempModel;
5304                }
5305                if (solution) {
5306                  memcpy(model2->primalColumnSolution(),solution,