source: branches/sandbox/Cbc/src/CbcSolver.cpp @ 1390

Last change on this file since 1390 was 1390, checked in by lou, 10 years ago

Remove CPX_KEEP_RESULT. Use CBC_OTHER_SOLVER = 1.

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