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

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

adding stuff for heuristics

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