Context Navigation

CbcSolver.cpp @ 822

Last change on this file since 822 was 822, checked in by forrest, 13 years ago
prepare to use cliques
File size: 328.0 KB

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