Context Navigation

CbcModel.cpp @ 2419

Last change on this file since 2419 was 2419, checked in by forrest, 2 years ago
more string parameters and some sos priorities
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 775.8 KB

Line
1	/* $Id: CbcModel.cpp 2419 2018-11-15 16:15:28Z forrest $ */
2	// Copyright (C) 2002, International Business Machines
3	// Corporation and others. All Rights Reserved.
4	// This code is licensed under the terms of the Eclipse Public License (EPL).
5
6	#if defined(_MSC_VER)
7	// Turn off compiler warning about long names
8	# pragma warning(disable:4786)
9	#endif
10
11	#include "CbcConfig.h"
12
13	#include <string>
14	//#define CBC_DEBUG 1
15	//#define CHECK_CUT_COUNTS
16	//#define CHECK_NODE
17	//#define CHECK_NODE_FULL
18	//#define NODE_LOG
19	//#define GLOBAL_CUTS_JUST_POINTERS
20	#ifdef CGL_DEBUG_GOMORY
21	extern int gomory_try;
22	#endif
23	#include <cassert>
24	#include <cmath>
25	#include <cfloat>
26	#ifdef COIN_HAS_CLP
27	// include Presolve from Clp
28	#include "ClpPresolve.hpp"
29	#include "OsiClpSolverInterface.hpp"
30	#include "ClpNode.hpp"
31	#include "ClpDualRowDantzig.hpp"
32	#include "ClpSimplexPrimal.hpp"
33	#endif
34
35	#include "CbcEventHandler.hpp"
36
37	#include "OsiSolverInterface.hpp"
38	#include "OsiAuxInfo.hpp"
39	#include "OsiSolverBranch.hpp"
40	#include "OsiChooseVariable.hpp"
41	#include "CoinWarmStartBasis.hpp"
42	#include "CoinPackedMatrix.hpp"
43	#include "CoinHelperFunctions.hpp"
44	#include "CbcBranchActual.hpp"
45	#include "CbcBranchDynamic.hpp"
46	#include "CbcHeuristic.hpp"
47	#include "CbcHeuristicFPump.hpp"
48	#include "CbcHeuristicRINS.hpp"
49	#include "CbcHeuristicDive.hpp"
50	#include "CbcModel.hpp"
51	#include "CbcTreeLocal.hpp"
52	#include "CbcStatistics.hpp"
53	#include "CbcStrategy.hpp"
54	#include "CbcMessage.hpp"
55	#include "OsiRowCut.hpp"
56	#include "OsiColCut.hpp"
57	#include "OsiRowCutDebugger.hpp"
58	#include "OsiCuts.hpp"
59	#include "CbcCountRowCut.hpp"
60	#include "CbcCutGenerator.hpp"
61	#include "CbcFeasibilityBase.hpp"
62	#include "CbcFathom.hpp"
63	#include "CbcFullNodeInfo.hpp"
64	#ifdef COIN_HAS_NTY
65	#include "CbcSymmetry.hpp"
66	#endif
67	// include Probing
68	#include "CglProbing.hpp"
69	#include "CglGomory.hpp"
70	#include "CglTwomir.hpp"
71	// include preprocessing
72	#include "CglPreProcess.hpp"
73	#include "CglDuplicateRow.hpp"
74	#include "CglStored.hpp"
75	#include "CglClique.hpp"
76	#include "CglKnapsackCover.hpp"
77
78	#include "CoinTime.hpp"
79	#include "CoinMpsIO.hpp"
80
81	#include "CbcCompareActual.hpp"
82	#include "CbcTree.hpp"
83	// This may be dummy
84	#include "CbcThread.hpp"
85	/* Various functions local to CbcModel.cpp */
86
87	typedef struct {
88	double useCutoff;
89	CbcModel * model;
90	int switches;
91	} rootBundle;
92	static void * doRootCbcThread(void * voidInfo);
93
94	namespace {
95
96	//-------------------------------------------------------------------
97	// Returns the greatest common denominator of two
98	// positive integers, a and b, found using Euclid's algorithm
99	//-------------------------------------------------------------------
100	static int gcd(int a, int b)
101	{
102	int remainder = -1;
103	// make sure a<=b (will always remain so)
104	if (a > b) {
105	// Swap a and b
106	int temp = a;
107	a = b;
108	b = temp;
109	}
110	// if zero then gcd is nonzero (zero may occur in rhs of packed)
111	if (!a) {
112	if (b) {
113	return b;
114	} else {
115	printf("**** gcd given two zeros!!\n");
116	abort();
117	}
118	}
119	while (remainder) {
120	remainder = b % a;
121	b = a;
122	a = remainder;
123	}
124	return b;
125	}
126
127
128
129	#ifdef CHECK_NODE_FULL
130
131	/*
132	Routine to verify that tree linkage is correct. The invariant that is tested
133	is
134
135	reference count = (number of actual references) + (number of branches left)
136
137	The routine builds a set of paired arrays, info and count, by traversing the
138	tree. Each CbcNodeInfo is recorded in info, and the number of times it is
139	referenced (via the parent field) is recorded in count. Then a final check is
140	made to see if the numberPointingToThis_ field agrees.
141	*/
142
143	void verifyTreeNodes (const CbcTree * branchingTree, const CbcModel &model)
144
145	{
146	if (model.getNodeCount() == 661) return;
147	printf("*** CHECKING tree after %d nodes\n", model.getNodeCount()) ;
148
149	int j ;
150	int nNodes = branchingTree->size() ;
151	# define MAXINFO 1000
152	int *count = new int [MAXINFO] ;
153	CbcNodeInfo *info = new CbcNodeInfo[MAXINFO] ;
154	int nInfo = 0 ;
155	/*
156	Collect all CbcNodeInfo objects in info, by starting from each live node and
157	traversing back to the root. Nodes in the live set should have unexplored
158	branches remaining.
159
160	TODO: The `while (nodeInfo)' loop could be made to break on reaching a
161	common ancester (nodeInfo is found in info[k]). Alternatively, the
162	check could change to signal an error if nodeInfo is not found above a
163	common ancestor.
164	*/
165	for (j = 0 ; j < nNodes ; j++) {
166	CbcNode *node = branchingTree->nodePointer(j) ;
167	if (!node)
168	continue;
169	CbcNodeInfo *nodeInfo = node->nodeInfo() ;
170	int change = node->nodeInfo()->numberBranchesLeft() ;
171	assert(change) ;
172	while (nodeInfo) {
173	int k ;
174	for (k = 0 ; k < nInfo ; k++) {
175	if (nodeInfo == info[k]) break ;
176	}
177	if (k == nInfo) {
178	assert(nInfo < MAXINFO) ;
179	nInfo++ ;
180	info[k] = nodeInfo ;
181	count[k] = 0 ;
182	}
183	nodeInfo = nodeInfo->parent() ;
184	}
185	}
186	/*
187	Walk the info array. For each nodeInfo, look up its parent in info and
188	increment the corresponding count.
189	*/
190	for (j = 0 ; j < nInfo ; j++) {
191	CbcNodeInfo *nodeInfo = info[j] ;
192	nodeInfo = nodeInfo->parent() ;
193	if (nodeInfo) {
194	int k ;
195	for (k = 0 ; k < nInfo ; k++) {
196	if (nodeInfo == info[k]) break ;
197	}
198	assert (k < nInfo) ;
199	count[k]++ ;
200	}
201	}
202	/*
203	Walk the info array one more time and check that the invariant holds. The
204	number of references (numberPointingToThis()) should equal the sum of the
205	number of actual references (held in count[]) plus the number of potential
206	references (unexplored branches, numberBranchesLeft()).
207	*/
208	for (j = 0; j < nInfo; j++) {
209	CbcNodeInfo * nodeInfo = info[j] ;
210	if (nodeInfo) {
211	int k ;
212	for (k = 0; k < nInfo; k++)
213	if (nodeInfo == info[k])
214	break ;
215	printf("Nodeinfo %x - %d left, %d count\n",
216	nodeInfo,
217	nodeInfo->numberBranchesLeft(),
218	nodeInfo->numberPointingToThis()) ;
219	assert(nodeInfo->numberPointingToThis() ==
220	count[k] + nodeInfo->numberBranchesLeft()) ;
221	}
222	}
223
224	delete [] count ;
225	delete [] info ;
226
227	return ;
228	}
229
230	#endif /* CHECK_NODE_FULL */
231
232
233
234	#ifdef CHECK_CUT_COUNTS
235
236	/*
237	Routine to verify that cut reference counts are correct.
238	*/
239	void verifyCutCounts (const CbcTree * branchingTree, CbcModel &model)
240
241	{
242	printf("*** CHECKING cuts after %d nodes\n", model.getNodeCount()) ;
243
244	int j ;
245	int nNodes = branchingTree->size() ;
246
247	/*
248	cut.tempNumber_ exists for the purpose of doing this verification. Clear it
249	in all cuts. We traverse the tree by starting from each live node and working
250	back to the root. At each CbcNodeInfo, check for cuts.
251	*/
252	for (j = 0 ; j < nNodes ; j++) {
253	CbcNode *node = branchingTree->nodePointer(j) ;
254	CbcNodeInfo * nodeInfo = node->nodeInfo() ;
255	assert (node->nodeInfo()->numberBranchesLeft()) ;
256	while (nodeInfo) {
257	int k ;
258	for (k = 0 ; k < nodeInfo->numberCuts() ; k++) {
259	CbcCountRowCut *cut = nodeInfo->cuts()[k] ;
260	if (cut) cut->tempNumber_ = 0;
261	}
262	nodeInfo = nodeInfo->parent() ;
263	}
264	}
265	/*
266	Walk the live set again, this time collecting the list of cuts in use at each
267	node. addCuts1 will collect the cuts in model.addedCuts_. Take into account
268	that when we recreate the basis for a node, we compress out the slack cuts.
269	*/
270	for (j = 0 ; j < nNodes ; j++) {
271	CoinWarmStartBasis *debugws = model.getEmptyBasis() ;
272	CbcNode *node = branchingTree->nodePointer(j) ;
273	CbcNodeInfo *nodeInfo = node->nodeInfo();
274	int change = node->nodeInfo()->numberBranchesLeft() ;
275	printf("Node %d %x (info %x) var %d way %d obj %g", j, node,
276	node->nodeInfo(), node->columnNumber(), node->way(),
277	node->objectiveValue()) ;
278
279	model.addCuts1(node, debugws) ;
280
281	int i ;
282	int numberRowsAtContinuous = model.numberRowsAtContinuous() ;
283	CbcCountRowCut **addedCuts = model.addedCuts() ;
284	for (i = 0 ; i < model.currentNumberCuts() ; i++) {
285	CoinWarmStartBasis::Status status =
286	debugws->getArtifStatus(i + numberRowsAtContinuous) ;
287	if (status != CoinWarmStartBasis::basic && addedCuts[i]) {
288	addedCuts[i]->tempNumber_ += change ;
289	}
290	}
291
292	while (nodeInfo) {
293	nodeInfo = nodeInfo->parent() ;
294	if (nodeInfo) printf(" -> %x", nodeInfo);
295	}
296	printf("\n") ;
297	delete debugws ;
298	}
299	/*
300	The moment of truth: We've tallied up the references by direct scan of the search tree. Check for agreement with the count in the cut.
301
302	TODO: Rewrite to check and print mismatch only when tempNumber_ == 0?
303	*/
304	for (j = 0 ; j < nNodes ; j++) {
305	CbcNode *node = branchingTree->nodePointer(j) ;
306	CbcNodeInfo *nodeInfo = node->nodeInfo();
307	while (nodeInfo) {
308	int k ;
309	for (k = 0 ; k < nodeInfo->numberCuts() ; k++) {
310	CbcCountRowCut *cut = nodeInfo->cuts()[k] ;
311	if (cut && cut->tempNumber_ >= 0) {
312	if (cut->tempNumber_ != cut->numberPointingToThis())
313	printf("mismatch %x %d %x %d %d\n", nodeInfo, k,
314	cut, cut->tempNumber_, cut->numberPointingToThis()) ;
315	else
316	printf(" match %x %d %x %d %d\n", nodeInfo, k,
317	cut, cut->tempNumber_, cut->numberPointingToThis()) ;
318	cut->tempNumber_ = -1 ;
319	}
320	}
321	nodeInfo = nodeInfo->parent() ;
322	}
323	}
324
325	return ;
326	}
327
328	#endif /* CHECK_CUT_COUNTS */
329
330
331	#ifdef CHECK_CUT_SIZE
332
333	/*
334	Routine to verify that cut reference counts are correct.
335	*/
336	void verifyCutSize (const CbcTree * branchingTree, CbcModel &model)
337	{
338
339	int j ;
340	int nNodes = branchingTree->size() ;
341	int totalCuts = 0;
342
343	/*
344	cut.tempNumber_ exists for the purpose of doing this verification. Clear it
345	in all cuts. We traverse the tree by starting from each live node and working
346	back to the root. At each CbcNodeInfo, check for cuts.
347	*/
348	for (j = 0 ; j < nNodes ; j++) {
349	CbcNode *node = branchingTree->nodePointer(j) ;
350	CbcNodeInfo * nodeInfo = node->nodeInfo() ;
351	assert (node->nodeInfo()->numberBranchesLeft()) ;
352	while (nodeInfo) {
353	totalCuts += nodeInfo->numberCuts();
354	nodeInfo = nodeInfo->parent() ;
355	}
356	}
357	printf("*** CHECKING cuts (size) after %d nodes - %d cuts\n", model.getNodeCount(), totalCuts) ;
358	return ;
359	}
360
361	#endif /* CHECK_CUT_SIZE */
362
363	}
364
365	/* End unnamed namespace for CbcModel.cpp */
366
367
368	void
369	CbcModel::analyzeObjective ()
370	/*
371	Try to find a minimum change in the objective function. The first scan
372	checks that there are no continuous variables with non-zero coefficients,
373	and grabs the largest objective coefficient associated with an unfixed
374	integer variable. The second scan attempts to scale up the objective
375	coefficients to a point where they are sufficiently close to integer that
376	we can pretend they are integer, and calculate a gcd over the coefficients
377	of interest. This will be the minimum increment for the scaled coefficients.
378	The final action is to scale the increment back for the original coefficients
379	and install it, if it's better than the existing value.
380
381	John's note: We could do better than this.
382
383	John's second note - apologies for changing s to z
384	*/
385	{
386	const double *objective = getObjCoefficients() ;
387	const double *lower = getColLower() ;
388	const double *upper = getColUpper() ;
389	/*
390	Scan continuous and integer variables to see if continuous
391	are cover or network with integral rhs.
392	*/
393	double continuousMultiplier = 1.0;
394	double * coeffMultiplier = NULL;
395	double largestObj = 0.0;
396	double smallestObj = COIN_DBL_MAX;
397	{
398	const double *rowLower = getRowLower() ;
399	const double *rowUpper = getRowUpper() ;
400	int numberRows = solver_->getNumRows() ;
401	double * rhs = new double [numberRows];
402	memset(rhs, 0, numberRows*sizeof(double));
403	int iColumn;
404	int numberColumns = solver_->getNumCols() ;
405	// Column copy of matrix
406	int problemType = -1;
407	const double * element = solver_->getMatrixByCol()->getElements();
408	const int * row = solver_->getMatrixByCol()->getIndices();
409	const CoinBigIndex * columnStart = solver_->getMatrixByCol()->getVectorStarts();
410	const int * columnLength = solver_->getMatrixByCol()->getVectorLengths();
411	int numberInteger = 0;
412	int numberIntegerObj = 0;
413	int numberGeneralIntegerObj = 0;
414	int numberIntegerWeight = 0;
415	int numberContinuousObj = 0;
416	double cost = COIN_DBL_MAX;
417	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
418	if (upper[iColumn] == lower[iColumn]) {
419	CoinBigIndex start = columnStart[iColumn];
420	CoinBigIndex end = start + columnLength[iColumn];
421	for (CoinBigIndex j = start; j < end; j++) {
422	int iRow = row[j];
423	rhs[iRow] += lower[iColumn] * element[j];
424	}
425	} else {
426	double objValue = objective[iColumn];
427	if (solver_->isInteger(iColumn))
428	numberInteger++;
429	if (objValue) {
430	if (!solver_->isInteger(iColumn)) {
431	numberContinuousObj++;
432	} else {
433	largestObj = CoinMax(largestObj, fabs(objValue));
434	smallestObj = CoinMin(smallestObj, fabs(objValue));
435	numberIntegerObj++;
436	if (cost == COIN_DBL_MAX)
437	cost = objValue;
438	else if (cost != objValue)
439	cost = -COIN_DBL_MAX;
440	int gap = static_cast<int> (upper[iColumn] - lower[iColumn]);
441	if (gap > 1) {
442	numberGeneralIntegerObj++;
443	numberIntegerWeight += gap;
444	}
445	}
446	}
447	}
448	}
449	int iType = 0;
450	if (!numberContinuousObj && numberIntegerObj <= 5 && numberIntegerWeight <= 100 &&
451	numberIntegerObj*3 < numberObjects_ && !parentModel_ && solver_->getNumRows() > 100)
452	iType = 1 + 4 + (((moreSpecialOptions_&536870912)==0) ? 2 : 0);
453	else if (!numberContinuousObj && numberIntegerObj <= 100 &&
454	numberIntegerObj*5 < numberObjects_ && numberIntegerWeight <= 100 &&
455	!parentModel_ &&
456	solver_->getNumRows() > 100 && cost != -COIN_DBL_MAX)
457	iType = 4 + (((moreSpecialOptions_&536870912)==0) ? 2 : 0);
458	else if (!numberContinuousObj && numberIntegerObj <= 100 &&
459	numberIntegerObj*5 < numberObjects_ &&
460	!parentModel_ &&
461	solver_->getNumRows() > 100 && cost != -COIN_DBL_MAX)
462	iType = 8;
463	int iTest = getMaximumNodes();
464	if (iTest >= 987654320 && iTest < 987654330 && numberObjects_ && !parentModel_) {
465	iType = iTest - 987654320;
466	printf("Testing %d integer variables out of %d objects (%d integer) have cost of %g - %d continuous\n",
467	numberIntegerObj, numberObjects_, numberInteger, cost, numberContinuousObj);
468	if (iType == 9)
469	exit(77);
470	if (numberContinuousObj)
471	iType = 0;
472	}
473
474	//if (!numberContinuousObj&&(numberIntegerObj<=5\|\|cost!=-COIN_DBL_MAX)&&
475	//numberIntegerObj*3<numberObjects_&&!parentModel_&&solver_->getNumRows()>100) {
476	if (iType) {
477	/*
478	A) put high priority on (if none)
479	B) create artificial objective (if clp)
480	*/
481	int iPriority = -1;
482	for (int i = 0; i < numberObjects_; i++) {
483	int k = object_[i]->priority();
484	if (iPriority == -1)
485	iPriority = k;
486	else if (iPriority != k)
487	iPriority = -2;
488	}
489	bool branchOnSatisfied = ((iType & 1) != 0);
490	bool createFake = ((iType & 2) != 0);
491	bool randomCost = ((iType & 4) != 0);
492	if (iPriority >= 0) {
493	char general[200];
494	if (cost == -COIN_DBL_MAX) {
495	sprintf(general, "%d integer variables out of %d objects (%d integer) have costs - high priority",
496	numberIntegerObj, numberObjects_, numberInteger);
497	} else if (cost == COIN_DBL_MAX) {
498	sprintf(general, "No integer variables out of %d objects (%d integer) have costs",
499	numberObjects_, numberInteger);
500	branchOnSatisfied = false;
501	} else {
502	sprintf(general, "%d integer variables out of %d objects (%d integer) have cost of %g - high priority",
503	numberIntegerObj, numberObjects_, numberInteger, cost);
504	}
505	messageHandler()->message(CBC_GENERAL,
506	messages())
507	<< general << CoinMessageEol ;
508	sprintf(general, "branch on satisfied %c create fake objective %c random cost %c",
509	branchOnSatisfied ? 'Y' : 'N',
510	createFake ? 'Y' : 'N',
511	randomCost ? 'Y' : 'N');
512	messageHandler()->message(CBC_GENERAL,
513	messages())
514	<< general << CoinMessageEol ;
515	// switch off clp type branching
516	// no ? fastNodeDepth_ = -1;
517	int highPriority = (branchOnSatisfied) ? -999 : 100;
518	for (int i = 0; i < numberObjects_; i++) {
519	CbcSimpleInteger * thisOne = dynamic_cast <CbcSimpleInteger *> (object_[i]);
520	object_[i]->setPriority(1000);
521	if (thisOne) {
522	int iColumn = thisOne->columnNumber();
523	if (objective[iColumn])
524	thisOne->setPriority(highPriority);
525	}
526	}
527	}
528	#ifdef COIN_HAS_CLP
529	OsiClpSolverInterface * clpSolver
530	= dynamic_cast<OsiClpSolverInterface *> (solver_);
531	if (clpSolver && createFake) {
532	// Create artificial objective to be used when all else fixed
533	int numberColumns = clpSolver->getNumCols();
534	double * fakeObj = new double [numberColumns];
535	// Column copy
536	const CoinPackedMatrix * matrixByCol = clpSolver->getMatrixByCol();
537	//const double * element = matrixByCol.getElements();
538	//const int * row = matrixByCol.getIndices();
539	//const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
540	const int * columnLength = matrixByCol->getVectorLengths();
541	const double * solution = clpSolver->getColSolution();
542	#ifdef JJF_ZERO
543	int nAtBound = 0;
544	for (int i = 0; i < numberColumns; i++) {
545	double lowerValue = lower[i];
546	double upperValue = upper[i];
547	if (clpSolver->isInteger(i)) {
548	double lowerValue = lower[i];
549	double upperValue = upper[i];
550	double value = solution[i];
551	if (value < lowerValue + 1.0e-6 \|\|
552	value > upperValue - 1.0e-6)
553	nAtBound++;
554	}
555	}
556	#endif
557	/*
558	Generate a random objective function for problems where the given objective
559	function is not terribly useful. (Nearly feasible, single integer variable,
560	that sort of thing.
561	*/
562	CoinDrand48(true, 1234567);
563	for (int i = 0; i < numberColumns; i++) {
564	double lowerValue = lower[i];
565	double upperValue = upper[i];
566	double value = (randomCost) ? ceil((CoinDrand48() + 0.5) * 1000)
567	: i + 1 + columnLength[i] * 1000;
568	value *= 0.001;
569	//value += columnLength[i];
570	if (lowerValue > -1.0e5 \|\| upperValue < 1.0e5) {
571	if (fabs(lowerValue) > fabs(upperValue))
572	value = - value;
573	if (clpSolver->isInteger(i)) {
574	double solValue = solution[i];
575	// Better to add in 0.5 or 1.0??
576	if (solValue < lowerValue + 1.0e-6)
577	value = fabs(value) + 0.5; //fabs(value*1.5);
578	else if (solValue > upperValue - 1.0e-6)
579	value = -fabs(value) - 0.5; //-fabs(value*1.5);
580	}
581	} else {
582	value = 0.0;
583	}
584	fakeObj[i] = value;
585	}
586	// pass to solver
587	clpSolver->setFakeObjective(fakeObj);
588	delete [] fakeObj;
589	}
590	#endif
591	} else if (largestObj < smallestObj*5.0 && !parentModel_ &&
592	!numberContinuousObj &&
593	!numberGeneralIntegerObj &&
594	numberIntegerObj*2 < CoinMin(numberColumns,20)) {
595	// up priorities on costed
596	int iPriority = -1;
597	for (int i = 0; i < numberObjects_; i++) {
598	int k = object_[i]->priority();
599	if (iPriority == -1)
600	iPriority = k;
601	else if (iPriority != k)
602	iPriority = -2;
603	}
604	if (iPriority >= 100) {
605	#if CBC_USEFUL_PRINTING>1
606	printf("Setting variables with obj to high priority\n");
607	#endif
608	for (int i = 0; i < numberObjects_; i++) {
609	CbcSimpleInteger * obj =
610	dynamic_cast <CbcSimpleInteger *>(object_[i]) ;
611	if (obj) {
612	int iColumn = obj->columnNumber();
613	if (objective[iColumn])
614	object_[i]->setPriority(iPriority - 1);
615	}
616	}
617	}
618	}
619	int iRow;
620	for (iRow = 0; iRow < numberRows; iRow++) {
621	if (rowLower[iRow] > -1.0e20 &&
622	fabs(rowLower[iRow] - rhs[iRow] - floor(rowLower[iRow] - rhs[iRow] + 0.5)) > 1.0e-10) {
623	continuousMultiplier = 0.0;
624	break;
625	}
626	if (rowUpper[iRow] < 1.0e20 &&
627	fabs(rowUpper[iRow] - rhs[iRow] - floor(rowUpper[iRow] - rhs[iRow] + 0.5)) > 1.0e-10) {
628	continuousMultiplier = 0.0;
629	break;
630	}
631	// set rhs to limiting value
632	if (rowLower[iRow] != rowUpper[iRow]) {
633	if (rowLower[iRow] > -1.0e20) {
634	if (rowUpper[iRow] < 1.0e20) {
635	// no good
636	continuousMultiplier = 0.0;
637	break;
638	} else {
639	rhs[iRow] = rowLower[iRow] - rhs[iRow];
640	if (problemType < 0)
641	problemType = 3; // set cover
642	else if (problemType != 3)
643	problemType = 4;
644	}
645	} else {
646	rhs[iRow] = rowUpper[iRow] - rhs[iRow];
647	if (problemType < 0)
648	problemType = 1; // set partitioning <=
649	else if (problemType != 1)
650	problemType = 4;
651	}
652	} else {
653	rhs[iRow] = rowUpper[iRow] - rhs[iRow];
654	if (problemType < 0)
655	problemType = 3; // set partitioning ==
656	else if (problemType != 2)
657	problemType = 2;
658	}
659	if (fabs(rhs[iRow] - 1.0) > 1.0e-12)
660	problemType = 4;
661	}
662	if (continuousMultiplier) {
663	// 1 network, 2 cover, 4 negative cover
664	int possible = 7;
665	bool unitRhs = true;
666	// See which rows could be set cover
667	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
668	if (upper[iColumn] > lower[iColumn] + 1.0e-8) {
669	CoinBigIndex start = columnStart[iColumn];
670	CoinBigIndex end = start + columnLength[iColumn];
671	for (CoinBigIndex j = start; j < end; j++) {
672	double value = element[j];
673	if (value == 1.0) {
674	} else if (value == -1.0) {
675	rhs[row[j]] = -0.5;
676	} else {
677	rhs[row[j]] = -COIN_DBL_MAX;
678	}
679	}
680	}
681	}
682	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
683	if (upper[iColumn] > lower[iColumn] + 1.0e-8) {
684	if (!isInteger(iColumn)) {
685	CoinBigIndex start = columnStart[iColumn];
686	CoinBigIndex end = start + columnLength[iColumn];
687	double rhsValue = 0.0;
688	// 1 all ones, -1 all -1s, 2 all +- 1, 3 no good
689	int type = 0;
690	for (CoinBigIndex j = start; j < end; j++) {
691	double value = element[j];
692	if (fabs(value) != 1.0) {
693	type = 3;
694	break;
695	} else if (value == 1.0) {
696	if (!type)
697	type = 1;
698	else if (type != 1)
699	type = 2;
700	} else {
701	if (!type)
702	type = -1;
703	else if (type != -1)
704	type = 2;
705	}
706	int iRow = row[j];
707	if (rhs[iRow] == -COIN_DBL_MAX) {
708	type = 3;
709	break;
710	} else if (rhs[iRow] == -0.5) {
711	// different values
712	unitRhs = false;
713	} else if (rhsValue) {
714	if (rhsValue != rhs[iRow])
715	unitRhs = false;
716	} else {
717	rhsValue = rhs[iRow];
718	}
719	}
720	// if no elements OK
721	if (type == 3) {
722	// no good
723	possible = 0;
724	break;
725	} else if (type == 2) {
726	if (end - start > 2) {
727	// no good
728	possible = 0;
729	break;
730	} else {
731	// only network
732	possible &= 1;
733	if (!possible)
734	break;
735	}
736	} else if (type == 1) {
737	// only cover
738	possible &= 2;
739	if (!possible)
740	break;
741	} else if (type == -1) {
742	// only negative cover
743	possible &= 4;
744	if (!possible)
745	break;
746	}
747	}
748	}
749	}
750	if ((possible == 2 \|\| possible == 4) && !unitRhs) {
751	#if COIN_DEVELOP>1
752	printf("XXXXXX Continuous all +1 but different rhs\n");
753	#endif
754	possible = 0;
755	}
756	// may be all integer
757	if (possible != 7) {
758	if (!possible)
759	continuousMultiplier = 0.0;
760	else if (possible == 1)
761	continuousMultiplier = 1.0;
762	else
763	continuousMultiplier = 0.0; // 0.5 was incorrect;
764	#if COIN_DEVELOP>1
765	if (continuousMultiplier)
766	printf("XXXXXX multiplier of %g\n", continuousMultiplier);
767	#endif
768	if (continuousMultiplier == 0.5) {
769	coeffMultiplier = new double [numberColumns];
770	bool allOne = true;
771	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
772	coeffMultiplier[iColumn] = 1.0;
773	if (upper[iColumn] > lower[iColumn] + 1.0e-8) {
774	if (!isInteger(iColumn)) {
775	CoinBigIndex start = columnStart[iColumn];
776	int iRow = row[start];
777	double value = rhs[iRow];
778	assert (value >= 0.0);
779	if (value != 0.0 && value != 1.0)
780	allOne = false;
781	coeffMultiplier[iColumn] = 0.5 * value;
782	}
783	}
784	}
785	if (allOne) {
786	// back to old way
787	delete [] coeffMultiplier;
788	coeffMultiplier = NULL;
789	}
790	}
791	} else {
792	// all integer
793	problemType_ = problemType;
794	#if COIN_DEVELOP>1
795	printf("Problem type is %d\n", problemType_);
796	#endif
797	}
798	}
799
800	// But try again
801	if (continuousMultiplier < 1.0) {
802	memset(rhs, 0, numberRows*sizeof(double));
803	int * count = new int [numberRows];
804	memset(count, 0, numberRows*sizeof(int));
805	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
806	CoinBigIndex start = columnStart[iColumn];
807	CoinBigIndex end = start + columnLength[iColumn];
808	if (upper[iColumn] == lower[iColumn]) {
809	for (CoinBigIndex j = start; j < end; j++) {
810	int iRow = row[j];
811	rhs[iRow] += lower[iColumn] * element[j];
812	}
813	} else if (solver_->isInteger(iColumn)) {
814	for (CoinBigIndex j = start; j < end; j++) {
815	int iRow = row[j];
816	if (fabs(element[j] - floor(element[j] + 0.5)) > 1.0e-10)
817	rhs[iRow] = COIN_DBL_MAX;
818	}
819	} else {
820	for (CoinBigIndex j = start; j < end; j++) {
821	int iRow = row[j];
822	count[iRow]++;
823	if (fabs(element[j]) != 1.0)
824	rhs[iRow] = COIN_DBL_MAX;
825	}
826	}
827	}
828	// now look at continuous
829	bool allGood = true;
830	double direction = solver_->getObjSense() ;
831	int numberObj = 0;
832	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
833	if (upper[iColumn] > lower[iColumn]) {
834	double objValue = objective[iColumn] * direction;
835	if (objValue && !solver_->isInteger(iColumn)) {
836	numberObj++;
837	CoinBigIndex start = columnStart[iColumn];
838	CoinBigIndex end = start + columnLength[iColumn];
839	if (objValue > 0.0) {
840	// wants to be as low as possible
841	if (lower[iColumn] < -1.0e10 \|\| fabs(lower[iColumn] - floor(lower[iColumn] + 0.5)) > 1.0e-10) {
842	allGood = false;
843	break;
844	} else if (upper[iColumn] < 1.0e10 && fabs(upper[iColumn] - floor(upper[iColumn] + 0.5)) > 1.0e-10) {
845	allGood = false;
846	break;
847	}
848	bool singletonRow = true;
849	bool equality = false;
850	for (CoinBigIndex j = start; j < end; j++) {
851	int iRow = row[j];
852	if (count[iRow] > 1)
853	singletonRow = false;
854	else if (rowLower[iRow] == rowUpper[iRow])
855	equality = true;
856	double rhsValue = rhs[iRow];
857	double lowerValue = rowLower[iRow];
858	double upperValue = rowUpper[iRow];
859	if (rhsValue < 1.0e20) {
860	if (lowerValue > -1.0e20)
861	lowerValue -= rhsValue;
862	if (upperValue < 1.0e20)
863	upperValue -= rhsValue;
864	}
865	if (fabs(rhsValue) > 1.0e20 \|\| fabs(rhsValue - floor(rhsValue + 0.5)) > 1.0e-10
866	\|\| fabs(element[j]) != 1.0) {
867	// no good
868	allGood = false;
869	break;
870	}
871	if (element[j] > 0.0) {
872	if (lowerValue > -1.0e20 && fabs(lowerValue - floor(lowerValue + 0.5)) > 1.0e-10) {
873	// no good
874	allGood = false;
875	break;
876	}
877	} else {
878	if (upperValue < 1.0e20 && fabs(upperValue - floor(upperValue + 0.5)) > 1.0e-10) {
879	// no good
880	allGood = false;
881	break;
882	}
883	}
884	}
885	if (!singletonRow && end > start + 1 && !equality)
886	allGood = false;
887	if (!allGood)
888	break;
889	} else {
890	// wants to be as high as possible
891	if (upper[iColumn] > 1.0e10 \|\| fabs(upper[iColumn] - floor(upper[iColumn] + 0.5)) > 1.0e-10) {
892	allGood = false;
893	break;
894	} else if (lower[iColumn] > -1.0e10 && fabs(lower[iColumn] - floor(lower[iColumn] + 0.5)) > 1.0e-10) {
895	allGood = false;
896	break;
897	}
898	bool singletonRow = true;
899	bool equality = false;
900	for (CoinBigIndex j = start; j < end; j++) {
901	int iRow = row[j];
902	if (count[iRow] > 1)
903	singletonRow = false;
904	else if (rowLower[iRow] == rowUpper[iRow])
905	equality = true;
906	double rhsValue = rhs[iRow];
907	double lowerValue = rowLower[iRow];
908	double upperValue = rowUpper[iRow];
909	if (rhsValue < 1.0e20) {
910	if (lowerValue > -1.0e20)
911	lowerValue -= rhsValue;
912	if (upperValue < 1.0e20)
913	upperValue -= rhsValue;
914	}
915	if (fabs(rhsValue) > 1.0e20 \|\| fabs(rhsValue - floor(rhsValue + 0.5)) > 1.0e-10
916	\|\| fabs(element[j]) != 1.0) {
917	// no good
918	allGood = false;
919	break;
920	}
921	if (element[j] < 0.0) {
922	if (lowerValue > -1.0e20 && fabs(lowerValue - floor(lowerValue + 0.5)) > 1.0e-10) {
923	// no good
924	allGood = false;
925	break;
926	}
927	} else {
928	if (upperValue < 1.0e20 && fabs(upperValue - floor(upperValue + 0.5)) > 1.0e-10) {
929	// no good
930	allGood = false;
931	break;
932	}
933	}
934	}
935	if (!singletonRow && end > start + 1 && !equality)
936	allGood = false;
937	if (!allGood)
938	break;
939	}
940	}
941	}
942	}
943	delete [] count;
944	if (allGood) {
945	#if COIN_DEVELOP>1
946	if (numberObj)
947	printf("YYYY analysis says all continuous with costs will be integer\n");
948	#endif
949	continuousMultiplier = 1.0;
950	}
951	}
952	delete [] rhs;
953	}
954	/*
955	Take a first scan to see if there are unfixed continuous variables in the
956	objective. If so, the minimum objective change could be arbitrarily small.
957	Also pick off the maximum coefficient of an unfixed integer variable.
958
959	If the objective is found to contain only integer variables, set the
960	fathoming discipline to strict.
961	*/
962	double maximumCost = 0.0 ;
963	//double trueIncrement=0.0;
964	int iColumn ;
965	int numberColumns = getNumCols() ;
966	double scaleFactor = 1.0; // due to rhs etc
967	/*
968	Original model did not have integer bounds.
969	*/
970	if ((specialOptions_&65536) == 0) {
971	/* be on safe side (later look carefully as may be able to
972	to get 0.5 say if bounds are multiples of 0.5 */
973	for (iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
974	if (upper[iColumn] > lower[iColumn] + 1.0e-8) {
975	double value;
976	value = fabs(lower[iColumn]);
977	if (floor(value + 0.5) != value) {
978	scaleFactor = CoinMin(scaleFactor, 0.5);
979	if (floor(2.0value + 0.5) != 2.0value) {
980	scaleFactor = CoinMin(scaleFactor, 0.25);
981	if (floor(4.0value + 0.5) != 4.0value) {
982	scaleFactor = 0.0;
983	}
984	}
985	}
986	value = fabs(upper[iColumn]);
987	if (floor(value + 0.5) != value) {
988	scaleFactor = CoinMin(scaleFactor, 0.5);
989	if (floor(2.0value + 0.5) != 2.0value) {
990	scaleFactor = CoinMin(scaleFactor, 0.25);
991	if (floor(4.0value + 0.5) != 4.0value) {
992	scaleFactor = 0.0;
993	}
994	}
995	}
996	}
997	}
998	}
999	bool possibleMultiple = continuousMultiplier != 0.0 && scaleFactor != 0.0 ;
1000	if (possibleMultiple) {
1001	for (iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
1002	if (upper[iColumn] > lower[iColumn] + 1.0e-8) {
1003	maximumCost = CoinMax(maximumCost, fabs(objective[iColumn])) ;
1004	}
1005	}
1006	}
1007	setIntParam(CbcModel::CbcFathomDiscipline, possibleMultiple) ;
1008	/*
1009	If a nontrivial increment is possible, try and figure it out. We're looking
1010	for gcd(c<j>) for all c<j> that are coefficients of unfixed integer
1011	variables. Since the c<j> might not be integers, try and inflate them
1012	sufficiently that they look like integers (and we'll deflate the gcd
1013	later).
1014
1015	2520.0 is used as it is a nice multiple of 2,3,5,7
1016	*/
1017	if (possibleMultiple && maximumCost) {
1018	int increment = 0 ;
1019	double multiplier = 2520.0 ;
1020	while (10.0multipliermaximumCost < 1.0e8)
1021	multiplier *= 10.0 ;
1022	int bigIntegers = 0; // Count of large costs which are integer
1023	for (iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
1024	if (upper[iColumn] > lower[iColumn] + 1.0e-8) {
1025	double objValue = fabs(objective[iColumn]);
1026	if (!isInteger(iColumn)) {
1027	if (!coeffMultiplier)
1028	objValue *= continuousMultiplier;
1029	else
1030	objValue *= coeffMultiplier[iColumn];
1031	}
1032	if (objValue) {
1033	double value = objValue * multiplier ;
1034	if (value < 2.1e9) {
1035	int nearest = static_cast<int> (floor(value + 0.5)) ;
1036	if (fabs(value - floor(value + 0.5)) > 1.0e-8) {
1037	increment = 0 ;
1038	break ;
1039	} else if (!increment) {
1040	increment = nearest ;
1041	} else {
1042	increment = gcd(increment, nearest) ;
1043	}
1044	} else {
1045	// large value - may still be multiple of 1.0
1046	if (fabs(objValue - floor(objValue + 0.5)) > 1.0e-8) {
1047	increment = 0;
1048	break;
1049	} else {
1050	bigIntegers++;
1051	}
1052	}
1053	}
1054	}
1055	}
1056	delete [] coeffMultiplier;
1057	/*
1058	If the increment beats the current value for objective change, install it.
1059	*/
1060	if (increment) {
1061	double value = increment ;
1062	double cutoff = getDblParam(CbcModel::CbcCutoffIncrement) ;
1063	if (bigIntegers) {
1064	// allow for 1.0
1065	increment = gcd(increment, static_cast<int> (multiplier));
1066	value = increment;
1067	}
1068	value /= multiplier ;
1069	value *= scaleFactor;
1070	//trueIncrement=CoinMax(cutoff,value);;
1071	if (value*0.999 > cutoff) {
1072	messageHandler()->message(CBC_INTEGERINCREMENT,
1073	messages())
1074	<< value << CoinMessageEol ;
1075	setDblParam(CbcModel::CbcCutoffIncrement, CoinMax(value*0.999,value-1.0e-4)) ;
1076	}
1077	}
1078	}
1079
1080	return ;
1081	}
1082
1083	/*
1084	saveModel called (carved out of) BranchandBound
1085	*/
1086	void CbcModel::saveModel(OsiSolverInterface * saveSolver, double * checkCutoffForRestart, bool * feasible)
1087	{
1088	if (saveSolver && (specialOptions_&32768) != 0) {
1089	// See if worth trying reduction
1090	*checkCutoffForRestart = getCutoff();
1091	bool tryNewSearch = solverCharacteristics_->reducedCostsAccurate() &&
1092	(*checkCutoffForRestart < 1.0e20);
1093	int numberColumns = getNumCols();
1094	if (tryNewSearch) {
1095	#if CBC_USEFUL_PRINTING>1
1096	printf("after %d nodes, cutoff %g - looking\n",
1097	numberNodes_, getCutoff());
1098	#endif
1099	saveSolver->resolve();
1100	double direction = saveSolver->getObjSense() ;
1101	double gap = checkCutoffForRestart - saveSolver->getObjValue() direction ;
1102	double tolerance;
1103	saveSolver->getDblParam(OsiDualTolerance, tolerance) ;
1104	if (gap <= 0.0)
1105	gap = tolerance;
1106	gap += 100.0 * tolerance;
1107	double integerTolerance = getDblParam(CbcIntegerTolerance) ;
1108
1109	const double *lower = saveSolver->getColLower() ;
1110	const double *upper = saveSolver->getColUpper() ;
1111	const double *solution = saveSolver->getColSolution() ;
1112	const double *reducedCost = saveSolver->getReducedCost() ;
1113
1114	int numberFixed = 0 ;
1115	int numberFixed2 = 0;
1116	for (int i = 0 ; i < numberIntegers_ ; i++) {
1117	int iColumn = integerVariable_[i] ;
1118	double djValue = direction * reducedCost[iColumn] ;
1119	if (upper[iColumn] - lower[iColumn] > integerTolerance) {
1120	if (solution[iColumn] < lower[iColumn] + integerTolerance && djValue > gap) {
1121	saveSolver->setColUpper(iColumn, lower[iColumn]) ;
1122	numberFixed++ ;
1123	} else if (solution[iColumn] > upper[iColumn] - integerTolerance && -djValue > gap) {
1124	saveSolver->setColLower(iColumn, upper[iColumn]) ;
1125	numberFixed++ ;
1126	}
1127	} else {
1128	numberFixed2++;
1129	}
1130	}
1131	#ifdef COIN_DEVELOP
1132	/*
1133	We're debugging. (specialOptions 1)
1134	*/
1135	if ((specialOptions_&1) != 0) {
1136	const OsiRowCutDebugger *debugger = saveSolver->getRowCutDebugger() ;
1137	if (debugger) {
1138	printf("Contains optimal\n") ;
1139	OsiSolverInterface * temp = saveSolver->clone();
1140	const double * solution = debugger->optimalSolution();
1141	const double *lower = temp->getColLower() ;
1142	const double *upper = temp->getColUpper() ;
1143	int n = temp->getNumCols();
1144	for (int i = 0; i < n; i++) {
1145	if (temp->isInteger(i)) {
1146	double value = floor(solution[i] + 0.5);
1147	assert (value >= lower[i] && value <= upper[i]);
1148	temp->setColLower(i, value);
1149	temp->setColUpper(i, value);
1150	}
1151	}
1152	temp->writeMps("reduced_fix");
1153	delete temp;
1154	saveSolver->writeMps("reduced");
1155	} else {
1156	abort();
1157	}
1158	}
1159	printf("Restart could fix %d integers (%d already fixed)\n",
1160	numberFixed + numberFixed2, numberFixed2);
1161	#endif
1162	numberFixed += numberFixed2;
1163	if (numberFixed*20 < numberColumns)
1164	tryNewSearch = false;
1165	}
1166	if (tryNewSearch) {
1167	// back to solver without cuts?
1168	OsiSolverInterface * solver2 = continuousSolver_->clone();
1169	const double *lower = saveSolver->getColLower() ;
1170	const double *upper = saveSolver->getColUpper() ;
1171	for (int i = 0 ; i < numberIntegers_ ; i++) {
1172	int iColumn = integerVariable_[i] ;
1173	solver2->setColLower(iColumn, lower[iColumn]);
1174	solver2->setColUpper(iColumn, upper[iColumn]);
1175	}
1176	// swap
1177	delete saveSolver;
1178	saveSolver = solver2;
1179	double * newSolution = new double[numberColumns];
1180	double objectiveValue = *checkCutoffForRestart;
1181	CbcSerendipity heuristic(*this);
1182	if (bestSolution_)
1183	heuristic.setInputSolution(bestSolution_, bestObjective_);
1184	heuristic.setFractionSmall(0.9);
1185	heuristic.setFeasibilityPumpOptions(1008013);
1186	// Use numberNodes to say how many are original rows
1187	heuristic.setNumberNodes(continuousSolver_->getNumRows());
1188	#ifdef COIN_DEVELOP
1189	if (continuousSolver_->getNumRows() <
1190	saveSolver->getNumRows())
1191	printf("%d rows added ZZZZZ\n",
1192	solver_->getNumRows() - continuousSolver_->getNumRows());
1193	#endif
1194	int returnCode = heuristic.smallBranchAndBound(saveSolver,
1195	-1, newSolution,
1196	objectiveValue,
1197	*checkCutoffForRestart, "Reduce");
1198	if (returnCode < 0) {
1199	#ifdef COIN_DEVELOP
1200	printf("Restart - not small enough to do search after fixing\n");
1201	#endif
1202	delete [] newSolution;
1203	} else {
1204	if ((returnCode&1) != 0) {
1205	// increment number of solutions so other heuristics can test
1206	numberSolutions_++;
1207	numberHeuristicSolutions_++;
1208	lastHeuristic_ = NULL;
1209	setBestSolution(CBC_ROUNDING, objectiveValue, newSolution) ;
1210	}
1211	delete [] newSolution;
1212	*feasible = false; // stop search
1213	}
1214	#if 0 // probably not needed def CBC_THREAD
1215	if (master_) {
1216	lockThread();
1217	if (parallelMode() > 0) {
1218	while (master_->waitForThreadsInTree(0)) {
1219	lockThread();
1220	double dummyBest;
1221	tree_->cleanTree(this, -COIN_DBL_MAX, dummyBest) ;
1222	//unlockThread();
1223	}
1224	}
1225	master_->waitForThreadsInTree(2);
1226	delete master_;
1227	master_ = NULL;
1228	masterThread_ = NULL;
1229	}
1230	#endif
1231	}
1232	}
1233	}
1234	/*
1235	Adds integers, called from BranchandBound()
1236	*/
1237	void CbcModel::AddIntegers()
1238	{
1239	int numberColumns = continuousSolver_->getNumCols();
1240	int numberRows = continuousSolver_->getNumRows();
1241	int numberOriginalIntegers = numberIntegers_;
1242	int * del = new int [CoinMax(numberColumns, numberRows)];
1243	int * original = new int [numberColumns];
1244	char * possibleRow = new char [numberRows];
1245	{
1246	const CoinPackedMatrix * rowCopy = continuousSolver_->getMatrixByRow();
1247	const int * column = rowCopy->getIndices();
1248	const int * rowLength = rowCopy->getVectorLengths();
1249	const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
1250	const double * rowLower = continuousSolver_->getRowLower();
1251	const double * rowUpper = continuousSolver_->getRowUpper();
1252	const double * element = rowCopy->getElements();
1253	for (int i = 0; i < numberRows; i++) {
1254	int nLeft = 0;
1255	bool possible = false;
1256	if (rowLower[i] < -1.0e20) {
1257	double value = rowUpper[i];
1258	if (fabs(value - floor(value + 0.5)) < 1.0e-8)
1259	possible = true;
1260	} else if (rowUpper[i] > 1.0e20) {
1261	double value = rowLower[i];
1262	if (fabs(value - floor(value + 0.5)) < 1.0e-8)
1263	possible = true;
1264	} else {
1265	double value = rowUpper[i];
1266	if (rowLower[i] == rowUpper[i] &&
1267	fabs(value - floor(value + 0.5)) < 1.0e-8)
1268	possible = true;
1269	}
1270	double allSame = (possible) ? 0.0 : -1.0;
1271	for (CoinBigIndex j = rowStart[i];
1272	j < rowStart[i] + rowLength[i]; j++) {
1273	int iColumn = column[j];
1274	if (continuousSolver_->isInteger(iColumn)) {
1275	if (fabs(element[j]) != 1.0)
1276	possible = false;
1277	} else {
1278	nLeft++;
1279	if (!allSame) {
1280	allSame = fabs(element[j]);
1281	} else if (allSame>0.0) {
1282	if (allSame!=fabs(element[j]))
1283	allSame = -1.0;
1284	}
1285	}
1286	}
1287	if (nLeft == rowLength[i] && allSame > 0.0)
1288	possibleRow[i] = 2;
1289	else if (possible \|\| !nLeft)
1290	possibleRow[i] = 1;
1291	else
1292	possibleRow[i] = 0;
1293	}
1294	}
1295	int nDel = 0;
1296	for (int i = 0; i < numberColumns; i++) {
1297	original[i] = i;
1298	if (continuousSolver_->isInteger(i))
1299	del[nDel++] = i;
1300	}
1301	{
1302	// we must not exclude current best solution (rounding errors)
1303	// also not if large values
1304	const int * row = continuousSolver_->getMatrixByCol()->getIndices();
1305	const CoinBigIndex * columnStart = continuousSolver_->getMatrixByCol()->getVectorStarts();
1306	const int * columnLength = continuousSolver_->getMatrixByCol()->getVectorLengths();
1307	const double * solution = continuousSolver_->getColSolution();
1308	for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
1309	if (!continuousSolver_->isInteger(iColumn)) {
1310	double value = bestSolution_ ? bestSolution_[iColumn] : 0.0;
1311	double value2 = solution[iColumn];
1312	if (fabs(value-floor(value+0.5))>1.0e-8 \|\|
1313	fabs(value2)>1.0e3) {
1314	CoinBigIndex start = columnStart[iColumn];
1315	CoinBigIndex end = start + columnLength[iColumn];
1316	for (CoinBigIndex j = start; j < end; j++) {
1317	int iRow = row[j];
1318	possibleRow[iRow]=0;
1319	}
1320	}
1321	}
1322	}
1323	}
1324	int nExtra = 0;
1325	OsiSolverInterface * copy1 = continuousSolver_->clone();
1326	int nPass = 0;
1327	while (nDel && nPass < 10) {
1328	nPass++;
1329	OsiSolverInterface * copy2 = copy1->clone();
1330	int nLeft = 0;
1331	for (int i = 0; i < nDel; i++)
1332	original[del[i]] = -1;
1333	for (int i = 0; i < numberColumns; i++) {
1334	int kOrig = original[i];
1335	if (kOrig >= 0)
1336	original[nLeft++] = kOrig;
1337	}
1338	assert (nLeft == numberColumns - nDel);
1339	copy2->deleteCols(nDel, del);
1340	numberColumns = copy2->getNumCols();
1341	const CoinPackedMatrix * rowCopy = copy2->getMatrixByRow();
1342	numberRows = rowCopy->getNumRows();
1343	const int * column = rowCopy->getIndices();
1344	const int * rowLength = rowCopy->getVectorLengths();
1345	const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
1346	const double * rowLower = copy2->getRowLower();
1347	const double * rowUpper = copy2->getRowUpper();
1348	const double * element = rowCopy->getElements();
1349	const CoinPackedMatrix * columnCopy = copy2->getMatrixByCol();
1350	const int * columnLength = columnCopy->getVectorLengths();
1351	nDel = 0;
1352	// Could do gcd stuff on ones with costs
1353	for (int i = 0; i < numberRows; i++) {
1354	if (!rowLength[i]) {
1355	del[nDel++] = i;
1356	} else if (possibleRow[i]) {
1357	if (rowLength[i] == 1) {
1358	CoinBigIndex k = rowStart[i];
1359	int iColumn = column[k];
1360	if (!copy2->isInteger(iColumn)) {
1361	double mult = 1.0 / fabs(element[k]);
1362	if (rowLower[i] < -1.0e20) {
1363	// treat rhs as multiple of 1 unless elements all same
1364	double value = ((possibleRow[i]==2) ? rowUpper[i] : 1.0) * mult;
1365	if (fabs(value - floor(value + 0.5)) < 1.0e-8) {
1366	del[nDel++] = i;
1367	if (columnLength[iColumn] == 1) {
1368	copy2->setInteger(iColumn);
1369	int kOrig = original[iColumn];
1370	setOptionalInteger(kOrig);
1371	}
1372	}
1373	} else if (rowUpper[i] > 1.0e20) {
1374	// treat rhs as multiple of 1 unless elements all same
1375	double value = ((possibleRow[i]==2) ? rowLower[i] : 1.0) * mult;
1376	if (fabs(value - floor(value + 0.5)) < 1.0e-8) {
1377	del[nDel++] = i;
1378	if (columnLength[iColumn] == 1) {
1379	copy2->setInteger(iColumn);
1380	int kOrig = original[iColumn];
1381	setOptionalInteger(kOrig);
1382	}
1383	}
1384	} else {
1385	// treat rhs as multiple of 1 unless elements all same
1386	double value = ((possibleRow[i]==2) ? rowUpper[i] : 1.0) * mult;
1387	if (rowLower[i] == rowUpper[i] &&
1388	fabs(value - floor(value + 0.5)) < 1.0e-8) {
1389	del[nDel++] = i;
1390	copy2->setInteger(iColumn);
1391	int kOrig = original[iColumn];
1392	setOptionalInteger(kOrig);
1393	}
1394	}
1395	}
1396	} else {
1397	// only if all singletons
1398	bool possible = false;
1399	if (rowLower[i] < -1.0e20) {
1400	double value = rowUpper[i];
1401	if (fabs(value - floor(value + 0.5)) < 1.0e-8)
1402	possible = true;
1403	} else if (rowUpper[i] > 1.0e20) {
1404	double value = rowLower[i];
1405	if (fabs(value - floor(value + 0.5)) < 1.0e-8)
1406	possible = true;
1407	} else {
1408	double value = rowUpper[i];
1409	if (rowLower[i] == rowUpper[i] &&
1410	fabs(value - floor(value + 0.5)) < 1.0e-8)
1411	possible = true;
1412	}
1413	if (possible) {
1414	for (CoinBigIndex j = rowStart[i];
1415	j < rowStart[i] + rowLength[i]; j++) {
1416	int iColumn = column[j];
1417	if (columnLength[iColumn] != 1 \|\| fabs(element[j]) != 1.0) {
1418	possible = false;
1419	break;
1420	}
1421	}
1422	if (possible) {
1423	for (CoinBigIndex j = rowStart[i];
1424	j < rowStart[i] + rowLength[i]; j++) {
1425	int iColumn = column[j];
1426	if (!copy2->isInteger(iColumn)) {
1427	copy2->setInteger(iColumn);
1428	int kOrig = original[iColumn];
1429	setOptionalInteger(kOrig);
1430	}
1431	}
1432	del[nDel++] = i;
1433	}
1434	}
1435	}
1436	}
1437	}
1438	if (nDel) {
1439	copy2->deleteRows(nDel, del);
1440	// pack down possible
1441	int n=0;
1442	for (int i=0;i<nDel;i++)
1443	possibleRow[del[i]]=-1;
1444	for (int i=0;i<numberRows;i++) {
1445	if (possibleRow[i]>=0)
1446	possibleRow[n++]=possibleRow[i];
1447	}
1448	}
1449	if (nDel != numberRows) {
1450	nDel = 0;
1451	for (int i = 0; i < numberColumns; i++) {
1452	if (copy2->isInteger(i)) {
1453	del[nDel++] = i;
1454	nExtra++;
1455	}
1456	}
1457	} else {
1458	nDel = 0;
1459	}
1460	delete copy1;
1461	copy1 = copy2->clone();
1462	delete copy2;
1463	}
1464	// See if what's left is a network
1465	bool couldBeNetwork = false;
1466	if (copy1->getNumRows() && copy1->getNumCols()) {
1467	#ifdef COIN_HAS_CLP
1468	OsiClpSolverInterface * clpSolver
1469	= dynamic_cast<OsiClpSolverInterface *> (copy1);
1470	if (false && clpSolver) {
1471	numberRows = clpSolver->getNumRows();
1472	char * rotate = new char[numberRows];
1473	int n = clpSolver->getModelPtr()->findNetwork(rotate, 1.0);
1474	delete [] rotate;
1475	#if CBC_USEFUL_PRINTING>1
1476	printf("INTA network %d rows out of %d\n", n, numberRows);
1477	#endif
1478	if (CoinAbs(n) == numberRows) {
1479	couldBeNetwork = true;
1480	for (int i = 0; i < numberRows; i++) {
1481	if (!possibleRow[i]) {
1482	couldBeNetwork = false;
1483	#if CBC_USEFUL_PRINTING>1
1484	printf("but row %d is bad\n", i);
1485	#endif
1486	break;
1487	}
1488	}
1489	}
1490	} else
1491	#endif
1492	{
1493	numberColumns = copy1->getNumCols();
1494	numberRows = copy1->getNumRows();
1495	const double * rowLower = copy1->getRowLower();
1496	const double * rowUpper = copy1->getRowUpper();
1497	couldBeNetwork = true;
1498	for (int i = 0; i < numberRows; i++) {
1499	if (rowLower[i] > -1.0e20 && fabs(rowLower[i] - floor(rowLower[i] + 0.5)) > 1.0e-12) {
1500	couldBeNetwork = false;
1501	break;
1502	}
1503	if (rowUpper[i] < 1.0e20 && fabs(rowUpper[i] - floor(rowUpper[i] + 0.5)) > 1.0e-12) {
1504	couldBeNetwork = false;
1505	break;
1506	}
1507	if (possibleRow[i]==0) {
1508	couldBeNetwork = false;
1509	break;
1510	}
1511	}
1512	if (couldBeNetwork) {
1513	const CoinPackedMatrix * matrixByCol = copy1->getMatrixByCol();
1514	const double * element = matrixByCol->getElements();
1515	//const int * row = matrixByCol->getIndices();
1516	const CoinBigIndex * columnStart = matrixByCol->getVectorStarts();
1517	const int * columnLength = matrixByCol->getVectorLengths();
1518	for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
1519	CoinBigIndex start = columnStart[iColumn];
1520	CoinBigIndex end = start + columnLength[iColumn];
1521	if (end > start + 2) {
1522	couldBeNetwork = false;
1523	break;
1524	}
1525	int type = 0;
1526	for (CoinBigIndex j = start; j < end; j++) {
1527	double value = element[j];
1528	if (fabs(value) != 1.0) {
1529	couldBeNetwork = false;
1530	break;
1531	} else if (value == 1.0) {
1532	if ((type&1) == 0)
1533	type \|= 1;
1534	else
1535	type = 7;
1536	} else if (value == -1.0) {
1537	if ((type&2) == 0)
1538	type \|= 2;
1539	else
1540	type = 7;
1541	}
1542	}
1543	if (type > 3) {
1544	couldBeNetwork = false;
1545	break;
1546	}
1547	}
1548	}
1549	}
1550	}
1551	if (couldBeNetwork) {
1552	for (int i = 0; i < numberColumns; i++)
1553	setOptionalInteger(original[i]);
1554	}
1555	if (nExtra \|\| couldBeNetwork) {
1556	numberColumns = copy1->getNumCols();
1557	numberRows = copy1->getNumRows();
1558	if (!numberColumns \|\| !numberRows) {
1559	int numberColumns = solver_->getNumCols();
1560	for (int i = 0; i < numberColumns; i++)
1561	assert(solver_->isInteger(i));
1562	}
1563	#if CBC_USEFUL_PRINTING>1
1564	if (couldBeNetwork \|\| nExtra)
1565	printf("INTA %d extra integers, %d left%s\n", nExtra,
1566	numberColumns,
1567	couldBeNetwork ? ", all network" : "");
1568	#endif
1569	findIntegers(true, 2);
1570	convertToDynamic();
1571	}
1572	#if CBC_USEFUL_PRINTING>1
1573	if (!couldBeNetwork && copy1->getNumCols() &&
1574	copy1->getNumRows()) {
1575	printf("INTA %d rows and %d columns remain\n",
1576	copy1->getNumRows(), copy1->getNumCols());
1577	if (copy1->getNumCols() < 200) {
1578	copy1->writeMps("moreint");
1579	printf("INTA Written remainder to moreint.mps.gz %d rows %d cols\n",
1580	copy1->getNumRows(), copy1->getNumCols());
1581	}
1582	}
1583	#endif
1584	delete copy1;
1585	delete [] del;
1586	delete [] original;
1587	delete [] possibleRow;
1588	// double check increment
1589	analyzeObjective();
1590	// If any changes - tell code
1591	if(numberOriginalIntegers<numberIntegers_)
1592	synchronizeModel();
1593	}
1594	/**
1595	\todo
1596	Normally, it looks like we enter here from command dispatch in the main
1597	routine, after calling the solver for an initial solution
1598	(CbcModel::initialSolve, which simply calls the solver's initialSolve
1599	routine.) The first thing we do is call resolve. Presumably there are
1600	circumstances where this is nontrivial? There's also a call from
1601	CbcModel::originalModel (tied up with integer presolve), which should be
1602	checked.
1603
1604	*/
1605
1606	#ifdef CONFLICT_CUTS
1607	#if PRINT_CONFLICT==1
1608	static int numberConflictCuts=0;
1609	static int lastNumberConflictCuts=0;
1610	static double lengthConflictCuts=0.0;
1611	#endif
1612	#endif
1613	/*
1614	The overall flow can be divided into three stages:
1615	* Prep: Check that the lp relaxation remains feasible at the root. If so,
1616	do all the setup for B&C.
1617	* Process the root node: Generate cuts, apply heuristics, and in general do
1618	the best we can to resolve the problem without B&C.
1619	* Do B&C search until we hit a limit or exhaust the search tree.
1620
1621	Keep in mind that in general there is no node in the search tree that
1622	corresponds to the active subproblem. The active subproblem is represented
1623	by the current state of the model, of the solver, and of the constraint
1624	system held by the solver.
1625	*/
1626	#ifdef COIN_HAS_CPX
1627	#include "OsiCpxSolverInterface.hpp"
1628	#include "cplex.h"
1629	#endif
1630	void CbcModel::branchAndBound(int doStatistics)
1631
1632	{
1633	if (!parentModel_) {
1634	/*
1635	Capture a time stamp before we start (unless set).
1636	*/
1637	if (!dblParam_[CbcStartSeconds]) {
1638	if (!useElapsedTime())
1639	dblParam_[CbcStartSeconds] = CoinCpuTime();
1640	else
1641	dblParam_[CbcStartSeconds] = CoinGetTimeOfDay();
1642	}
1643	}
1644	dblParam_[CbcSmallestChange] = COIN_DBL_MAX;
1645	dblParam_[CbcSumChange] = 0.0;
1646	dblParam_[CbcLargestChange] = 0.0;
1647	intParam_[CbcNumberBranches] = 0;
1648	double lastBestPossibleObjective=-COIN_DBL_MAX;
1649	// when to check for restart
1650	int nextCheckRestart=50;
1651	// Force minimization !!!!
1652	bool flipObjective = (solver_->getObjSense()<0.0);
1653	if (flipObjective)
1654	flipModel();
1655	dblParam_[CbcOptimizationDirection] = 1.0; // was solver_->getObjSense();
1656	strongInfo_[0] = 0;
1657	strongInfo_[1] = 0;
1658	strongInfo_[2] = 0;
1659	strongInfo_[3] = 0;
1660	strongInfo_[4] = 0;
1661	strongInfo_[5] = 0;
1662	strongInfo_[6] = 0;
1663	numberStrongIterations_ = 0;
1664	currentNode_ = NULL;
1665	// See if should do cuts old way
1666	if (parallelMode() < 0) {
1667	specialOptions_ \|= 4096 + 8192;
1668	} else if (parallelMode() > 0) {
1669	specialOptions_ \|= 4096;
1670	}
1671	int saveMoreSpecialOptions = moreSpecialOptions_;
1672	if (dynamic_cast<CbcTreeLocal *> (tree_))
1673	specialOptions_ \|= 4096 + 8192;
1674	#ifdef COIN_HAS_CLP
1675	{
1676	OsiClpSolverInterface * clpSolver
1677	= dynamic_cast<OsiClpSolverInterface *> (solver_);
1678	if (clpSolver) {
1679	// pass in disaster handler
1680	CbcDisasterHandler handler(this);
1681	clpSolver->passInDisasterHandler(&handler);
1682	// Initialise solvers seed (unless users says not)
1683	if ((specialOptions_&4194304)==0)
1684	clpSolver->getModelPtr()->setRandomSeed(1234567);
1685	#ifdef JJF_ZERO
1686	// reduce factorization frequency
1687	int frequency = clpSolver->getModelPtr()->factorizationFrequency();
1688	clpSolver->getModelPtr()->setFactorizationFrequency(CoinMin(frequency, 120));
1689	#endif
1690	}
1691	}
1692	#endif
1693	// original solver (only set if pre-processing)
1694	OsiSolverInterface * originalSolver = NULL;
1695	int numberOriginalObjects = numberObjects_;
1696	OsiObject ** originalObject = NULL;
1697	// Save whether there were any objects
1698	bool noObjects = (numberObjects_ == 0);
1699	// Set up strategies
1700	/*
1701	See if the user has supplied a strategy object and deal with it if present.
1702	The call to setupOther will set numberStrong_ and numberBeforeTrust_, and
1703	perform integer preprocessing, if requested.
1704
1705	We need to hang on to a pointer to solver_. setupOther will assign a
1706	preprocessed solver to model, but will instruct assignSolver not to trash the
1707	existing one.
1708	*/
1709	if (strategy_) {
1710	// May do preprocessing
1711	originalSolver = solver_;
1712	strategy_->setupOther(*this);
1713	if (strategy_->preProcessState()) {
1714	// pre-processing done
1715	if (strategy_->preProcessState() < 0) {
1716	// infeasible (or unbounded)
1717	status_ = 0 ;
1718	if (!solver_->isProvenDualInfeasible()) {
1719	handler_->message(CBC_INFEAS, messages_) << CoinMessageEol ;
1720	secondaryStatus_ = 1;
1721	} else {
1722	handler_->message(CBC_UNBOUNDED,
1723	messages_) << CoinMessageEol ;
1724	secondaryStatus_ = 7;
1725	}
1726	originalContinuousObjective_ = COIN_DBL_MAX;
1727	if (flipObjective)
1728	flipModel();
1729	return ;
1730	} else if (numberObjects_ && object_) {
1731	numberOriginalObjects = numberObjects_;
1732	// redo sequence
1733	numberIntegers_ = 0;
1734	int numberColumns = getNumCols();
1735	int nOrig = originalSolver->getNumCols();
1736	CglPreProcess * process = strategy_->process();
1737	assert (process);
1738	const int * originalColumns = process->originalColumns();
1739	// allow for cliques etc
1740	nOrig = CoinMax(nOrig, originalColumns[numberColumns-1] + 1);
1741	// try and redo debugger
1742	OsiRowCutDebugger * debugger = const_cast<OsiRowCutDebugger *> (solver_->getRowCutDebuggerAlways());
1743	if (debugger) {
1744	if (numberColumns<=debugger->numberColumns())
1745	debugger->redoSolution(numberColumns, originalColumns);
1746	else
1747	debugger=NULL; // no idea how to handle (SOS?)
1748	}
1749	// User-provided solution might have been best. Synchronise.
1750	if (bestSolution_) {
1751	// need to redo - in case no better found in BAB
1752	// just get integer part right
1753	for (int i = 0; i < numberColumns; i++) {
1754	int jColumn = originalColumns[i];
1755	bestSolution_[i] = bestSolution_[jColumn];
1756	}
1757	}
1758	originalObject = object_;
1759	// object number or -1
1760	int * temp = new int[nOrig];
1761	int iColumn;
1762	for (iColumn = 0; iColumn < nOrig; iColumn++)
1763	temp[iColumn] = -1;
1764	int iObject;
1765	int nNonInt = 0;
1766	for (iObject = 0; iObject < numberOriginalObjects; iObject++) {
1767	iColumn = originalObject[iObject]->columnNumber();
1768	if (iColumn < 0) {
1769	nNonInt++;
1770	} else {
1771	temp[iColumn] = iObject;
1772	}
1773	}
1774	int numberNewIntegers = 0;
1775	int numberOldIntegers = 0;
1776	int numberOldOther = 0;
1777	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1778	int jColumn = originalColumns[iColumn];
1779	if (temp[jColumn] >= 0) {
1780	int iObject = temp[jColumn];
1781	CbcSimpleInteger * obj =
1782	dynamic_cast <CbcSimpleInteger *>(originalObject[iObject]) ;
1783	if (obj)
1784	numberOldIntegers++;
1785	else
1786	numberOldOther++;
1787	} else if (isInteger(iColumn)) {
1788	numberNewIntegers++;
1789	}
1790	}
1791	/*
1792	Allocate an array to hold the indices of the integer variables.
1793	Make a large enough array for all objects
1794	*/
1795	numberObjects_ = numberNewIntegers + numberOldIntegers + numberOldOther + nNonInt;
1796	object_ = new OsiObject * [numberObjects_];
1797	delete [] integerVariable_;
1798	integerVariable_ = new int [numberNewIntegers+numberOldIntegers];
1799	/*
1800	Walk the variables again, filling in the indices and creating objects for
1801	the integer variables. Initially, the objects hold the index and upper &
1802	lower bounds.
1803	*/
1804	numberIntegers_ = 0;
1805	int n = originalColumns[numberColumns-1] + 1;
1806	int * backward = new int[n];
1807	int i;
1808	for ( i = 0; i < n; i++)
1809	backward[i] = -1;
1810	for (i = 0; i < numberColumns; i++)
1811	backward[originalColumns[i]] = i;
1812	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1813	int jColumn = originalColumns[iColumn];
1814	if (temp[jColumn] >= 0) {
1815	int iObject = temp[jColumn];
1816	CbcSimpleInteger * obj =
1817	dynamic_cast <CbcSimpleInteger *>(originalObject[iObject]) ;
1818	if (obj) {
1819	object_[numberIntegers_] = originalObject[iObject]->clone();
1820	// redo ids etc
1821	//object_[numberIntegers_]->resetSequenceEtc(numberColumns,originalColumns);
1822	object_[numberIntegers_]->resetSequenceEtc(numberColumns, backward);
1823	integerVariable_[numberIntegers_++] = iColumn;
1824	}
1825	} else if (isInteger(iColumn)) {
1826	object_[numberIntegers_] =
1827	new CbcSimpleInteger(this, iColumn);
1828	integerVariable_[numberIntegers_++] = iColumn;
1829	}
1830	}
1831	delete [] backward;
1832	numberObjects_ = numberIntegers_;
1833	// Now append other column stuff
1834	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1835	int jColumn = originalColumns[iColumn];
1836	if (temp[jColumn] >= 0) {
1837	int iObject = temp[jColumn];
1838	CbcSimpleInteger * obj =
1839	dynamic_cast <CbcSimpleInteger *>(originalObject[iObject]) ;
1840	if (!obj) {
1841	object_[numberObjects_] = originalObject[iObject]->clone();
1842	// redo ids etc
1843	CbcObject * obj =
1844	dynamic_cast <CbcObject *>(object_[numberObjects_]) ;
1845	assert (obj);
1846	obj->redoSequenceEtc(this, numberColumns, originalColumns);
1847	numberObjects_++;
1848	}
1849	}
1850	}
1851	// now append non column stuff
1852	for (iObject = 0; iObject < numberOriginalObjects; iObject++) {
1853	iColumn = originalObject[iObject]->columnNumber();
1854	if (iColumn < 0) {
1855	// already has column numbers changed
1856	object_[numberObjects_] = originalObject[iObject]->clone();
1857	#ifdef JJF_ZERO
1858	// redo ids etc
1859	CbcObject * obj =
1860	dynamic_cast <CbcObject *>(object_[numberObjects_]) ;
1861	assert (obj);
1862	obj->redoSequenceEtc(this, numberColumns, originalColumns);
1863	#endif
1864	numberObjects_++;
1865	}
1866	}
1867	delete [] temp;
1868	if (!numberObjects_)
1869	handler_->message(CBC_NOINT, messages_) << CoinMessageEol ;
1870	} else {
1871	int numberColumns = getNumCols();
1872	CglPreProcess * process = strategy_->process();
1873	assert (process);
1874	const int * originalColumns = process->originalColumns();
1875	// try and redo debugger
1876	OsiRowCutDebugger * debugger = const_cast<OsiRowCutDebugger *> (solver_->getRowCutDebuggerAlways());
1877	if (debugger)
1878	debugger->redoSolution(numberColumns, originalColumns);
1879	}
1880	} else {
1881	//no preprocessing
1882	originalSolver = NULL;
1883	}
1884	strategy_->setupCutGenerators(*this);
1885	strategy_->setupHeuristics(*this);
1886	// Set strategy print level to models
1887	strategy_->setupPrinting(*this, handler_->logLevel());
1888	}
1889	eventHappened_ = false;
1890	CbcEventHandler *eventHandler = getEventHandler() ;
1891	if (eventHandler)
1892	eventHandler->setModel(this);
1893	#define CLIQUE_ANALYSIS
1894	#ifdef CLIQUE_ANALYSIS
1895	// set up for probing
1896	// If we're doing clever stuff with cliques, additional info here.
1897	if (!parentModel_)
1898	probingInfo_ = new CglTreeProbingInfo(solver_);
1899	else
1900	probingInfo_ = NULL;
1901	#else
1902	probingInfo_ = NULL;
1903	#endif
1904
1905	// Try for dominated columns
1906	if ((specialOptions_&64) != 0) {
1907	CglDuplicateRow dupcuts(solver_);
1908	dupcuts.setMode(2);
1909	CglStored * storedCuts = dupcuts.outDuplicates(solver_);
1910	if (storedCuts) {
1911	COIN_DETAIL_PRINT(printf("adding dup cuts\n"));
1912	addCutGenerator(storedCuts, 1, "StoredCuts from dominated",
1913	true, false, false, -200);
1914	}
1915	}
1916	if (!nodeCompare_)
1917	nodeCompare_ = new CbcCompareDefault();;
1918	// See if hot start wanted
1919	CbcCompareBase * saveCompare = NULL;
1920	// User supplied hotstart. Adapt for preprocessing.
1921	if (hotstartSolution_) {
1922	if (strategy_ && strategy_->preProcessState() > 0) {
1923	CglPreProcess * process = strategy_->process();
1924	assert (process);
1925	int n = solver_->getNumCols();
1926	const int * originalColumns = process->originalColumns();
1927	// columns should be in order ... but
1928	double * tempS = new double[n];
1929	for (int i = 0; i < n; i++) {
1930	int iColumn = originalColumns[i];
1931	tempS[i] = hotstartSolution_[iColumn];
1932	}
1933	delete [] hotstartSolution_;
1934	hotstartSolution_ = tempS;
1935	if (hotstartPriorities_) {
1936	int * tempP = new int [n];
1937	for (int i = 0; i < n; i++) {
1938	int iColumn = originalColumns[i];
1939	tempP[i] = hotstartPriorities_[iColumn];
1940	}
1941	delete [] hotstartPriorities_;
1942	hotstartPriorities_ = tempP;
1943	}
1944	}
1945	saveCompare = nodeCompare_;
1946	// depth first
1947	nodeCompare_ = new CbcCompareDepth();
1948	}
1949	if (!problemFeasibility_)
1950	problemFeasibility_ = new CbcFeasibilityBase();
1951	# ifdef CBC_DEBUG
1952	std::string problemName ;
1953	solver_->getStrParam(OsiProbName, problemName) ;
1954	printf("Problem name - %s\n", problemName.c_str()) ;
1955	solver_->setHintParam(OsiDoReducePrint, false, OsiHintDo, 0) ;
1956	# endif
1957	/*
1958	Assume we're done, and see if we're proven wrong.
1959	*/
1960	status_ = 0 ;
1961	secondaryStatus_ = 0;
1962	phase_ = 0;
1963	/*
1964	Scan the variables, noting the integer variables. Create an
1965	CbcSimpleInteger object for each integer variable.
1966	*/
1967	findIntegers(false) ;
1968	// Say not dynamic pseudo costs
1969	ownership_ &= ~0x40000000;
1970	// If dynamic pseudo costs then do
1971	if (numberBeforeTrust_)
1972	convertToDynamic();
1973	// Set up char array to say if integer (speed)
1974	delete [] integerInfo_;
1975	{
1976	int n = solver_->getNumCols();
1977	integerInfo_ = new char [n];
1978	for (int i = 0; i < n; i++) {
1979	if (solver_->isInteger(i))
1980	integerInfo_[i] = 1;
1981	else
1982	integerInfo_[i] = 0;
1983	}
1984	}
1985	if (preferredWay_) {
1986	// set all unset ones
1987	for (int iObject = 0 ; iObject < numberObjects_ ; iObject++) {
1988	CbcObject * obj =
1989	dynamic_cast <CbcObject *>(object_[iObject]) ;
1990	if (obj && !obj->preferredWay())
1991	obj->setPreferredWay(preferredWay_);
1992	}
1993	}
1994	/*
1995	Ensure that objects on the lists of OsiObjects, heuristics, and cut
1996	generators attached to this model all refer to this model.
1997	*/
1998	synchronizeModel() ;
1999	if (!solverCharacteristics_) {
2000	OsiBabSolver * solverCharacteristics = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
2001	if (solverCharacteristics) {
2002	solverCharacteristics_ = solverCharacteristics;
2003	} else {
2004	// replace in solver
2005	OsiBabSolver defaultC;
2006	solver_->setAuxiliaryInfo(&defaultC);
2007	solverCharacteristics_ = dynamic_cast<OsiBabSolver *> (solver_->getAuxiliaryInfo());
2008	}
2009	}
2010
2011	solverCharacteristics_->setSolver(solver_);
2012	// Set so we can tell we are in initial phase in resolve
2013	continuousObjective_ = -COIN_DBL_MAX ;
2014	/*
2015	Solve the relaxation.
2016
2017	Apparently there are circumstances where this will be non-trivial --- i.e.,
2018	we've done something since initialSolve that's trashed the solution to the
2019	continuous relaxation.
2020	*/
2021	/* Tell solver we are in Branch and Cut
2022	Could use last parameter for subtle differences */
2023	solver_->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
2024	#ifdef COIN_HAS_CLP
2025	{
2026	OsiClpSolverInterface * clpSolver
2027	= dynamic_cast<OsiClpSolverInterface *> (solver_);
2028	if (clpSolver) {
2029	ClpSimplex * clpSimplex = clpSolver->getModelPtr();
2030	if ((specialOptions_&32) == 0) {
2031	// take off names (unless going to be saving)
2032	if (numberAnalyzeIterations_>=0\|\|(-numberAnalyzeIterations_&64)==0)
2033	clpSimplex->dropNames();
2034	}
2035	// no crunch if mostly continuous
2036	if ((clpSolver->specialOptions()&(1 + 8)) != (1 + 8)) {
2037	int numberColumns = solver_->getNumCols();
2038	if (numberColumns > 1000 && numberIntegers_*4 < numberColumns)
2039	clpSolver->setSpecialOptions(clpSolver->specialOptions()&(~1));
2040	}
2041	//#define NO_CRUNCH
2042	#ifdef NO_CRUNCH
2043	printf("TEMP switching off crunch\n");
2044	int iOpt = clpSolver->specialOptions();
2045	iOpt &= ~1;
2046	iOpt \|= 65536;
2047	clpSolver->setSpecialOptions(iOpt);
2048	#endif
2049	}
2050	}
2051	#endif
2052	bool feasible;
2053	numberSolves_ = 0 ;
2054	{
2055	// check
2056	int numberOdd = 0;
2057	for (int i = 0; i < numberObjects_; i++) {
2058	CbcSimpleInteger * obj =
2059	dynamic_cast <CbcSimpleInteger *>(object_[i]) ;
2060	if (!obj)
2061	numberOdd++;
2062	}
2063	if (numberOdd)
2064	moreSpecialOptions_ \|= 1073741824;
2065	}
2066	// If NLP then we assume already solved outside branchAndbound
2067	if (!solverCharacteristics_->solverType() \|\| solverCharacteristics_->solverType() == 4) {
2068	feasible = resolve(NULL, 0) != 0 ;
2069	} else {
2070	// pick up given status
2071	feasible = (solver_->isProvenOptimal() &&
2072	!solver_->isDualObjectiveLimitReached()) ;
2073	}
2074	if (problemFeasibility_->feasible(this, 0) < 0) {
2075	feasible = false; // pretend infeasible
2076	}
2077	numberSavedSolutions_ = 0;
2078	int saveNumberStrong = numberStrong_;
2079	int saveNumberBeforeTrust = numberBeforeTrust_;
2080	/*
2081	If the linear relaxation of the root is infeasible, bail out now. Otherwise,
2082	continue with processing the root node.
2083	*/
2084	if (!feasible) {
2085	status_ = 0 ;
2086	if (!solver_->isProvenDualInfeasible()) {
2087	handler_->message(CBC_INFEAS, messages_) << CoinMessageEol ;
2088	secondaryStatus_ = 1;
2089	} else {
2090	handler_->message(CBC_UNBOUNDED, messages_) << CoinMessageEol ;
2091	secondaryStatus_ = 7;
2092	}
2093	originalContinuousObjective_ = COIN_DBL_MAX;
2094	if (bestSolution_ &&
2095	((specialOptions_&8388608)==0\|\|(specialOptions_&2048)!=0)) {
2096	// best solution found by various heuristics - set solution
2097	char general[200];
2098	sprintf(general,"Solution of %g already found by heuristic",
2099	bestObjective_);
2100	messageHandler()->message(CBC_GENERAL,
2101	messages())
2102	<< general << CoinMessageEol ;
2103	setCutoff(1.0e50) ; // As best solution should be worse than cutoff
2104	// change cutoff as constraint if wanted
2105	if (cutoffRowNumber_>=0) {
2106	if (solver_->getNumRows()>cutoffRowNumber_)
2107	solver_->setRowUpper(cutoffRowNumber_,1.0e50);
2108	}
2109	// also in continuousSolver_
2110	if (continuousSolver_) {
2111	// Solvers know about direction
2112	double direction = solver_->getObjSense();
2113	continuousSolver_->setDblParam(OsiDualObjectiveLimit, 1.0e50*direction);
2114	} else {
2115	continuousSolver_ = solver_->clone();
2116	}
2117	phase_ = 5;
2118	double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
2119	if ((specialOptions_&4) == 0)
2120	bestObjective_ += 100.0 * increment + 1.0e-3; // only set if we are going to solve
2121	setBestSolution(CBC_END_SOLUTION, bestObjective_, bestSolution_, 1) ;
2122	continuousSolver_->resolve() ;
2123	if (!continuousSolver_->isProvenOptimal()) {
2124	continuousSolver_->messageHandler()->setLogLevel(2) ;
2125	continuousSolver_->initialSolve() ;
2126	}
2127	delete solver_ ;
2128	solverCharacteristics_=NULL;
2129	solver_ = continuousSolver_ ;
2130	setPointers(solver_);
2131	continuousSolver_ = NULL ;
2132	}
2133	solverCharacteristics_ = NULL;
2134	if (flipObjective)
2135	flipModel();
2136	return ;
2137	} else if (!numberObjects_) {
2138	// nothing to do
2139	// Undo preprocessing performed during BaB.
2140	if (strategy_ && strategy_->preProcessState() > 0) {
2141	// undo preprocessing
2142	CglPreProcess * process = strategy_->process();
2143	assert (process);
2144	int n = originalSolver->getNumCols();
2145	if (bestSolution_) {
2146	delete [] bestSolution_;
2147	bestSolution_ = new double [n];
2148	process->postProcess(*solver_);
2149	}
2150	strategy_->deletePreProcess();
2151	// Solution now back in originalSolver
2152	delete solver_;
2153	solver_ = originalSolver;
2154	if (bestSolution_) {
2155	bestObjective_ = solver_->getObjValue() * solver_->getObjSense();
2156	memcpy(bestSolution_, solver_->getColSolution(), n*sizeof(double));
2157	}
2158	// put back original objects if there were any
2159	if (originalObject) {
2160	int iColumn;
2161	assert (ownObjects_);
2162	for (iColumn = 0; iColumn < numberObjects_; iColumn++)
2163	delete object_[iColumn];
2164	delete [] object_;
2165	numberObjects_ = numberOriginalObjects;
2166	object_ = originalObject;
2167	delete [] integerVariable_;
2168	numberIntegers_ = 0;
2169	for (iColumn = 0; iColumn < n; iColumn++) {
2170	if (solver_->isInteger(iColumn))
2171	numberIntegers_++;
2172	}
2173	integerVariable_ = new int[numberIntegers_];
2174	numberIntegers_ = 0;
2175	for (iColumn = 0; iColumn < n; iColumn++) {
2176	if (solver_->isInteger(iColumn))
2177	integerVariable_[numberIntegers_++] = iColumn;
2178	}
2179	}
2180	}
2181	if (flipObjective)
2182	flipModel();
2183	solverCharacteristics_ = NULL;
2184	bestObjective_ = solver_->getObjValue() * solver_->getObjSense();
2185	int numberColumns = solver_->getNumCols();
2186	delete [] bestSolution_;
2187	bestSolution_ = new double[numberColumns];
2188	CoinCopyN(solver_->getColSolution(), numberColumns, bestSolution_);
2189	return ;
2190	}
2191	/*
2192	See if we're using the Osi side of the branching hierarchy. If so, either
2193	convert existing CbcObjects to OsiObjects, or generate them fresh. In the
2194	first case, CbcModel owns the objects on the object_ list. In the second
2195	case, the solver holds the objects and object_ simply points to the
2196	solver's list.
2197
2198	080417 The conversion code here (the block protected by `if (obj)') cannot
2199	possibly be correct. On the Osi side, descent is OsiObject -> OsiObject2 ->
2200	all other Osi object classes. On the Cbc side, it's OsiObject -> CbcObject
2201	-> all other Cbc object classes. It's structurally impossible for any Osi
2202	object to descend from CbcObject. The only thing I can see is that this is
2203	really dead code, and object detection is now handled from the Osi side.
2204	*/
2205	// Convert to Osi if wanted
2206	//OsiBranchingInformation * persistentInfo = NULL;
2207	if (branchingMethod_ && branchingMethod_->chooseMethod()) {
2208	//persistentInfo = new OsiBranchingInformation(solver_);
2209	if (numberOriginalObjects) {
2210	for (int iObject = 0 ; iObject < numberObjects_ ; iObject++) {
2211	CbcObject * obj =
2212	dynamic_cast <CbcObject *>(object_[iObject]) ;
2213	if (obj) {
2214	CbcSimpleInteger * obj2 =
2215	dynamic_cast <CbcSimpleInteger *>(obj) ;
2216	if (obj2) {
2217	// back to Osi land
2218	object_[iObject] = obj2->osiObject();
2219	delete obj;
2220	} else {
2221	OsiSimpleInteger * obj3 =
2222	dynamic_cast <OsiSimpleInteger *>(obj) ;
2223	if (!obj3) {
2224	OsiSOS * obj4 =
2225	dynamic_cast <OsiSOS *>(obj) ;
2226	if (!obj4) {
2227	CbcSOS * obj5 =
2228	dynamic_cast <CbcSOS *>(obj) ;
2229	if (obj5) {
2230	// back to Osi land
2231	object_[iObject] = obj5->osiObject(solver_);
2232	} else {
2233	printf("Code up CbcObject type in Osi land\n");
2234	abort();
2235	}
2236	}
2237	}
2238	}
2239	}
2240	}
2241	// and add to solver
2242	//if (!solver_->numberObjects()) {
2243	solver_->addObjects(numberObjects_, object_);
2244	//} else {
2245	//if (solver_->numberObjects()!=numberOriginalObjects) {
2246	//printf("should have trapped that solver has objects before\n");
2247	//abort();
2248	//}
2249	//}
2250	} else {
2251	/*
2252	As of 080104, findIntegersAndSOS is misleading --- the default OSI
2253	implementation finds only integers.
2254	*/
2255	// do from solver
2256	deleteObjects(false);
2257	solver_->findIntegersAndSOS(false);
2258	numberObjects_ = solver_->numberObjects();
2259	object_ = solver_->objects();
2260	ownObjects_ = false;
2261	}
2262	branchingMethod_->chooseMethod()->setSolver(solver_);
2263	}
2264	// take off heuristics if have to (some do not work with SOS, for example)
2265	// object should know what's safe.
2266	{
2267	int numberOdd = 0;
2268	int numberSOS = 0;
2269	for (int i = 0; i < numberObjects_; i++) {
2270	if (!object_[i]->canDoHeuristics())
2271	numberOdd++;
2272	CbcSOS * obj =
2273	dynamic_cast <CbcSOS *>(object_[i]) ;
2274	if (obj)
2275	numberSOS++;
2276	}
2277	if (numberOdd) {
2278	if (numberHeuristics_ && (specialOptions_&1024)==0 ) {
2279	int k = 0;
2280	for (int i = 0; i < numberHeuristics_; i++) {
2281	if (!heuristic_[i]->canDealWithOdd())
2282	delete heuristic_[i];
2283	else
2284	heuristic_[k++] = heuristic_[i];
2285	}
2286	if (!k) {
2287	delete [] heuristic_;
2288	heuristic_ = NULL;
2289	}
2290	numberHeuristics_ = k;
2291	handler_->message(CBC_HEURISTICS_OFF, messages_) << numberOdd << CoinMessageEol ;
2292	}
2293	// If odd switch off AddIntegers
2294	specialOptions_ &= ~65536;
2295	// switch off fast nodes for now
2296	fastNodeDepth_ = -1;
2297	moreSpecialOptions_ &= ~33554432; // no diving
2298	} else if (numberSOS) {
2299	specialOptions_ \|= 128; // say can do SOS in dynamic mode
2300	// switch off fast nodes for now
2301	fastNodeDepth_ = -1;
2302	moreSpecialOptions_ &= ~33554432; // no diving
2303	}
2304	if (numberThreads_ > 0) {
2305	/* switch off fast nodes for now
2306	Trouble is that by time mini bab finishes code is
2307	looking at a different node
2308	*/
2309	fastNodeDepth_ = -1;
2310	}
2311	}
2312	// Save objective (just so user can access it)
2313	originalContinuousObjective_ = solver_->getObjValue()* solver_->getObjSense();
2314	bestPossibleObjective_ = originalContinuousObjective_;
2315	sumChangeObjective1_ = 0.0;
2316	sumChangeObjective2_ = 0.0;
2317	/*
2318	OsiRowCutDebugger knows an optimal answer for a subset of MIP problems.
2319	Assuming it recognises the problem, when called upon it will check a cut to
2320	see if it cuts off the optimal answer.
2321	*/
2322	// If debugger exists set specialOptions_ bit
2323	if (solver_->getRowCutDebuggerAlways()) {
2324	specialOptions_ \|= 1;
2325	}
2326
2327	# ifdef CBC_DEBUG
2328	if ((specialOptions_&1) == 0)
2329	solver_->activateRowCutDebugger(problemName.c_str()) ;
2330	if (solver_->getRowCutDebuggerAlways())
2331	specialOptions_ \|= 1;
2332	# endif
2333
2334	/*
2335	Begin setup to process a feasible root node.
2336	*/
2337	bestObjective_ = CoinMin(bestObjective_, 1.0e50) ;
2338	if (!bestSolution_) {
2339	numberSolutions_ = 0 ;
2340	numberHeuristicSolutions_ = 0 ;
2341	}
2342	stateOfSearch_ = 0;
2343	// Everything is minimization
2344	{
2345	// needed to sync cutoffs
2346	double value ;
2347	solver_->getDblParam(OsiDualObjectiveLimit, value) ;
2348	dblParam_[CbcCurrentCutoff] = value * solver_->getObjSense();
2349	}
2350	double cutoff = getCutoff() ;
2351	double direction = solver_->getObjSense() ;
2352	dblParam_[CbcOptimizationDirection] = direction;
2353	if (cutoff < 1.0e20 && direction < 0.0)
2354	messageHandler()->message(CBC_CUTOFF_WARNING1,
2355	messages())
2356	<< cutoff << -cutoff << CoinMessageEol ;
2357	if (cutoff > bestObjective_)
2358	cutoff = bestObjective_ ;
2359	setCutoff(cutoff) ;
2360	/*
2361	We probably already have a current solution, but just in case ...
2362	*/
2363	int numberColumns = getNumCols() ;
2364	if (!currentSolution_)
2365	currentSolution_ = new double[numberColumns] ;
2366	testSolution_ = currentSolution_;
2367	/*
2368	Create a copy of the solver, thus capturing the original (root node)
2369	constraint system (aka the continuous system).
2370	*/
2371	delete continuousSolver_;
2372	continuousSolver_ = solver_->clone() ;
2373	#ifdef CONFLICT_CUTS
2374	if ((moreSpecialOptions_&4194304)!=0) {
2375	#ifdef COIN_HAS_CLP
2376	OsiClpSolverInterface * clpSolver
2377	= dynamic_cast<OsiClpSolverInterface *> (solver_);
2378	if (clpSolver) {
2379	int specialOptions=clpSolver->getModelPtr()->specialOptions();
2380	// 2097152 switches on rays in crunch
2381	if (!parentModel_)
2382	clpSolver->getModelPtr()->setSpecialOptions(specialOptions\|32\|2097152);
2383	else
2384	clpSolver->getModelPtr()->setSpecialOptions(specialOptions&~(32\|2097152));
2385	}
2386	}
2387	#endif
2388	#endif
2389	#ifdef COIN_HAS_NTY
2390	// maybe allow on fix and restart later
2391	if ((moreSpecialOptions2_&(128\|256))!=0&&!parentModel_) {
2392	symmetryInfo_ = new CbcSymmetry();
2393	symmetryInfo_->setupSymmetry(*continuousSolver_);
2394	int numberGenerators = symmetryInfo_->statsOrbits(this,0);
2395	if (!symmetryInfo_->numberUsefulOrbits()&&(moreSpecialOptions2_&(128\|256))!=(128\|256)) {
2396	delete symmetryInfo_;
2397	symmetryInfo_=NULL;
2398	moreSpecialOptions2_ &= ~(128\|256);
2399	}
2400	if ((moreSpecialOptions2_&(128\|256))==(128\|256)) {
2401	//moreSpecialOptions2_ &= ~256;
2402	}
2403	}
2404	#endif
2405
2406	// add cutoff as constraint if wanted
2407	if (cutoffRowNumber_==-2) {
2408	if (!parentModel_) {
2409	int numberColumns=solver_->getNumCols();
2410	double * obj = CoinCopyOfArray(solver_->getObjCoefficients(),numberColumns);
2411	int * indices = new int [numberColumns];
2412	int n=0;
2413	for (int i=0;i<numberColumns;i++) {
2414	if (obj[i]) {
2415	indices[n]=i;
2416	obj[n++]=obj[i];
2417	}
2418	}
2419	if (n) {
2420	double cutoff=getCutoff();
2421	// relax a little bit
2422	cutoff += 1.0e-4;
2423	double offset;
2424	solver_->getDblParam(OsiObjOffset, offset);
2425	cutoffRowNumber_ = solver_->getNumRows();
2426	solver_->addRow(n,indices,obj,-COIN_DBL_MAX,CoinMin(cutoff,1.0e25)+offset);
2427	} else {
2428	// no objective!
2429	cutoffRowNumber_ = -1;
2430	}
2431	delete [] indices;
2432	delete [] obj;
2433	} else {
2434	// switch off
2435	cutoffRowNumber_ = -1;
2436	}
2437	}
2438	numberRowsAtContinuous_ = getNumRows() ;
2439	solver_->saveBaseModel();
2440	/*
2441	Check the objective to see if we can deduce a nontrivial increment. If
2442	it's better than the current value for CbcCutoffIncrement, it'll be
2443	installed.
2444	*/
2445	if (solverCharacteristics_->reducedCostsAccurate())
2446	analyzeObjective() ;
2447	{
2448	// may be able to change cutoff now
2449	double cutoff = getCutoff();
2450	double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
2451	if (cutoff > bestObjective_ - increment) {
2452	cutoff = bestObjective_ - increment ;
2453	setCutoff(cutoff) ;
2454	}
2455	}
2456	#ifdef COIN_HAS_CLP
2457	// Possible save of pivot method
2458	ClpDualRowPivot * savePivotMethod = NULL;
2459	{
2460	// pass tolerance and increment to solver
2461	OsiClpSolverInterface * clpSolver
2462	= dynamic_cast<OsiClpSolverInterface *> (solver_);
2463	if (clpSolver)
2464	clpSolver->setStuff(getIntegerTolerance(), getCutoffIncrement());
2465	#ifdef CLP_RESOLVE
2466	if ((moreSpecialOptions_&1048576)!=0&&!parentModel_&&clpSolver) {
2467	resolveClp(clpSolver,0);
2468	}
2469	#endif
2470	}
2471	#endif
2472	/*
2473	Set up for cut generation. addedCuts_ holds the cuts which are relevant for
2474	the active subproblem. whichGenerator will be used to record the generator
2475	that produced a given cut.
2476	*/
2477	#define INITIAL_MAXIMUM_WHICH 1000
2478	maximumWhich_ = INITIAL_MAXIMUM_WHICH ;
2479	delete [] whichGenerator_;
2480	whichGenerator_ = new int[maximumWhich_] ;
2481	memset(whichGenerator_, 0, maximumWhich_*sizeof(int));
2482	maximumNumberCuts_ = 0 ;
2483	currentNumberCuts_ = 0 ;
2484	delete [] addedCuts_ ;
2485	addedCuts_ = NULL ;
2486	OsiObject ** saveObjects = NULL;
2487	maximumRows_ = numberRowsAtContinuous_;
2488	currentDepth_ = 0;
2489	workingBasis_.resize(maximumRows_, numberColumns);
2490	/*
2491	Set up an empty heap and associated data structures to hold the live set
2492	(problems which require further exploration).
2493	*/
2494	CbcCompareDefault * compareActual
2495	= dynamic_cast<CbcCompareDefault *> (nodeCompare_);
2496	if (compareActual) {
2497	compareActual->setBestPossible(direction*solver_->getObjValue());
2498	compareActual->setCutoff(getCutoff());
2499	#ifdef JJF_ZERO
2500	if (false && !numberThreads_ && !parentModel_) {
2501	printf("CbcTreeArray ? threads ? parentArray\n");
2502	// Setup new style tree
2503	delete tree_;
2504	tree_ = new CbcTreeArray();
2505	}
2506	#endif
2507	}
2508	tree_->setComparison(*nodeCompare_) ;
2509	/*
2510	Used to record the path from a node to the root of the search tree, so that
2511	we can then traverse from the root to the node when restoring a subproblem.
2512	*/
2513	maximumDepth_ = 10 ;
2514	delete [] walkback_ ;
2515	walkback_ = new CbcNodeInfo * [maximumDepth_] ;
2516	lastDepth_ = 0;
2517	delete [] lastNodeInfo_ ;
2518	lastNodeInfo_ = new CbcNodeInfo * [maximumDepth_] ;
2519	delete [] lastNumberCuts_ ;
2520	lastNumberCuts_ = new int [maximumDepth_] ;
2521	maximumCuts_ = 100;
2522	lastNumberCuts2_ = 0;
2523	delete [] lastCut_;
2524	lastCut_ = new const OsiRowCut * [maximumCuts_];
2525	/*
2526	Used to generate bound edits for CbcPartialNodeInfo.
2527	*/
2528	double * lowerBefore = new double [numberColumns] ;
2529	double * upperBefore = new double [numberColumns] ;
2530	/*
2531	Set up to run heuristics and generate cuts at the root node. The heavy
2532	lifting is hidden inside the calls to doHeuristicsAtRoot and solveWithCuts.
2533
2534	To start, tell cut generators they can be a bit more aggressive at the
2535	root node.
2536
2537	QUESTION: phase_ = 0 is documented as `initial solve', phase = 1 as `solve
2538	with cuts at root'. Is phase_ = 1 the correct indication when
2539	doHeurisiticsAtRoot is called to run heuristics outside of the main
2540	cut / heurisitc / reoptimise loop in solveWithCuts?
2541
2542	Generate cuts at the root node and reoptimise. solveWithCuts does the heavy
2543	lifting. It will iterate a generate/reoptimise loop (including reduced cost
2544	fixing) until no cuts are generated, the change in objective falls off, or
2545	the limit on the number of rounds of cut generation is exceeded.
2546
2547	At the end of all this, any cuts will be recorded in cuts and also
2548	installed in the solver's constraint system. We'll have reoptimised, and
2549	removed any slack cuts (numberOldActiveCuts_ and numberNewCuts_ have been
2550	adjusted accordingly).
2551
2552	Tell cut generators they can be a bit more aggressive at root node
2553
2554	TODO: Why don't we make a copy of the solution after solveWithCuts?
2555	TODO: If numberUnsatisfied == 0, don't we have a solution?
2556	*/
2557	phase_ = 1;
2558	int iCutGenerator;
2559	for (iCutGenerator = 0; iCutGenerator < numberCutGenerators_; iCutGenerator++) {
2560	// If parallel switch off global cuts
2561	if (numberThreads_) {
2562	generator_[iCutGenerator]->setGlobalCuts(false);
2563	generator_[iCutGenerator]->setGlobalCutsAtRoot(false);
2564	}
2565	CglCutGenerator * generator = generator_[iCutGenerator]->generator();
2566	generator->setAggressiveness(generator->getAggressiveness() + 100);
2567	if (!generator->canDoGlobalCuts())
2568	generator->setGlobalCuts(false);
2569	}
2570	OsiCuts cuts ;
2571	int anyAction = -1 ;
2572	numberOldActiveCuts_ = 0 ;
2573	numberNewCuts_ = 0 ;
2574	// Array to mark solution
2575	delete [] usedInSolution_;
2576	usedInSolution_ = new int[numberColumns];
2577	CoinZeroN(usedInSolution_, numberColumns);
2578	/*
2579	For printing totals and for CbcNode (numberNodes_)
2580	*/
2581	numberIterations_ = 0 ;
2582	numberNodes_ = 0 ;
2583	numberNodes2_ = 0 ;
2584	maximumStatistics_ = 0;
2585	maximumDepthActual_ = 0;
2586	numberDJFixed_ = 0.0;
2587	if (!parentModel_) {
2588	if ((specialOptions_&262144) != 0) {
2589	// create empty stored cuts
2590	//storedRowCuts_ = new CglStored(solver_->getNumCols());
2591	} else if ((specialOptions_&524288) != 0 && storedRowCuts_) {
2592	// tighten and set best solution
2593	// A) tight bounds on integer variables
2594	/*
2595	storedRowCuts_ are coming in from outside, probably for nonlinear.
2596	John was unsure about origin.
2597	*/
2598	const double * lower = solver_->getColLower();
2599	const double * upper = solver_->getColUpper();
2600	const double * tightLower = storedRowCuts_->tightLower();
2601	const double * tightUpper = storedRowCuts_->tightUpper();
2602	int nTightened = 0;
2603	for (int i = 0; i < numberIntegers_; i++) {
2604	int iColumn = integerVariable_[i];
2605	if (tightLower[iColumn] > lower[iColumn]) {
2606	nTightened++;
2607	solver_->setColLower(iColumn, tightLower[iColumn]);
2608	}
2609	if (tightUpper[iColumn] < upper[iColumn]) {
2610	nTightened++;
2611	solver_->setColUpper(iColumn, tightUpper[iColumn]);
2612	}
2613	}
2614	if (nTightened)
2615	COIN_DETAIL_PRINT(printf("%d tightened by alternate cuts\n", nTightened));
2616	if (storedRowCuts_->bestObjective() < bestObjective_) {
2617	// B) best solution
2618	double objValue = storedRowCuts_->bestObjective();
2619	setBestSolution(CBC_SOLUTION, objValue,
2620	storedRowCuts_->bestSolution()) ;
2621	// Do heuristics
2622	// Allow RINS
2623	for (int i = 0; i < numberHeuristics_; i++) {
2624	CbcHeuristicRINS * rins
2625	= dynamic_cast<CbcHeuristicRINS *> (heuristic_[i]);
2626	if (rins) {
2627	rins->setLastNode(-100);
2628	}
2629	}
2630	}
2631	}
2632	}
2633	#ifdef SWITCH_VARIABLES
2634	// see if any switching variables
2635	if (numberIntegers_<solver_->getNumCols())
2636	findSwitching();
2637	#endif
2638	/*
2639	Run heuristics at the root. This is the only opportunity to run FPump; it
2640	will be removed from the heuristics list by doHeuristicsAtRoot.
2641	*/
2642	// See if multiple runs wanted
2643	CbcModel ** rootModels=NULL;
2644	if (!parentModel_&&multipleRootTries_%100) {
2645	double rootTimeCpu=CoinCpuTime();
2646	double startTimeRoot=CoinGetTimeOfDay();
2647	int numberRootThreads=1;
2648	/* undocumented fine tuning
2649	aabbcc where cc is number of tries
2650	bb if nonzero is number of threads
2651	aa if nonzero just do heuristics
2652	*/
2653	int numberModels = multipleRootTries_%100;
2654	#ifdef CBC_THREAD
2655	numberRootThreads = (multipleRootTries_/100)%100;
2656	if (!numberRootThreads) {
2657	if (numberThreads_<2)
2658	numberRootThreads=numberModels;
2659	else
2660	numberRootThreads=CoinMin(numberThreads_,numberModels);
2661	}
2662	#endif
2663	int otherOptions = (multipleRootTries_/10000)%100;
2664	rootModels = new CbcModel * [numberModels];
2665	unsigned int newSeed = randomSeed_;
2666	if (newSeed==0) {
2667	double time = fabs(CoinGetTimeOfDay());
2668	while (time>=COIN_INT_MAX)
2669	time *= 0.5;
2670	newSeed = static_cast<unsigned int>(time);
2671	} else if (newSeed<0) {
2672	newSeed = 123456789;
2673	#ifdef COIN_HAS_CLP
2674	OsiClpSolverInterface * clpSolver
2675	= dynamic_cast<OsiClpSolverInterface *> (solver_);
2676	if (clpSolver) {
2677	newSeed += clpSolver->getModelPtr()->randomNumberGenerator()->getSeed();
2678	}
2679	#endif
2680	}
2681	CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis *> (solver_->getEmptyWarmStart());
2682	for (int i=0;i<numberModels;i++) {
2683	rootModels[i]=new CbcModel(*this);
2684	rootModels[i]->setNumberThreads(0);
2685	rootModels[i]->setMaximumNodes(otherOptions ? -1 : 0);
2686	rootModels[i]->setRandomSeed(newSeed+10000000*i);
2687	rootModels[i]->randomNumberGenerator()->setSeed(newSeed+50000000*i);
2688	rootModels[i]->setMultipleRootTries(0);
2689	#ifdef COIN_HAS_NTY
2690	rootModels[i]->zapSymmetry();
2691	rootModels[i]->moreSpecialOptions2_ &= ~(128\|256); // off nauty
2692	#endif
2693	// use seed
2694	rootModels[i]->setSpecialOptions(specialOptions_ \|(4194304\|8388608));
2695	rootModels[i]->setMoreSpecialOptions(moreSpecialOptions_ &
2696	(~(134217728\|4194304)));
2697	rootModels[i]->setMoreSpecialOptions2(moreSpecialOptions2_ &
2698	(~(128\|256)));
2699	rootModels[i]->solver_->setWarmStart(basis);
2700	#ifdef COIN_HAS_CLP
2701	OsiClpSolverInterface * clpSolver
2702	= dynamic_cast<OsiClpSolverInterface *> (rootModels[i]->solver_);
2703	#define NEW_RANDOM_BASIS
2704	#ifdef NEW_RANDOM_BASIS
2705	if (i==0)
2706	continue;
2707	#endif
2708	if (clpSolver) {
2709	ClpSimplex * simplex = clpSolver->getModelPtr();
2710	if (defaultHandler_)
2711	simplex->setDefaultMessageHandler();
2712	simplex->setRandomSeed(newSeed+20000000*i);
2713	simplex->allSlackBasis();
2714	int logLevel=simplex->logLevel();
2715	if (logLevel==1)
2716	simplex->setLogLevel(0);
2717	if (i!=0) {
2718	#ifdef NEW_RANDOM_BASIS
2719	int numberRows = simplex->numberRows();
2720	int throwOut=20;//2+numberRows/100;
2721	for (int iThrow=0;iThrow<throwOut;iThrow++) {
2722	double random = simplex->randomNumberGenerator()->randomDouble();
2723	int iStart=static_cast<int>(random*numberRows);
2724	for (int j=iStart;j<numberRows;j++) {
2725	if (simplex->getRowStatus(j)!=ClpSimplex::basic) {
2726	simplex->setRowStatus(j,ClpSimplex::basic);
2727	break;
2728	}
2729	}
2730	}
2731	clpSolver->setWarmStart(NULL);
2732	#else
2733	double random = simplex->randomNumberGenerator()->randomDouble();
2734	int bias = static_cast<int>(random*(numberIterations/4));
2735	simplex->setMaximumIterations(numberIterations/2+bias);
2736	simplex->primal();
2737	simplex->setMaximumIterations(COIN_INT_MAX);
2738	simplex->dual();
2739	#endif
2740	} else {
2741	#ifndef NEW_RANDOM_BASIS
2742	simplex->primal();
2743	#endif
2744	#endif
2745	}
2746	#ifdef NEW_RANDOM_BASIS
2747	simplex->setLogLevel(logLevel);
2748	clpSolver->setWarmStart(NULL);
2749	#endif
2750	}
2751	for (int j=0;j<numberHeuristics_;j++)
2752	rootModels[i]->heuristic_[j]->setSeed(rootModels[i]->heuristic_[j]->getSeed()+100000000*i);
2753	for (int j=0;j<numberCutGenerators_;j++)
2754	rootModels[i]->generator_[j]->generator()->refreshSolver(rootModels[i]->solver_);
2755	}
2756	delete basis;
2757	#ifdef CBC_THREAD
2758	if (numberRootThreads==1) {
2759	#endif
2760	for (int iModel=0;iModel<numberModels;iModel++) {
2761	doRootCbcThread(rootModels[iModel]);
2762	// see if solved at root node
2763	if (rootModels[iModel]->getMaximumNodes()) {
2764	feasible=false;
2765	break;
2766	}
2767	}
2768	#ifdef CBC_THREAD
2769	} else {
2770	Coin_pthread_t * threadId = new Coin_pthread_t [numberRootThreads];
2771	for (int kModel=0;kModel<numberModels;kModel+=numberRootThreads) {
2772	bool finished=false;
2773	for (int iModel=kModel;iModel<CoinMin(numberModels,kModel+numberRootThreads);iModel++) {
2774	pthread_create(&(threadId[iModel-kModel].thr), NULL,
2775	doRootCbcThread,
2776	rootModels[iModel]);
2777	}
2778	// wait
2779	for (int iModel=kModel;iModel<CoinMin(numberModels,kModel+numberRootThreads);iModel++) {
2780	pthread_join(threadId[iModel-kModel].thr, NULL);
2781	}
2782	// see if solved at root node
2783	for (int iModel=kModel;iModel<CoinMin(numberModels,kModel+numberRootThreads);iModel++) {
2784	if (rootModels[iModel]->getMaximumNodes())
2785	finished=true;
2786	}
2787	if (finished) {
2788	feasible=false;
2789	break;
2790	}
2791	}
2792	delete [] threadId;
2793	}
2794	#endif
2795	// sort solutions
2796	int * which = new int [numberModels];
2797	double * value = new double [numberModels];
2798	int numberSolutions=0;
2799	for (int iModel=0;iModel<numberModels;iModel++) {
2800	if (rootModels[iModel]->bestSolution()) {
2801	which[numberSolutions]=iModel;
2802	value[numberSolutions++]=
2803	-rootModels[iModel]->getMinimizationObjValue();
2804	}
2805	}
2806	char general[100];
2807	rootTimeCpu=CoinCpuTime()-rootTimeCpu;
2808	if (numberRootThreads==1)
2809	sprintf(general,"Multiple root solvers took a total of %.2f seconds\n",
2810	rootTimeCpu);
2811	else
2812	sprintf(general,"Multiple root solvers took a total of %.2f seconds (%.2f elapsed)\n",
2813	rootTimeCpu,CoinGetTimeOfDay()-startTimeRoot);
2814	messageHandler()->message(CBC_GENERAL,
2815	messages())
2816	<< general << CoinMessageEol ;
2817	CoinSort_2(value,value+numberSolutions,which);
2818	// to get name
2819	CbcHeuristicRINS dummyHeuristic;
2820	dummyHeuristic.setHeuristicName("Multiple root solvers");
2821	lastHeuristic_=&dummyHeuristic;
2822	for (int i=0;i<numberSolutions;i++) {
2823	double objValue = - value[i];
2824	if (objValue<getCutoff()) {
2825	int iModel=which[i];
2826	setBestSolution(CBC_ROUNDING,objValue,
2827	rootModels[iModel]->bestSolution());
2828	}
2829	}
2830	lastHeuristic_=NULL;
2831	delete [] which;
2832	delete [] value;
2833	}
2834	// Do heuristics
2835	if (numberObjects_&&!rootModels)
2836	doHeuristicsAtRoot();
2837	if (solverCharacteristics_->solutionAddsCuts()) {
2838	// With some heuristics solver needs a resolve here
2839	solver_->resolve();
2840	if(!isProvenOptimal()){
2841	solver_->initialSolve();
2842	}
2843	}
2844	/*
2845	Grepping through the code, it would appear that this is a command line
2846	debugging hook. There's no obvious place in the code where this is set to
2847	a negative value.
2848
2849	User hook, says John.
2850	*/
2851	if ( intParam_[CbcMaxNumNode] < 0
2852	\|\|numberSolutions_>=getMaximumSolutions())
2853	eventHappened_ = true; // stop as fast as possible
2854	stoppedOnGap_ = false ;
2855	// See if can stop on gap
2856	bestPossibleObjective_ = solver_->getObjValue() * solver_->getObjSense();
2857	if(canStopOnGap()) {
2858	if (bestPossibleObjective_ < getCutoff())
2859	stoppedOnGap_ = true ;
2860	feasible = false;
2861	//eventHappened_=true; // stop as fast as possible
2862	}
2863	/*
2864	Set up for statistics collection, if requested. Standard values are
2865	documented in CbcModel.hpp. The magic number 100 will trigger a dump of
2866	CbcSimpleIntegerDynamicPseudoCost objects (no others). Looks like another
2867	command line debugging hook.
2868	*/
2869	statistics_ = NULL;
2870	// Do on switch
2871	if (doStatistics > 0 && doStatistics <= 100) {
2872	maximumStatistics_ = 10000;
2873	statistics_ = new CbcStatistics * [maximumStatistics_];
2874	memset(statistics_, 0, maximumStatistics_sizeof(CbcStatistics ));
2875	}
2876	// See if we can add integers
2877	if (noObjects && numberIntegers_ < solver_->getNumCols() && (specialOptions_&65536) != 0 && !parentModel_ && false) {
2878	int numberIntegers1=0;
2879	int numberColumns = solver_->getNumCols();
2880	for (int i=0;i<numberColumns;i++) {
2881	if (solver_->isInteger(i))
2882	numberIntegers1++;
2883	}
2884	AddIntegers();
2885	// make sure in sync
2886	int numberIntegers2=0;
2887	for (int i=0;i<numberColumns;i++) {
2888	if (solver_->isInteger(i))
2889	numberIntegers2++;
2890	}
2891	if (numberIntegers1<numberIntegers2) {
2892	findIntegers(true,2);
2893	convertToDynamic();
2894	}
2895	}
2896
2897	/*
2898	Do an initial round of cut generation for the root node. Depending on the
2899	type of underlying solver, we may want to do this even if the initial query
2900	to the objects indicates they're satisfied.
2901
2902	solveWithCuts does the heavy lifting. It will iterate a generate/reoptimise
2903	loop (including reduced cost fixing) until no cuts are generated, the
2904	change in objective falls off, or the limit on the number of rounds of cut
2905	generation is exceeded.
2906
2907	At the end of all this, any cuts will be recorded in cuts and also
2908	installed in the solver's constraint system. We'll have reoptimised, and
2909	removed any slack cuts (numberOldActiveCuts_ and numberNewCuts_ have been
2910	adjusted accordingly).
2911	*/
2912	int iObject ;
2913	int numberUnsatisfied = 0 ;
2914	delete [] currentSolution_;
2915	currentSolution_ = new double [numberColumns];
2916	testSolution_ = currentSolution_;
2917	memcpy(currentSolution_, solver_->getColSolution(),
2918	numberColumns*sizeof(double)) ;
2919	// point to useful information
2920	OsiBranchingInformation usefulInfo = usefulInformation();
2921
2922	for (iObject = 0 ; iObject < numberObjects_ ; iObject++) {
2923	double infeasibility =
2924	object_[iObject]->checkInfeasibility(&usefulInfo) ;
2925	if (infeasibility ) numberUnsatisfied++ ;
2926	}
2927	// replace solverType
2928	double * tightBounds = NULL;
2929	if (solverCharacteristics_->tryCuts()) {
2930
2931	if (numberUnsatisfied) {
2932	// User event
2933	if (!eventHappened_ && feasible) {
2934	if (rootModels) {
2935	// for fixings
2936	int numberColumns=solver_->getNumCols();
2937	tightBounds = new double [2*numberColumns];
2938	{
2939	const double * lower = solver_->getColLower();
2940	const double * upper = solver_->getColUpper();
2941	for (int i=0;i<numberColumns;i++) {
2942	tightBounds[2*i+0]=lower[i];
2943	tightBounds[2*i+1]=upper[i];
2944	}
2945	}
2946	int numberModels = multipleRootTries_%100;
2947	const OsiSolverInterface ** solvers = new
2948	const OsiSolverInterface * [numberModels];
2949	int numberRows=continuousSolver_->getNumRows();
2950	int maxCuts=0;
2951	for (int i=0;i<numberModels;i++) {
2952	solvers[i]=rootModels[i]->solver();
2953	const double * lower = solvers[i]->getColLower();
2954	const double * upper = solvers[i]->getColUpper();
2955	for (int j=0;j<numberColumns;j++) {
2956	tightBounds[2j+0]=CoinMax(lower[j],tightBounds[2j+0]);
2957	tightBounds[2j+1]=CoinMin(upper[j],tightBounds[2j+1]);
2958	}
2959	int numberRows2=solvers[i]->getNumRows();
2960	assert (numberRows2>=numberRows);
2961	maxCuts += numberRows2-numberRows;
2962	// accumulate statistics
2963	for (int j=0;j<numberCutGenerators_;j++) {
2964	generator_[j]->addStatistics(rootModels[i]->cutGenerator(j));
2965	}
2966	}
2967	for (int j=0;j<numberCutGenerators_;j++) {
2968	generator_[j]->scaleBackStatistics(numberModels);
2969	}
2970	//CbcRowCuts rowCut(maxCuts);
2971	const OsiRowCutDebugger *debugger = NULL;
2972	if ((specialOptions_&1) != 0)
2973	debugger = solver_->getRowCutDebugger() ;
2974	for (int iModel=0;iModel<numberModels;iModel++) {
2975	int numberRows2=solvers[iModel]->getNumRows();
2976	const CoinPackedMatrix * rowCopy = solvers[iModel]->getMatrixByRow();
2977	const int * rowLength = rowCopy->getVectorLengths();
2978	const double * elements = rowCopy->getElements();
2979	const int * column = rowCopy->getIndices();
2980	const CoinBigIndex * rowStart = rowCopy->getVectorStarts();
2981	const double * rowLower = solvers[iModel]->getRowLower();
2982	const double * rowUpper = solvers[iModel]->getRowUpper();
2983	for (int iRow=numberRows;iRow<numberRows2;iRow++) {
2984	OsiRowCut rc;
2985	rc.setLb(rowLower[iRow]);
2986	rc.setUb(rowUpper[iRow]);
2987	CoinBigIndex start = rowStart[iRow];
2988	rc.setRow(rowLength[iRow],column+start,elements+start,false);
2989	if (debugger)
2990	CoinAssert (!debugger->invalidCut(rc));
2991	globalCuts_.addCutIfNotDuplicate(rc);
2992	}
2993	//int cutsAdded=globalCuts_.numberCuts()-numberCuts;
2994	//numberCuts += cutsAdded;
2995	//printf("Model %d gave %d cuts (out of %d possible)\n",
2996	// iModel,cutsAdded,numberRows2-numberRows);
2997	}
2998	// normally replace global cuts
2999	//if (!globalCuts_.())
3000	//globalCuts_=rowCutrowCut.addCuts(globalCuts_);
3001	//rowCut.addCuts(globalCuts_);
3002	int nTightened=0;
3003	assert(feasible);
3004	{
3005	double tolerance=1.0e-5;
3006	const double * lower = solver_->getColLower();
3007	const double * upper = solver_->getColUpper();
3008	for (int i=0;i<numberColumns;i++) {
3009	if (tightBounds[2i+0]>tightBounds[2i+1]+1.0e-9) {
3010	feasible=false;
3011	char general[200];
3012	sprintf(general,"Solvers give infeasible bounds on %d %g,%g was %g,%g - search finished\n",
3013	i,tightBounds[2i+0],tightBounds[2i+1],lower[i],upper[i]);
3014	messageHandler()->message(CBC_GENERAL,messages())
3015	<< general << CoinMessageEol ;
3016	break;
3017	}
3018	double oldLower=lower[i];
3019	double oldUpper=upper[i];
3020	if (tightBounds[2*i+0]>oldLower+tolerance) {
3021	nTightened++;
3022	solver_->setColLower(i,tightBounds[2*i+0]);
3023	}
3024	if (tightBounds[2*i+1]<oldUpper-tolerance) {
3025	nTightened++;
3026	solver_->setColUpper(i,tightBounds[2*i+1]);
3027	}
3028	}
3029	}
3030	delete [] tightBounds;
3031	tightBounds=NULL;
3032	char printBuffer[200];
3033	sprintf(printBuffer,"%d solvers added %d different cuts out of pool of %d",
3034	numberModels,globalCuts_.sizeRowCuts(),maxCuts);
3035	messageHandler()->message(CBC_GENERAL, messages())
3036	<< printBuffer << CoinMessageEol ;
3037	if (nTightened) {
3038	sprintf(printBuffer,"%d bounds were tightened",
3039	nTightened);
3040	messageHandler()->message(CBC_GENERAL, messages())
3041	<< printBuffer << CoinMessageEol ;
3042	}
3043	delete [] solvers;
3044	}
3045	if (!parentModel_&&(moreSpecialOptions_&67108864) != 0) {
3046	// load cuts from file
3047	FILE * fp = fopen("global.cuts","rb");
3048	if (fp) {
3049	size_t nRead;
3050	int numberColumns=solver_->getNumCols();
3051	int numCols;
3052	nRead = fread(&numCols, sizeof(int), 1, fp);
3053	if (nRead != 1)
3054	throw("Error in fread");
3055	if (numberColumns!=numCols) {
3056	printf("Mismatch on columns %d %d\n",numberColumns,numCols);
3057	fclose(fp);
3058	} else {
3059	// If rootModel just do some
3060	double threshold=-1.0;
3061	if (!multipleRootTries_)
3062	threshold=0.5;
3063	int initialCuts=0;
3064	int initialGlobal = globalCuts_.sizeRowCuts();
3065	double * elements = new double [numberColumns+2];
3066	int * indices = new int [numberColumns];
3067	int numberEntries=1;
3068	while (numberEntries>0) {
3069	nRead = fread(&numberEntries, sizeof(int), 1, fp);
3070	if (nRead != 1)
3071	throw("Error in fread");
3072	double randomNumber=randomNumberGenerator_.randomDouble();
3073	if (numberEntries>0) {
3074	initialCuts++;
3075	nRead = fread(elements, sizeof(double), numberEntries+2, fp);
3076	if (nRead != static_cast<size_t>(numberEntries+2))
3077	throw("Error in fread");
3078	nRead = fread(indices, sizeof(int), numberEntries, fp);
3079	if (nRead != static_cast<size_t>(numberEntries))
3080	throw("Error in fread");
3081	if (randomNumber>threshold) {
3082	OsiRowCut rc;
3083	rc.setLb(elements[numberEntries]);
3084	rc.setUb(elements[numberEntries+1]);
3085	rc.setRow(numberEntries,indices,elements,
3086	false);
3087	rc.setGloballyValidAsInteger(2);
3088	globalCuts_.addCutIfNotDuplicate(rc) ;
3089	}
3090	}
3091	}
3092	fclose(fp);
3093	// fixes
3094	int nTightened=0;
3095	fp = fopen("global.fix","rb");
3096	if (fp) {
3097	nRead = fread(indices, sizeof(int), 2, fp);
3098	if (nRead != 2)
3099	throw("Error in fread");
3100	if (numberColumns!=indices[0]) {
3101	printf("Mismatch on columns %d %d\n",numberColumns,
3102	indices[0]);
3103	} else {
3104	indices[0]=1;
3105	while (indices[0]>=0) {
3106	nRead = fread(indices, sizeof(int), 2, fp);
3107	if (nRead != 2)
3108	throw("Error in fread");
3109	int iColumn=indices[0];
3110	if (iColumn>=0) {
3111	nTightened++;
3112	nRead = fread(elements, sizeof(double), 4, fp);
3113	if (nRead != 4)
3114	throw("Error in fread");
3115	solver_->setColLower(iColumn,elements[0]);
3116	solver_->setColUpper(iColumn,elements[1]);
3117	}
3118	}
3119	}
3120	}
3121	if (fp)
3122	fclose(fp);
3123	char printBuffer[200];
3124	sprintf(printBuffer,"%d cuts read in of which %d were unique, %d bounds tightened",
3125	initialCuts,
3126	globalCuts_.sizeRowCuts()-initialGlobal,nTightened);
3127	messageHandler()->message(CBC_GENERAL, messages())
3128	<< printBuffer << CoinMessageEol ;
3129	delete [] elements;
3130	delete [] indices;
3131	}
3132	}
3133	}
3134	if (feasible)
3135	feasible = solveWithCuts(cuts, maximumCutPassesAtRoot_,
3136	NULL);
3137	if (multipleRootTries_&&
3138	(moreSpecialOptions_&134217728)!=0) {
3139	FILE * fp=NULL;
3140	size_t nRead;
3141	int numberColumns=solver_->getNumCols();
3142	int initialCuts=0;
3143	if ((moreSpecialOptions_&134217728)!=0) {
3144	// append so go down to end
3145	fp = fopen("global.cuts","r+b");
3146	if (fp) {
3147	int numCols;
3148	nRead = fread(&numCols, sizeof(int), 1, fp);
3149	if (nRead != 1)
3150	throw("Error in fread");
3151	if (numberColumns!=numCols) {
3152	printf("Mismatch on columns %d %d\n",numberColumns,numCols);
3153	fclose(fp);
3154	fp=NULL;
3155	}
3156	}
3157	}
3158	double * elements = new double [numberColumns+2];
3159	int * indices = new int [numberColumns];
3160	if (fp) {
3161	int numberEntries=1;
3162	while (numberEntries>0) {
3163	fpos_t position;
3164	fgetpos(fp, &position);
3165	nRead = fread(&numberEntries, sizeof(int), 1, fp);
3166	if (nRead != 1)
3167	throw("Error in fread");
3168	if (numberEntries>0) {
3169	initialCuts++;
3170	nRead = fread(elements, sizeof(double), numberEntries+2, fp);
3171	if (nRead != static_cast<size_t>(numberEntries+2))
3172	throw("Error in fread");
3173	nRead = fread(indices, sizeof(int), numberEntries, fp);
3174	if (nRead != static_cast<size_t>(numberEntries))
3175	throw("Error in fread");
3176	} else {
3177	// end
3178	fsetpos(fp, &position);
3179	}
3180	}
3181	} else {
3182	fp = fopen("global.cuts","wb");
3183	size_t nWrite;
3184	nWrite=fwrite(&numberColumns,sizeof(int),1,fp);
3185	if (nWrite != 1)
3186	throw("Error in fwrite");
3187	}
3188	size_t nWrite;
3189	// now append binding cuts
3190	int numberC=continuousSolver_->getNumRows();
3191	int numberRows=solver_->getNumRows();
3192	printf("Saving %d cuts (up from %d)\n",
3193	initialCuts+numberRows-numberC,initialCuts);
3194	const double * rowLower = solver_->getRowLower();
3195	const double * rowUpper = solver_->getRowUpper();
3196	// Row copy
3197	CoinPackedMatrix matrixByRow(*solver_->getMatrixByRow());
3198	const double * elementByRow = matrixByRow.getElements();
3199	const int * column = matrixByRow.getIndices();
3200	const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
3201	const int * rowLength = matrixByRow.getVectorLengths();
3202	for (int iRow=numberC;iRow<numberRows;iRow++) {
3203	int n=rowLength[iRow];
3204	assert (n);
3205	CoinBigIndex start=rowStart[iRow];
3206	memcpy(elements,elementByRow+start,n*sizeof(double));
3207	memcpy(indices,column+start,n*sizeof(int));
3208	elements[n]=rowLower[iRow];
3209	elements[n+1]=rowUpper[iRow];
3210	nWrite=fwrite(&n,sizeof(int),1,fp);
3211	if (nWrite != 1)
3212	throw("Error in fwrite");
3213	nWrite=fwrite(elements,sizeof(double),n+2,fp);
3214	if (nWrite != static_cast<size_t>(n+2))
3215	throw("Error in fwrite");
3216	nWrite=fwrite(indices,sizeof(int),n,fp);
3217	if (nWrite != static_cast<size_t>(n))
3218	throw("Error in fwrite");
3219	}
3220	// eof marker
3221	int eofMarker=-1;
3222	nWrite=fwrite(&eofMarker,sizeof(int),1,fp);
3223	if (nWrite != 1)
3224	throw("Error in fwrite");
3225	fclose(fp);
3226	// do tighter bounds (? later extra to original columns)
3227	int nTightened=0;
3228	const double * lower = solver_->getColLower();
3229	const double * upper = solver_->getColUpper();
3230	const double * originalLower = continuousSolver_->getColLower();
3231	const double * originalUpper = continuousSolver_->getColUpper();
3232	double tolerance=1.0e-5;
3233	for (int i=0;i<numberColumns;i++) {
3234	if (lower[i]>originalLower[i]+tolerance) {
3235	nTightened++;
3236	}
3237	if (upper[i]<originalUpper[i]-tolerance) {
3238	nTightened++;
3239	}
3240	}
3241	if (nTightened) {
3242	fp = fopen("global.fix","wb");
3243	size_t nWrite;
3244	indices[0]=numberColumns;
3245	if (originalColumns_)
3246	indices[1]=COIN_INT_MAX;
3247	else
3248	indices[1]=-1;
3249	nWrite=fwrite(indices,sizeof(int),2,fp);
3250	if (nWrite != 2)
3251	throw("Error in fwrite");
3252	for (int i=0;i<numberColumns;i++) {
3253	int nTightened=0;
3254	if (lower[i]>originalLower[i]+tolerance) {
3255	nTightened++;
3256	}
3257	if (upper[i]<originalUpper[i]-tolerance) {
3258	nTightened++;
3259	}
3260	if (nTightened) {
3261	indices[0]=i;
3262	if (originalColumns_)
3263	indices[1]=originalColumns_[i];
3264	elements[0]=lower[i];
3265	elements[1]=upper[i];
3266	elements[2]=originalLower[i];
3267	elements[3]=originalUpper[i];
3268	nWrite=fwrite(indices,sizeof(int),2,fp);
3269	if (nWrite != 2)
3270	throw("Error in fwrite");
3271	nWrite=fwrite(elements,sizeof(double),4,fp);
3272	if (nWrite != 4)
3273	throw("Error in fwrite");
3274	}
3275	}
3276	// eof marker
3277	indices[0]=-1;
3278	nWrite=fwrite(indices,sizeof(int),2,fp);
3279	if (nWrite != 2)
3280	throw("Error in fwrite");
3281	fclose(fp);
3282	}
3283	delete [] elements;
3284	delete [] indices;
3285	}
3286	if ((specialOptions_&524288) != 0 && !parentModel_
3287	&& storedRowCuts_) {
3288	if (feasible) {
3289	/* pick up stuff and try again
3290	add cuts, maybe keep around
3291	do best solution and if so new heuristics
3292	obviously tighten bounds
3293	*/
3294	// A and B probably done on entry
3295	// A) tight bounds on integer variables
3296	const double * lower = solver_->getColLower();
3297	const double * upper = solver_->getColUpper();
3298	const double * tightLower = storedRowCuts_->tightLower();
3299	const double * tightUpper = storedRowCuts_->tightUpper();
3300	int nTightened = 0;
3301	for (int i = 0; i < numberIntegers_; i++) {
3302	int iColumn = integerVariable_[i];
3303	if (tightLower[iColumn] > lower[iColumn]) {
3304	nTightened++;
3305	solver_->setColLower(iColumn, tightLower[iColumn]);
3306	}
3307	if (tightUpper[iColumn] < upper[iColumn]) {
3308	nTightened++;
3309	solver_->setColUpper(iColumn, tightUpper[iColumn]);
3310	}
3311	}
3312	if (nTightened)
3313	COIN_DETAIL_PRINT(printf("%d tightened by alternate cuts\n", nTightened));
3314	if (storedRowCuts_->bestObjective() < bestObjective_) {
3315	// B) best solution
3316	double objValue = storedRowCuts_->bestObjective();
3317	setBestSolution(CBC_SOLUTION, objValue,
3318	storedRowCuts_->bestSolution()) ;
3319	// Do heuristics
3320	// Allow RINS
3321	for (int i = 0; i < numberHeuristics_; i++) {
3322	CbcHeuristicRINS * rins
3323	= dynamic_cast<CbcHeuristicRINS *> (heuristic_[i]);
3324	if (rins) {
3325	rins->setLastNode(-100);
3326	}
3327	}
3328	doHeuristicsAtRoot();
3329	}
3330	#ifdef JJF_ZERO
3331	int nCuts = storedRowCuts_->sizeRowCuts();
3332	// add to global list
3333	for (int i = 0; i < nCuts; i++) {
3334	OsiRowCut newCut(*storedRowCuts_->rowCutPointer(i));
3335	newCut.setGloballyValidAsInteger(2);
3336	newCut.mutableRow().setTestForDuplicateIndex(false);
3337	globalCuts_.insert(newCut) ;
3338	}
3339	#else
3340	addCutGenerator(storedRowCuts_, -99, "Stored from previous run",
3341	true, false, false, -200);
3342	#endif
3343	// Set cuts as active
3344	delete [] addedCuts_ ;
3345	maximumNumberCuts_ = cuts.sizeRowCuts();
3346	if (maximumNumberCuts_) {
3347	addedCuts_ = new CbcCountRowCut * [maximumNumberCuts_];
3348	} else {
3349	addedCuts_ = NULL;
3350	}
3351	for (int i = 0; i < maximumNumberCuts_; i++)
3352	addedCuts_[i] = new CbcCountRowCut(*cuts.rowCutPtr(i),
3353	NULL, -1, -1, 2);
3354	COIN_DETAIL_PRINT(printf("size %d\n", cuts.sizeRowCuts()));
3355	cuts = OsiCuts();
3356	currentNumberCuts_ = maximumNumberCuts_;
3357	feasible = solveWithCuts(cuts, maximumCutPassesAtRoot_,
3358	NULL);
3359	for (int i = 0; i < maximumNumberCuts_; i++)
3360	delete addedCuts_[i];
3361	}
3362	delete storedRowCuts_;
3363	storedRowCuts_ = NULL;
3364	}
3365	} else {
3366	feasible = false;
3367	}
3368	} else if (solverCharacteristics_->solutionAddsCuts() \|\|
3369	solverCharacteristics_->alwaysTryCutsAtRootNode()) {
3370	// may generate cuts and turn the solution
3371	//to an infeasible one
3372	feasible = solveWithCuts(cuts, 2,
3373	NULL);
3374	}
3375	}
3376	if (rootModels) {
3377	int numberModels = multipleRootTries_%100;
3378	for (int i=0;i<numberModels;i++)
3379	delete rootModels[i];
3380	delete [] rootModels;
3381	}
3382	// check extra info on feasibility
3383	if (!solverCharacteristics_->mipFeasible())
3384	feasible = false;
3385	// If max nodes==0 - don't do strong branching
3386	if (!getMaximumNodes()) {
3387	if (feasible)
3388	feasible=false;
3389	else
3390	setMaximumNodes(1); //allow to stop on success
3391	}
3392	topOfTree_=NULL;
3393	#ifdef CLP_RESOLVE
3394	if ((moreSpecialOptions_&2097152)!=0&&!parentModel_&&feasible) {
3395	OsiClpSolverInterface * clpSolver
3396	= dynamic_cast<OsiClpSolverInterface *> (solver_);
3397	if (clpSolver)
3398	resolveClp(clpSolver,0);
3399	}
3400	#endif
3401	// make cut generators less aggressive
3402	for (iCutGenerator = 0; iCutGenerator < numberCutGenerators_; iCutGenerator++) {
3403	CglCutGenerator * generator = generator_[iCutGenerator]->generator();
3404	generator->setAggressiveness(generator->getAggressiveness() - 100);
3405	}
3406	currentNumberCuts_ = numberNewCuts_ ;
3407	if (solverCharacteristics_->solutionAddsCuts()) {
3408	// With some heuristics solver needs a resolve here (don't know if this is bug in heuristics)
3409	solver_->resolve();
3410	if(!isProvenOptimal()){
3411	solver_->initialSolve();
3412	}
3413	}
3414	// See if can stop on gap
3415	bestPossibleObjective_ = solver_->getObjValue() * solver_->getObjSense();
3416	if(canStopOnGap()) {
3417	if (bestPossibleObjective_ < getCutoff())
3418	stoppedOnGap_ = true ;
3419	feasible = false;
3420	}
3421	// User event
3422	if (eventHappened_)
3423	feasible = false;
3424	#if defined(COIN_HAS_CLP)&&defined(COIN_HAS_CPX)
3425	/*
3426	This is the notion of using Cbc stuff to get going, then calling cplex to
3427	finish off.
3428	*/
3429	if (feasible && (specialOptions_&16384) != 0 && fastNodeDepth_ == -2 && !parentModel_) {
3430	// Use Cplex to do search!
3431	double time1 = CoinCpuTime();
3432	OsiClpSolverInterface * clpSolver
3433	= dynamic_cast<OsiClpSolverInterface *> (solver_);
3434	OsiCpxSolverInterface cpxSolver;
3435	double direction = clpSolver->getObjSense();
3436	cpxSolver.setObjSense(direction);
3437	// load up cplex
3438	const CoinPackedMatrix * matrix = continuousSolver_->getMatrixByCol();
3439	const double * rowLower = continuousSolver_->getRowLower();
3440	const double * rowUpper = continuousSolver_->getRowUpper();
3441	const double * columnLower = continuousSolver_->getColLower();
3442	const double * columnUpper = continuousSolver_->getColUpper();
3443	const double * objective = continuousSolver_->getObjCoefficients();
3444	cpxSolver.loadProblem(*matrix, columnLower, columnUpper,
3445	objective, rowLower, rowUpper);
3446	double * setSol = new double [numberIntegers_];
3447	int * setVar = new int [numberIntegers_];
3448	// cplex doesn't know about objective offset
3449	double offset = clpSolver->getModelPtr()->objectiveOffset();
3450	for (int i = 0; i < numberIntegers_; i++) {
3451	int iColumn = integerVariable_[i];
3452	cpxSolver.setInteger(iColumn);
3453	if (bestSolution_) {
3454	setSol[i] = bestSolution_[iColumn];
3455	setVar[i] = iColumn;
3456	}
3457	}
3458	CPXENVptr env = cpxSolver.getEnvironmentPtr();
3459	CPXLPptr lpPtr = cpxSolver.getLpPtr(OsiCpxSolverInterface::KEEPCACHED_ALL);
3460	cpxSolver.switchToMIP();
3461	if (bestSolution_) {
3462	#if 0
3463	CPXcopymipstart(env, lpPtr, numberIntegers_, setVar, setSol);
3464	#else
3465	int zero = 0;
3466	CPXaddmipstarts(env, lpPtr, 1, numberIntegers_, &zero, setVar, setSol, NULL, NULL);
3467	#endif
3468	}
3469	if (clpSolver->getNumRows() > continuousSolver_->getNumRows() && false) {
3470	// add cuts
3471	const CoinPackedMatrix * matrix = clpSolver->getMatrixByRow();
3472	const double * rhs = clpSolver->getRightHandSide();
3473	const char * rowSense = clpSolver->getRowSense();
3474	const double * elementByRow = matrix->getElements();
3475	const int * column = matrix->getIndices();
3476	const CoinBigIndex * rowStart = matrix->getVectorStarts();
3477	const int * rowLength = matrix->getVectorLengths();
3478	int nStart = continuousSolver_->getNumRows();
3479	int nRows = clpSolver->getNumRows();
3480	int size = rowStart[nRows-1] + rowLength[nRows-1] -
3481	rowStart[nStart];
3482	int nAdd = 0;
3483	double * rmatval = new double [size];
3484	int * rmatind = new int [size];
3485	int * rmatbeg = new int [nRows-nStart+1];
3486	size = 0;
3487	rmatbeg[0] = 0;
3488	for (int i = nStart; i < nRows; i++) {
3489	for (int k = rowStart[i]; k < rowStart[i] + rowLength[i]; k++) {
3490	rmatind[size] = column[k];
3491	rmatval[size++] = elementByRow[k];
3492	}
3493	nAdd++;
3494	rmatbeg[nAdd] = size;
3495	}
3496	CPXaddlazyconstraints(env, lpPtr, nAdd, size,
3497	rhs, rowSense, rmatbeg,
3498	rmatind, rmatval, NULL);
3499	CPXsetintparam( env, CPX_PARAM_REDUCE,
3500	// CPX_PREREDUCE_NOPRIMALORDUAL (0)
3501	CPX_PREREDUCE_PRIMALONLY);
3502	}
3503	if (getCutoff() < 1.0e50) {
3504	double useCutoff = getCutoff() + offset;
3505	if (bestObjective_ < 1.0e50)
3506	useCutoff = bestObjective_ + offset + 1.0e-7;
3507	cpxSolver.setDblParam(OsiDualObjectiveLimit, useCutoff*
3508	direction);
3509	if ( direction > 0.0 )
3510	CPXsetdblparam( env, CPX_PARAM_CUTUP, useCutoff ) ; // min
3511	else
3512	CPXsetdblparam( env, CPX_PARAM_CUTLO, useCutoff ) ; // max
3513	}
3514	CPXsetdblparam(env, CPX_PARAM_EPGAP, dblParam_[CbcAllowableFractionGap]);
3515	delete [] setSol;
3516	delete [] setVar;
3517	char printBuffer[200];
3518	if (offset) {
3519	sprintf(printBuffer, "Add %g to all Cplex messages for true objective",
3520	-offset);
3521	messageHandler()->message(CBC_GENERAL, messages())
3522	<< printBuffer << CoinMessageEol ;
3523	cpxSolver.setDblParam(OsiObjOffset, offset);
3524	}
3525	cpxSolver.branchAndBound();
3526	double timeTaken = CoinCpuTime() - time1;
3527	sprintf(printBuffer, "Cplex took %g seconds",
3528	timeTaken);
3529	messageHandler()->message(CBC_GENERAL, messages())
3530	<< printBuffer << CoinMessageEol ;
3531	numberExtraNodes_ = CPXgetnodecnt(env, lpPtr);
3532	numberExtraIterations_ = CPXgetmipitcnt(env, lpPtr);
3533	double value = cpxSolver.getObjValue() * direction;
3534	if (cpxSolver.isProvenOptimal() && value <= getCutoff()) {
3535	feasible = true;
3536	clpSolver->setWarmStart(NULL);
3537	// try and do solution
3538	double * newSolution =
3539	CoinCopyOfArray(cpxSolver.getColSolution(),
3540	getNumCols());
3541	setBestSolution(CBC_STRONGSOL, value, newSolution) ;
3542	delete [] newSolution;
3543	}
3544	feasible = false;
3545	}
3546	#endif
3547	if (!parentModel_&&(moreSpecialOptions_&268435456) != 0) {
3548	// try idiotic idea
3549	CbcObject * obj = new CbcIdiotBranch(this);
3550	obj->setPriority(1); // temp
3551	addObjects(1, &obj);
3552	delete obj;
3553	}
3554
3555	/*
3556	A hook to use clp to quickly explore some part of the tree.
3557	*/
3558	if (fastNodeDepth_ == 1000 &&/!parentModel_/(specialOptions_&2048) == 0) {
3559	fastNodeDepth_ = -1;
3560	CbcObject * obj =
3561	new CbcFollowOn(this);
3562	obj->setPriority(1);
3563	addObjects(1, &obj);
3564	delete obj;
3565	}
3566	int saveNumberSolves = numberSolves_;
3567	int saveNumberIterations = numberIterations_;
3568	if ((fastNodeDepth_ >= 0\|\|(moreSpecialOptions_&33554432)!=0)
3569	&&/!parentModel_/(specialOptions_&2048) == 0) {
3570	// add in a general depth object doClp
3571	int type = (fastNodeDepth_ <= 100) ? fastNodeDepth_ : -(fastNodeDepth_ - 100);
3572	if ((moreSpecialOptions_&33554432)!=0)
3573	type=12;
3574	else
3575	fastNodeDepth_ += 1000000; // mark as done
3576	CbcObject * obj =
3577	new CbcGeneralDepth(this, type);
3578	addObjects(1, &obj);
3579	delete obj;
3580	// fake number of objects
3581	numberObjects_--;
3582	if (parallelMode() < -1) {
3583	// But make sure position is correct
3584	OsiObject * obj2 = object_[numberObjects_];
3585	obj = dynamic_cast<CbcObject *> (obj2);
3586	assert (obj);
3587	obj->setPosition(numberObjects_);
3588	}
3589	}
3590	#ifdef COIN_HAS_CLP
3591	#ifdef NO_CRUNCH
3592	if (true) {
3593	OsiClpSolverInterface * clpSolver
3594	= dynamic_cast<OsiClpSolverInterface *> (solver_);
3595	if (clpSolver && !parentModel_) {
3596	ClpSimplex * clpSimplex = clpSolver->getModelPtr();
3597	clpSimplex->setSpecialOptions(clpSimplex->specialOptions() \| 131072);
3598	//clpSimplex->startPermanentArrays();
3599	clpSimplex->setPersistenceFlag(2);
3600	}
3601	}
3602	#endif
3603	#endif
3604	// Save copy of solver
3605	OsiSolverInterface * saveSolver = NULL;
3606	if (!parentModel_ && (specialOptions_&(512 + 32768)) != 0)
3607	saveSolver = solver_->clone();
3608	double checkCutoffForRestart = 1.0e100;
3609	saveModel(saveSolver, &checkCutoffForRestart, &feasible);
3610	if ((specialOptions_&262144) != 0 && !parentModel_) {
3611	// Save stuff and return!
3612	storedRowCuts_->saveStuff(bestObjective_, bestSolution_,
3613	solver_->getColLower(),
3614	solver_->getColUpper());
3615	delete [] lowerBefore;
3616	delete [] upperBefore;
3617	delete saveSolver;
3618	if (flipObjective)
3619	flipModel();
3620	return;
3621	}
3622	/*
3623	We've taken the continuous relaxation as far as we can. Time to branch.
3624	The first order of business is to actually create a node. chooseBranch
3625	currently uses strong branching to evaluate branch object candidates,
3626	unless forced back to simple branching. If chooseBranch concludes that a
3627	branching candidate is monotone (anyAction == -1) or infeasible (anyAction
3628	== -2) when forced to integer values, it returns here immediately.
3629
3630	Monotone variables trigger a call to resolve(). If the problem remains
3631	feasible, try again to choose a branching variable. At the end of the loop,
3632	resolved == true indicates that some variables were fixed.
3633
3634	Loss of feasibility will result in the deletion of newNode.
3635	*/
3636
3637	bool resolved = false ;
3638	CbcNode *newNode = NULL ;
3639	numberFixedAtRoot_ = 0;
3640	numberFixedNow_ = 0;
3641	if (!parentModel_&&(moreSpecialOptions2_&2)!=0) {
3642	#ifdef COIN_HAS_CLP
3643	OsiClpSolverInterface * clpSolver
3644	= dynamic_cast<OsiClpSolverInterface *> (solver_);
3645	if (clpSolver) {
3646	if (getCutoff()>1.0e20) {
3647	printf("Zapping costs\n");
3648	int numberColumns=solver_->getNumCols();
3649	double * zeroCost = new double [numberColumns];
3650	// could make small random
3651	memset(zeroCost,0,numberColumns*sizeof(double));
3652	solver_->setObjective(zeroCost);
3653	double objValue = solver_->getObjValue();
3654	solver_->setDblParam(OsiObjOffset,-objValue);
3655	clpSolver->getModelPtr()->setObjectiveValue(objValue);
3656	delete [] zeroCost;
3657	} else {
3658	moreSpecialOptions2_ &= ~2;
3659	}
3660	} else {
3661	#endif
3662	moreSpecialOptions2_ &= ~2;
3663	#ifdef COIN_HAS_CLP
3664	}
3665	#endif
3666	}
3667	int numberIterationsAtContinuous = numberIterations_;
3668	//solverCharacteristics_->setSolver(solver_);
3669	if (feasible) {
3670	// mark all cuts as globally valid
3671	int numberCuts=cuts.sizeRowCuts();
3672	resizeWhichGenerator(0,numberCuts);
3673	for (int i=0;i<numberCuts;i++) {
3674	cuts.rowCutPtr(i)->setGloballyValid();
3675	whichGenerator_[i]=20000+(whichGenerator_[i]%10000);
3676	}
3677	#define HOTSTART -1
3678	#if HOTSTART<0
3679	if (bestSolution_ && !parentModel_ && !hotstartSolution_ &&
3680	(moreSpecialOptions_&1024) != 0 && (specialOptions_&2048) == 0) {
3681	// Set priorities so only branch on ones we need to
3682	// use djs and see if only few branches needed
3683	#ifndef NDEBUG
3684	double integerTolerance = getIntegerTolerance() ;
3685	#endif
3686	bool possible = true;
3687	const double * saveLower = continuousSolver_->getColLower();
3688	const double * saveUpper = continuousSolver_->getColUpper();
3689	for (int i = 0; i < numberObjects_; i++) {
3690	const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object_[i]);
3691	if (thisOne) {
3692	int iColumn = thisOne->columnNumber();
3693	if (saveUpper[iColumn] > saveLower[iColumn] + 1.5) {
3694	possible = false;
3695	break;
3696	}
3697	} else {
3698	possible = false;
3699	break;
3700	}
3701	}
3702	if (possible) {
3703	OsiSolverInterface * solver = continuousSolver_->clone();
3704	int numberColumns = solver->getNumCols();
3705	for (int iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
3706	double value = bestSolution_[iColumn] ;
3707	value = CoinMax(value, saveLower[iColumn]) ;
3708	value = CoinMin(value, saveUpper[iColumn]) ;
3709	value = floor(value + 0.5);
3710	if (solver->isInteger(iColumn)) {
3711	solver->setColLower(iColumn, value);
3712	solver->setColUpper(iColumn, value);
3713	}
3714	}
3715	solver->setHintParam(OsiDoDualInResolve, false, OsiHintTry);
3716	// objlim and all slack
3717	double direction = solver->getObjSense();
3718	solver->setDblParam(OsiDualObjectiveLimit, 1.0e50*direction);
3719	CoinWarmStartBasis * basis = dynamic_cast<CoinWarmStartBasis *> (solver->getEmptyWarmStart());
3720	solver->setWarmStart(basis);
3721	delete basis;
3722	bool changed = true;
3723	hotstartPriorities_ = new int [numberColumns];
3724	for (int iColumn = 0; iColumn < numberColumns; iColumn++)
3725	hotstartPriorities_[iColumn] = 1;
3726	while (changed) {
3727	changed = false;
3728	solver->resolve();
3729	if (!solver->isProvenOptimal()) {
3730	possible = false;
3731	break;
3732	}
3733	const double * dj = solver->getReducedCost();
3734	const double * colLower = solver->getColLower();
3735	const double * colUpper = solver->getColUpper();
3736	const double * solution = solver->getColSolution();
3737	int nAtLbNatural = 0;
3738	int nAtUbNatural = 0;
3739	int nZeroDj = 0;
3740	int nForced = 0;
3741	for (int iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
3742	double value = solution[iColumn] ;
3743	value = CoinMax(value, saveLower[iColumn]) ;
3744	value = CoinMin(value, saveUpper[iColumn]) ;
3745	if (solver->isInteger(iColumn)) {
3746	assert(fabs(value - solution[iColumn]) <= integerTolerance) ;
3747	if (hotstartPriorities_[iColumn] == 1) {
3748	if (dj[iColumn] < -1.0e-6) {
3749	// negative dj
3750	if (saveUpper[iColumn] == colUpper[iColumn]) {
3751	nAtUbNatural++;
3752	hotstartPriorities_[iColumn] = 2;
3753	solver->setColLower(iColumn, saveLower[iColumn]);
3754	solver->setColUpper(iColumn, saveUpper[iColumn]);
3755	} else {
3756	nForced++;
3757	}
3758	} else if (dj[iColumn] > 1.0e-6) {
3759	// positive dj
3760	if (saveLower[iColumn] == colLower[iColumn]) {
3761	nAtLbNatural++;
3762	hotstartPriorities_[iColumn] = 2;
3763	solver->setColLower(iColumn, saveLower[iColumn]);
3764	solver->setColUpper(iColumn, saveUpper[iColumn]);
3765	} else {
3766	nForced++;
3767	}
3768	} else {
3769	// zero dj
3770	nZeroDj++;
3771	}
3772	}
3773	}
3774	}
3775	#if CBC_USEFUL_PRINTING>1
3776	printf("%d forced, %d naturally at lower, %d at upper - %d zero dj\n",
3777	nForced, nAtLbNatural, nAtUbNatural, nZeroDj);
3778	#endif
3779	if (nAtLbNatural \|\| nAtUbNatural) {
3780	changed = true;
3781	} else {
3782	if (nForced + nZeroDj > 5000 \|\|
3783	(nForced + nZeroDj)*2 > numberIntegers_)
3784	possible = false;
3785	}
3786	}
3787	delete solver;
3788	}
3789	if (possible) {
3790	setHotstartSolution(bestSolution_);
3791	if (!saveCompare) {
3792	// create depth first comparison
3793	saveCompare = nodeCompare_;
3794	// depth first
3795	nodeCompare_ = new CbcCompareDepth();
3796	tree_->setComparison(*nodeCompare_) ;
3797	}
3798	} else {
3799	delete [] hotstartPriorities_;
3800	hotstartPriorities_ = NULL;
3801	}
3802	}
3803	#endif
3804	#if HOTSTART>0
3805	if (hotstartSolution_ && !hotstartPriorities_) {
3806	// Set up hot start
3807	OsiSolverInterface * solver = solver_->clone();
3808	double direction = solver_->getObjSense() ;
3809	int numberColumns = solver->getNumCols();
3810	double * saveLower = CoinCopyOfArray(solver->getColLower(), numberColumns);
3811	double * saveUpper = CoinCopyOfArray(solver->getColUpper(), numberColumns);
3812	// move solution
3813	solver->setColSolution(hotstartSolution_);
3814	// point to useful information
3815	const double * saveSolution = testSolution_;
3816	testSolution_ = solver->getColSolution();
3817	OsiBranchingInformation usefulInfo = usefulInformation();
3818	testSolution_ = saveSolution;
3819	/*
3820	Run through the objects and use feasibleRegion() to set variable bounds
3821	so as to fix the variables specified in the objects at their value in this
3822	solution. Since the object list contains (at least) one object for every
3823	integer variable, this has the effect of fixing all integer variables.
3824	*/
3825	for (int i = 0; i < numberObjects_; i++)
3826	object_[i]->feasibleRegion(solver, &usefulInfo);
3827	solver->resolve();
3828	assert (solver->isProvenOptimal());
3829	double gap = CoinMax((solver->getObjValue() - solver_->getObjValue()) * direction, 0.0) ;
3830	const double * dj = solver->getReducedCost();
3831	const double * colLower = solver->getColLower();
3832	const double * colUpper = solver->getColUpper();
3833	const double * solution = solver->getColSolution();
3834	int nAtLbNatural = 0;
3835	int nAtUbNatural = 0;
3836	int nAtLbNaturalZero = 0;
3837	int nAtUbNaturalZero = 0;
3838	int nAtLbFixed = 0;
3839	int nAtUbFixed = 0;
3840	int nAtOther = 0;
3841	int nAtOtherNatural = 0;
3842	int nNotNeeded = 0;
3843	delete [] hotstartSolution_;
3844	hotstartSolution_ = new double [numberColumns];
3845	delete [] hotstartPriorities_;
3846	hotstartPriorities_ = new int [numberColumns];
3847	int * order = (int *) saveUpper;
3848	int nFix = 0;
3849	double bestRatio = COIN_DBL_MAX;
3850	for (int iColumn = 0 ; iColumn < numberColumns ; iColumn++) {
3851	double value = solution[iColumn] ;
3852	value = CoinMax(value, saveLower[iColumn]) ;
3853	value = CoinMin(value, saveUpper[iColumn]) ;
3854	double sortValue = COIN_DBL_MAX;
3855	if (solver->isInteger(iColumn)) {
3856	assert(fabs(value - solution[iColumn]) <= 1.0e-5) ;
3857	double value2 = floor(value + 0.5);
3858	if (dj[iColumn] < -1.0e-6) {
3859	// negative dj
3860	//assert (value2==colUpper[iColumn]);
3861	if (saveUpper[iColumn] == colUpper[iColumn]) {
3862	nAtUbNatural++;
3863	sortValue = 0.0;
3864	double value = -dj[iColumn];
3865	if (value > gap)
3866	nFix++;
3867	else if (gap < value*bestRatio)
3868	bestRatio = gap / value;
3869	if (saveLower[iColumn] != colLower[iColumn]) {
3870	nNotNeeded++;
3871	sortValue = 1.0e20;
3872	}
3873	} else if (saveLower[iColumn] == colUpper[iColumn]) {
3874	nAtLbFixed++;
3875	sortValue = dj[iColumn];
3876	} else {
3877	nAtOther++;
3878	sortValue = 0.0;
3879	if (saveLower[iColumn] != colLower[iColumn] &&
3880	saveUpper[iColumn] != colUpper[iColumn]) {
3881	nNotNeeded++;
3882	sortValue = 1.0e20;
3883	}
3884	}
3885	} else if (dj[iColumn] > 1.0e-6) {
3886	// positive dj
3887	//assert (value2==colLower[iColumn]);
3888	if (saveLower[iColumn] == colLower[iColumn]) {
3889	nAtLbNatural++;
3890	sortValue = 0.0;
3891	double value = dj[iColumn];
3892	if (value > gap)
3893	nFix++;
3894	else if (gap < value*bestRatio)
3895	bestRatio = gap / value;
3896	if (saveUpper[iColumn] != colUpper[iColumn]) {
3897	nNotNeeded++;
3898	sortValue = 1.0e20;
3899	}
3900	} else if (saveUpper[iColumn] == colLower[iColumn]) {
3901	nAtUbFixed++;
3902	sortValue = -dj[iColumn];
3903	} else {
3904	nAtOther++;
3905	sortValue = 0.0;
3906	if (saveLower[iColumn] != colLower[iColumn] &&
3907	saveUpper[iColumn] != colUpper[iColumn]) {
3908	nNotNeeded++;
3909	sortValue = 1.0e20;
3910	}
3911	}
3912	} else {
3913	// zero dj
3914	if (value2 == saveUpper[iColumn]) {
3915	nAtUbNaturalZero++;
3916	sortValue = 0.0;
3917	if (saveLower[iColumn] != colLower[iColumn]) {
3918	nNotNeeded++;
3919	sortValue = 1.0e20;
3920	}
3921	} else if (value2 == saveLower[iColumn]) {
3922	nAtLbNaturalZero++;
3923	sortValue = 0.0;
3924	} else {
3925	nAtOtherNatural++;
3926	sortValue = 0.0;
3927	if (saveLower[iColumn] != colLower[iColumn] &&
3928	saveUpper[iColumn] != colUpper[iColumn]) {
3929	nNotNeeded++;
3930	sortValue = 1.0e20;
3931	}
3932	}
3933	}
3934	#if HOTSTART==3
3935	sortValue = -fabs(dj[iColumn]);
3936	#endif
3937	}
3938	hotstartSolution_[iColumn] = value ;
3939	saveLower[iColumn] = sortValue;
3940	order[iColumn] = iColumn;
3941	}
3942	COIN_DETAIL_PRINT(printf("** can fix %d columns - best ratio for others is %g on gap of %g\n",
3943	nFix, bestRatio, gap));
3944	int nNeg = 0;
3945	CoinSort_2(saveLower, saveLower + numberColumns, order);
3946	for (int i = 0; i < numberColumns; i++) {
3947	if (saveLower[i] < 0.0) {
3948	nNeg++;
3949	#if HOTSTART==2\|\|HOTSTART==3
3950	// swap sign ?
3951	saveLower[i] = -saveLower[i];
3952	#endif
3953	}
3954	}
3955	CoinSort_2(saveLower, saveLower + nNeg, order);
3956	for (int i = 0; i < numberColumns; i++) {
3957	#if HOTSTART==1
3958	hotstartPriorities_[order[i]] = 100;
3959	#else
3960	hotstartPriorities_[order[i]] = -(i + 1);
3961	#endif
3962	}
3963	COIN_DETAIL_PRINT(printf("nAtLbNat %d,nAtUbNat %d,nAtLbNatZero %d,nAtUbNatZero %d,nAtLbFixed %d,nAtUbFixed %d,nAtOther %d,nAtOtherNat %d, useless %d %d\n",
3964	nAtLbNatural,
3965	nAtUbNatural,
3966	nAtLbNaturalZero,
3967	nAtUbNaturalZero,
3968	nAtLbFixed,
3969	nAtUbFixed,
3970	nAtOther,
3971	nAtOtherNatural, nNotNeeded, nNeg));
3972	delete [] saveLower;
3973	delete [] saveUpper;
3974	if (!saveCompare) {
3975	// create depth first comparison
3976	saveCompare = nodeCompare_;
3977	// depth first
3978	nodeCompare_ = new CbcCompareDepth();
3979	tree_->setComparison(*nodeCompare_) ;
3980	}
3981	}
3982	#endif
3983	newNode = new CbcNode ;
3984	// Set objective value (not so obvious if NLP etc)
3985	setObjectiveValue(newNode, NULL);
3986	anyAction = -1 ;
3987	// To make depth available we may need a fake node
3988	CbcNode fakeNode;
3989	if (!currentNode_) {
3990	// Not true if sub trees assert (!numberNodes_);
3991	currentNode_ = &fakeNode;
3992	}
3993	phase_ = 3;
3994	// only allow 1000 passes
3995	int numberPassesLeft = 1000;
3996	// This is first crude step
3997	if (numberAnalyzeIterations_ && !parentModel_) {
3998	delete [] analyzeResults_;
3999	//int numberColumns = solver_->getNumCols();
4000	analyzeResults_ = new double [5*numberIntegers_];
4001	numberFixedAtRoot_ = newNode->analyze(this, analyzeResults_);
4002	if (numberFixedAtRoot_ > 0) {
4003	COIN_DETAIL_PRINT(printf("%d fixed by analysis\n", numberFixedAtRoot_));
4004	setPointers(solver_);
4005	numberFixedNow_ = numberFixedAtRoot_;
4006	} else if (numberFixedAtRoot_ < 0) {
4007	COIN_DETAIL_PRINT(printf("analysis found to be infeasible\n"));
4008	anyAction = -2;
4009	delete newNode ;
4010	newNode = NULL ;
4011	feasible = false ;
4012	}
4013	}
4014	OsiSolverBranch * branches = NULL;
4015	if (feasible)
4016	anyAction = chooseBranch(newNode, numberPassesLeft, NULL, cuts, resolved,
4017	NULL, NULL, NULL, branches);
4018	if (anyAction == -2 \|\| newNode->objectiveValue() >= cutoff) {
4019	if (anyAction != -2) {
4020	// zap parent nodeInfo
4021	#ifdef COIN_DEVELOP
4022	printf("zapping CbcNodeInfo %x\n", newNode->nodeInfo()->parent());
4023	#endif
4024	if (newNode->nodeInfo())
4025	newNode->nodeInfo()->nullParent();
4026	}
4027	delete newNode ;
4028	newNode = NULL ;
4029	feasible = false ;
4030	}
4031	}
4032	if (newNode && probingInfo_) {
4033	int number01 = probingInfo_->numberIntegers();
4034	//const fixEntry * entry = probingInfo_->fixEntries();
4035	const int * toZero = probingInfo_->toZero();
4036	//const int * toOne = probingInfo_->toOne();
4037	//const int * integerVariable = probingInfo_->integerVariable();
4038	if (toZero[number01]) {
4039	CglTreeProbingInfo info(*probingInfo_);
4040	if ((moreSpecialOptions2_&64)!=0&&!parentModel_) {
4041	/*
4042	Marginal idea. Further exploration probably good. Build some extra
4043	cliques from probing info. Not quite worth the effort?
4044	*/
4045	CglProbing generator1;
4046	generator1.setUsingObjective(false);
4047	generator1.setMaxPass(1);
4048	generator1.setMaxPassRoot(1);
4049	generator1.setMaxLook(100);
4050	generator1.setRowCuts(3);
4051	generator1.setMaxElements(300);
4052	generator1.setMaxProbeRoot(solver_->getNumCols());
4053	CoinThreadRandom randomGenerator;
4054	//CglTreeProbingInfo info(solver_);
4055	info.level = 0;
4056	info.formulation_rows = solver_->getNumRows();
4057	info.inTree = false;
4058	if (parentModel_) {
4059	info.parentSolver=parentModel_->continuousSolver();
4060	// indicate if doing full search
4061	info.hasParent = ((specialOptions_&67108864)==0) ? 1 : 2;
4062	} else {
4063	info.hasParent=0;
4064	info.parentSolver=NULL;
4065	}
4066	info.originalColumns=originalColumns();
4067	info.randomNumberGenerator=&randomGenerator;
4068	info.pass=4;
4069	generator1.setMode(8);
4070	OsiCuts cs;
4071	generator1.generateCutsAndModify(*solver_,cs,&info);
4072	// very clunky
4073	OsiSolverInterface * temp = generator1.cliqueModel(solver_,2);
4074	CglPreProcess dummy;
4075	OsiSolverInterface * newSolver=dummy.cliqueIt(*temp,0.0001);
4076	delete temp;
4077	OsiSolverInterface * fake = NULL;
4078	if (newSolver) {
4079	#if 0
4080	int numberCliques = generator1.numberCliques();
4081	cliqueEntry * entry = generator1.cliqueEntry();
4082	cliqueType * type = new cliqueType [numberCliques];
4083	int * start = new int [numberCliques+1];
4084	start[numberCliques]=2*numberCliques;
4085	int n=0;
4086	for (int i=0;i<numberCliques;i++) {
4087	start[i]=2*i;
4088	setOneFixesInCliqueEntry(entry[2*i],true);
4089	setOneFixesInCliqueEntry(entry[2*i+1],true);
4090	type[i]=0;
4091	}
4092	fake = info.analyze(*solver_, 1,numberCliques,start,
4093	entry,type);
4094	delete [] type;
4095	delete [] entry;
4096	#else
4097	fake = info.analyze(*newSolver, 1,-1);
4098	#endif
4099	delete newSolver;
4100	} else {
4101	fake = info.analyze(*solver_, 1);
4102	}
4103	if (fake) {
4104	//fake->writeMps("fake");
4105	CglFakeClique cliqueGen(fake);
4106	cliqueGen.setStarCliqueReport(false);
4107	cliqueGen.setRowCliqueReport(false);
4108	cliqueGen.setMinViolation(0.1);
4109	addCutGenerator(&cliqueGen, 1, "Fake cliques", true, false, false, -200);
4110	generator_[numberCutGenerators_-1]->setTiming(true);
4111	for (int i = 0; i < numberCutGenerators_; i++) {
4112	CglKnapsackCover * cutGen =
4113	dynamic_cast<CglKnapsackCover *>(generator_[i]->generator());
4114	if (cutGen) {
4115	cutGen->createCliques(*fake,2,200,false);
4116	}
4117	}
4118	}
4119	}
4120	if (probingInfo_->packDown()) {
4121	#if CBC_USEFUL_PRINTING>1
4122	printf("%d implications on %d 0-1\n", toZero[number01], number01);
4123	#endif
4124	// Create a cut generator that remembers implications discovered at root.
4125	CglImplication implication(probingInfo_);
4126	addCutGenerator(&implication, 1, "ImplicationCuts", true, false, false, -200);
4127	generator_[numberCutGenerators_-1]->setGlobalCuts(true);
4128	generator_[numberCutGenerators_-1]->setTiming(true);
4129	} else {
4130	delete probingInfo_;
4131	probingInfo_ = NULL;
4132	}
4133	} else {
4134	delete probingInfo_;
4135
4136	probingInfo_ = NULL;
4137	}
4138	}
4139	/*
4140	At this point, the root subproblem is infeasible or fathomed by bound
4141	(newNode == NULL), or we're live with an objective value that satisfies the
4142	current objective cutoff.
4143	*/
4144	assert (!newNode \|\| newNode->objectiveValue() <= cutoff) ;
4145	// Save address of root node as we don't want to delete it
4146	/*
4147	The common case is that the lp relaxation is feasible but doesn't satisfy
4148	integrality (i.e., newNode->branchingObject(), indicating we've been able to
4149	select a branching variable). Remove any cuts that have gone slack due to
4150	forcing monotone variables. Then tack on an CbcFullNodeInfo object and full
4151	basis (via createInfo()) and stash the new cuts in the nodeInfo (via
4152	addCuts()). If, by some miracle, we have an integral solution at the root
4153	(newNode->branchingObject() is NULL), takeOffCuts() will ensure that the solver holds
4154	a valid solution for use by setBestSolution().
4155	*/
4156	CoinWarmStartBasis *lastws = NULL ;
4157	if (feasible && newNode->branchingObject()) {
4158	if (resolved) {
4159	takeOffCuts(cuts, false, NULL) ;
4160	# ifdef CHECK_CUT_COUNTS
4161	{ printf("Number of rows after chooseBranch fix (root)"
4162	"(active only) %d\n",
4163	numberRowsAtContinuous_ + numberNewCuts_ + numberOldActiveCuts_) ;
4164	const CoinWarmStartBasis* debugws =
4165	dynamic_cast <const CoinWarmStartBasis*>(solver_->getWarmStart()) ;
4166	debugws->print() ;
4167	delete debugws ;
4168	}
4169	# endif
4170	}
4171	//newNode->createInfo(this,NULL,NULL,NULL,NULL,0,0) ;
4172	//newNode->nodeInfo()->addCuts(cuts,
4173	// newNode->numberBranches(),whichGenerator_) ;
4174	if (lastws) delete lastws ;
4175	lastws = dynamic_cast<CoinWarmStartBasis*>(solver_->getWarmStart()) ;
4176	}
4177	/*
4178	Continuous data to be used later
4179	*/
4180	continuousObjective_ = solver_->getObjValue() * solver_->getObjSense();
4181	continuousInfeasibilities_ = 0 ;
4182	if (newNode) {
4183	continuousObjective_ = newNode->objectiveValue() ;
4184	delete [] continuousSolution_;
4185	continuousSolution_ = CoinCopyOfArray(solver_->getColSolution(),
4186	numberColumns);
4187	continuousInfeasibilities_ = newNode->numberUnsatisfied() ;
4188	}
4189	/*
4190	Bound may have changed so reset in objects
4191	*/
4192	{
4193	int i ;
4194	for (i = 0; i < numberObjects_; i++)
4195	object_[i]->resetBounds(solver_) ;
4196	}
4197	/*
4198	Feasible? Then we should have either a live node prepped for future
4199	expansion (indicated by variable() >= 0), or (miracle of miracles) an
4200	integral solution at the root node.
4201
4202	initializeInfo sets the reference counts in the nodeInfo object. Since
4203	this node is still live, push it onto the heap that holds the live set.
4204	*/
4205	if (newNode) {
4206	if (newNode->branchingObject()) {
4207	newNode->initializeInfo() ;
4208	tree_->push(newNode) ;
4209	// save pointer to root node - so can pick up bounds
4210	if (!topOfTree_)
4211	topOfTree_ = dynamic_cast<CbcFullNodeInfo *>(newNode->nodeInfo()) ;
4212	if (statistics_) {
4213	if (numberNodes2_ == maximumStatistics_) {
4214	maximumStatistics_ = 2 * maximumStatistics_;
4215	CbcStatistics ** temp = new CbcStatistics * [maximumStatistics_];
4216	memset(temp, 0, maximumStatistics_sizeof(CbcStatistics ));
4217	memcpy(temp, statistics_, numberNodes2_sizeof(CbcStatistics ));
4218	delete [] statistics_;
4219	statistics_ = temp;
4220	}
4221	assert (!statistics_[numberNodes2_]);
4222	statistics_[numberNodes2_] = new CbcStatistics(newNode, this);
4223	}
4224	numberNodes2_++;
4225	# ifdef CHECK_NODE
4226	printf("Node %x on tree\n", newNode) ;
4227	# endif
4228	} else {
4229	// continuous is integer
4230	double objectiveValue = newNode->objectiveValue();
4231	setBestSolution(CBC_SOLUTION, objectiveValue,
4232	solver_->getColSolution()) ;
4233	if (eventHandler) {
4234	// we are stopping anyway so no need to test return code
4235	eventHandler->event(CbcEventHandler::solution);
4236	}
4237	delete newNode ;
4238	newNode = NULL ;
4239	}
4240	}
4241
4242	if (printFrequency_ <= 0) {
4243	printFrequency_ = 1000 ;
4244	if (getNumCols() > 2000)
4245	printFrequency_ = 100 ;
4246	}
4247	/*
4248	It is possible that strong branching fixes one variable and then the code goes round
4249	again and again. This can take too long. So we need to warn user - just once.
4250	*/
4251	numberLongStrong_ = 0;
4252	CbcNode * createdNode = NULL;
4253	#ifdef CBC_THREAD
4254	if ((specialOptions_&2048) != 0)
4255	numberThreads_ = 0;
4256	if (numberThreads_ ) {
4257	nodeCompare_->sayThreaded(); // need to use addresses
4258	master_ = new CbcBaseModel(*this,
4259	(parallelMode() < -1) ? 1 : 0);
4260	masterThread_ = master_->masterThread();
4261	}
4262	#endif
4263	#ifdef COIN_HAS_CLP
4264	{
4265	OsiClpSolverInterface * clpSolver
4266	= dynamic_cast<OsiClpSolverInterface *> (solver_);
4267	if (clpSolver && !parentModel_) {
4268	clpSolver->computeLargestAway();
4269	}
4270	}
4271	#endif
4272	/*
4273	At last, the actual branch-and-cut search loop, which will iterate until
4274	the live set is empty or we hit some limit (integrality gap, time, node
4275	count, etc.). The overall flow is to rebuild a subproblem, reoptimise using
4276	solveWithCuts(), choose a branching pattern with chooseBranch(), and finally
4277	add the node to the live set.
4278
4279	The first action is to winnow the live set to remove nodes which are worse
4280	than the current objective cutoff.
4281	*/
4282	if (solver_->getRowCutDebuggerAlways()) {
4283	OsiRowCutDebugger * debuggerX = const_cast<OsiRowCutDebugger *> (solver_->getRowCutDebuggerAlways());
4284	const OsiRowCutDebugger *debugger = solver_->getRowCutDebugger() ;
4285	if (!debugger) {
4286	// infeasible!!
4287	printf("before search\n");
4288	const double * lower = solver_->getColLower();
4289	const double * upper = solver_->getColUpper();
4290	const double * solution = debuggerX->optimalSolution();
4291	int numberColumns = solver_->getNumCols();
4292	for (int i = 0; i < numberColumns; i++) {
4293	if (solver_->isInteger(i)) {
4294	if (solution[i] < lower[i] - 1.0e-6 \|\| solution[i] > upper[i] + 1.0e-6)
4295	printf("**** ");
4296	printf("%d %g <= %g <= %g\n",
4297	i, lower[i], solution[i], upper[i]);
4298	}
4299	}
4300	//abort();
4301	}
4302	}
4303	{
4304	// may be able to change cutoff now
4305	double cutoff = getCutoff();
4306	double increment = getDblParam(CbcModel::CbcCutoffIncrement) ;
4307	if (cutoff > bestObjective_ - increment) {
4308	cutoff = bestObjective_ - increment ;
4309	setCutoff(cutoff) ;
4310	}
4311	}
4312	#ifdef CBC_THREAD
4313	bool goneParallel = false;
4314	#endif
4315	#define MAX_DEL_NODE 1
4316	CbcNode * delNode[MAX_DEL_NODE+1];
4317	int nDeleteNode = 0;
4318	// For Printing etc when parallel
4319	int lastEvery1000 = 0;
4320	int lastPrintEvery = 0;
4321	int numberConsecutiveInfeasible = 0;
4322	#define PERTURB_IN_FATHOM
4323	#ifdef PERTURB_IN_FATHOM
4324	// allow in fathom
4325	if ((moreSpecialOptions_& 262144) != 0)
4326	specialOptions_ \|= 131072;
4327	#endif
4328	while (true) {
4329	lockThread();
4330	#ifdef COIN_HAS_CLP
4331	// See if we want dantzig row choice
4332	goToDantzig(100, savePivotMethod);
4333	#endif
4334	//#define REPORT_DYNAMIC 2
4335	#if REPORT_DYNAMIC
4336	if (numberNodes_&&!parentModel_&&(tree_->empty()\|\|(numberNodes_%10000)==0)) {
4337	// Get average up and down costs
4338	double averageUp = 0.0;
4339	double averageDown = 0.0;
4340	int numberUp = 0;
4341	int numberDown = 0;
4342	int minTimesDown = COIN_INT_MAX;
4343	int maxTimesDown = 0;
4344	int neverBranchedDown = 0;
4345	int infeasibleTimesDown = 0;
4346	int minTimesUp = COIN_INT_MAX;
4347	int maxTimesUp = 0;
4348	int infeasibleTimesUp = 0;
4349	int neverBranchedUp = 0;
4350	int neverBranched = 0;
4351	int i;
4352	int numberInts=0;
4353	bool endOfSearch = tree_->empty();
4354	int numberUp2 = 0;
4355	int numberDown2 = 0;
4356	for ( i = 0; i < numberObjects_; i++) {
4357	OsiObject * object = object_[i];
4358	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
4359	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
4360	if (dynamicObject) {
4361	numberInts++;
4362	if (dynamicObject->numberTimesUp()\|\|
4363	dynamicObject->numberTimesDown()) {
4364	int nUp = 0;
4365	int nDown = 0;
4366	double up = 0.0;
4367	double down = 0.0;
4368	if (dynamicObject->numberTimesUp()) {
4369	numberUp++;
4370	nUp = dynamicObject->numberTimesUp();
4371	minTimesUp = CoinMin(minTimesUp,nUp);
4372	maxTimesUp = CoinMax(maxTimesUp,nUp);
4373	up = dynamicObject->upDynamicPseudoCost();
4374	averageUp += up;
4375	numberUp2 += nUp;
4376	infeasibleTimesUp += dynamicObject->numberTimesUpInfeasible();
4377	} else {
4378	neverBranchedUp++;
4379	}
4380	if (dynamicObject->numberTimesDown()) {
4381	numberDown++;
4382	nDown = dynamicObject->numberTimesDown();
4383	minTimesDown = CoinMin(minTimesDown,nDown);
4384	maxTimesDown = CoinMax(maxTimesDown,nDown);
4385	down = dynamicObject->downDynamicPseudoCost();
4386	averageDown += down;
4387	numberDown2 += dynamicObject->numberTimesDown();
4388	infeasibleTimesDown += dynamicObject->numberTimesDownInfeasible();
4389	} else {
4390	neverBranchedDown++;
4391	}
4392	#if REPORT_DYNAMIC > 1
4393	#if REPORT_DYNAMIC == 2
4394	if (endOfSearch&&numberIntegers_<400) {
4395	#elif REPORT_DYNAMIC == 3
4396	if (endOfSearch) {
4397	#else
4398	{
4399	#endif
4400	dynamicObject->print(0,0.0);
4401	}
4402	#endif
4403	} else {
4404	neverBranched++;
4405	#if REPORT_DYNAMIC > 2
4406	#if REPORT_DYNAMIC == 3
4407	if (endOfSearch&&numberIntegers_<400) {
4408	#elif REPORT_DYNAMIC == 4
4409	if (endOfSearch) {
4410	#else
4411	{
4412	#endif
4413	printf("col %d - never branched on\n",dynamicObject->columnNumber());
4414	}
4415	#endif
4416	}
4417	}
4418	}
4419	if (numberUp)
4420	averageUp /= static_cast<double> (numberUp);
4421	else
4422	averageUp = 0.0;
4423	if (numberDown)
4424	averageDown /= static_cast<double> (numberDown);
4425	else
4426	averageDown = 0.0;
4427	printf("Report for %d variables (%d never branched on) after %d nodes - total solves down %d up %d\n",
4428	numberInts,neverBranched,numberNodes_,numberDown2,numberUp2);
4429	if ((neverBranchedDown\|\|neverBranchedUp)&&endOfSearch)
4430	printf("odd %d never branched down and %d never branched up\n",
4431	neverBranchedDown,neverBranchedUp);
4432	printf("down average %g times (%d infeasible) average increase %g min/max times (%d,%d)\n",
4433	static_cast<double>(numberDown2)/numberDown,infeasibleTimesDown,averageDown,
4434	minTimesDown,maxTimesDown);
4435	printf("up average %g times (%d infeasible) average increase %g min/max times (%d,%d)\n",
4436	static_cast<double>(numberUp2)/numberUp,infeasibleTimesUp,averageUp,
4437	minTimesUp,maxTimesUp);
4438	}
4439	#endif
4440	if (tree_->empty()) {
4441	#ifdef CBC_THREAD
4442	if (parallelMode() > 0 && master_) {
4443	int anyLeft = master_->waitForThreadsInTree(0);
4444	if (!anyLeft) {
4445	master_->stopThreads(-1);
4446	break;
4447	}
4448	} else {
4449	break;
4450	}
4451	#else
4452	break;
4453	#endif
4454	} else {
4455	unlockThread();
4456	}
4457	// If done 50/100 nodes see if worth trying reduction
4458	if (numberNodes_ >= nextCheckRestart) {
4459	if (nextCheckRestart<100)
4460	nextCheckRestart=100;
4461	else
4462	nextCheckRestart=COIN_INT_MAX;
4463	#ifdef COIN_HAS_CLP
4464	OsiClpSolverInterface * clpSolver
4465	= dynamic_cast<OsiClpSolverInterface *> (solver_);
4466	if (clpSolver && ((specialOptions_&131072) == 0) && true) {
4467	ClpSimplex * simplex = clpSolver->getModelPtr();
4468	int perturbation = simplex->perturbation();
4469	#if CBC_USEFUL_PRINTING>1
4470	printf("Testing its n,s %d %d solves n,s %d %d - pert %d\n",
4471	numberIterations_, saveNumberIterations,
4472	numberSolves_, saveNumberSolves, perturbation);
4473	#endif
4474	if (perturbation == 50 && (numberIterations_ - saveNumberIterations) <
4475	8*(numberSolves_ - saveNumberSolves)) {
4476	// switch off perturbation
4477	simplex->setPerturbation(100);
4478	#if CBC_USEFUL_PRINTING>1
4479	printf("Perturbation switched off\n");
4480	#endif
4481	}
4482	}
4483	#endif
4484	/*
4485	Decide if we want to do a restart.
4486	*/
4487	if (saveSolver && (specialOptions_&(512 + 32768)) != 0) {
4488	bool tryNewSearch = solverCharacteristics_->reducedCostsAccurate() &&
4489	(getCutoff() < 1.0e20 && getCutoff() < checkCutoffForRestart);
4490	int numberColumns = getNumCols();
4491	if (tryNewSearch) {
4492	// adding increment back allows current best - tiny bit weaker
4493	checkCutoffForRestart = getCutoff() + getCutoffIncrement() ;
4494	#if CBC_USEFUL_PRINTING>1
4495	printf("after %d nodes, cutoff %g - looking\n",
4496	numberNodes_, getCutoff());
4497	#endif
4498	saveSolver->resolve();
4499	double direction = saveSolver->getObjSense() ;
4500	double gap = checkCutoffForRestart - saveSolver->getObjValue() * direction ;
4501	double tolerance;
4502	saveSolver->getDblParam(OsiDualTolerance, tolerance) ;
4503	if (gap <= 0.0)
4504	gap = tolerance;
4505	gap += 100.0 * tolerance;
4506	double integerTolerance = getDblParam(CbcIntegerTolerance) ;
4507
4508	const double *lower = saveSolver->getColLower() ;
4509	const double *upper = saveSolver->getColUpper() ;
4510	const double *solution = saveSolver->getColSolution() ;
4511	const double *reducedCost = saveSolver->getReducedCost() ;
4512
4513	int numberFixed = 0 ;
4514	int numberFixed2 = 0;
4515	#ifdef COIN_DEVELOP
4516	printf("gap %g\n", gap);
4517	#endif
4518	for (int i = 0 ; i < numberIntegers_ ; i++) {
4519	int iColumn = integerVariable_[i] ;
4520	double djValue = direction * reducedCost[iColumn] ;
4521	if (upper[iColumn] - lower[iColumn] > integerTolerance) {
4522	if (solution[iColumn] < lower[iColumn] + integerTolerance && djValue > gap) {
4523	//printf("%d to lb on dj of %g - bounds %g %g\n",
4524	// iColumn,djValue,lower[iColumn],upper[iColumn]);
4525	saveSolver->setColUpper(iColumn, lower[iColumn]) ;
4526	numberFixed++ ;
4527	} else if (solution[iColumn] > upper[iColumn] - integerTolerance && -djValue > gap) {
4528	//printf("%d to ub on dj of %g - bounds %g %g\n",
4529	// iColumn,djValue,lower[iColumn],upper[iColumn]);
4530	saveSolver->setColLower(iColumn, upper[iColumn]) ;
4531	numberFixed++ ;
4532	}
4533	} else {
4534	//printf("%d has dj of %g - already fixed to %g\n",
4535	// iColumn,djValue,lower[iColumn]);
4536	numberFixed2++;
4537	}
4538	}
4539	#ifdef COIN_DEVELOP
4540	if ((specialOptions_&1) != 0) {
4541	const OsiRowCutDebugger *debugger = saveSolver->getRowCutDebugger() ;
4542	if (debugger) {
4543	printf("Contains optimal\n") ;
4544	OsiSolverInterface * temp = saveSolver->clone();
4545	const double * solution = debugger->optimalSolution();
4546	const double *lower = temp->getColLower() ;
4547	const double *upper = temp->getColUpper() ;
4548	int n = temp->getNumCols();
4549	for (int i = 0; i < n; i++) {
4550	if (temp->isInteger(i)) {
4551	double value = floor(solution[i] + 0.5);
4552	assert (value >= lower[i] && value <= upper[i]);
4553	temp->setColLower(i, value);
4554	temp->setColUpper(i, value);
4555	}
4556	}
4557	temp->writeMps("reduced_fix");
4558	delete temp;
4559	saveSolver->writeMps("reduced");
4560	} else {
4561	abort();
4562	}
4563	}
4564	printf("Restart could fix %d integers (%d already fixed)\n",
4565	numberFixed + numberFixed2, numberFixed2);
4566	#endif
4567	numberFixed += numberFixed2;
4568	if (numberFixed10 < numberColumns && numberFixed4 < numberIntegers_)
4569	tryNewSearch = false;
4570	}
4571	#ifdef CONFLICT_CUTS
4572	// temporary
4573	//if ((moreSpecialOptions_&4194304)!=0)
4574	//tryNewSearch=false;
4575	#endif
4576	if (tryNewSearch) {
4577	// back to solver without cuts?
4578	OsiSolverInterface * solver2 = saveSolver->clone();
4579	const double *lower = saveSolver->getColLower() ;
4580	const double *upper = saveSolver->getColUpper() ;
4581	for (int i = 0 ; i < numberIntegers_ ; i++) {
4582	int iColumn = integerVariable_[i] ;
4583	solver2->setColLower(iColumn, lower[iColumn]);
4584	solver2->setColUpper(iColumn, upper[iColumn]);
4585	}
4586	// swap
4587	delete saveSolver;
4588	saveSolver = solver2;
4589	double * newSolution = new double[numberColumns];
4590	double objectiveValue = checkCutoffForRestart;
4591	// Save the best solution so far.
4592	CbcSerendipity heuristic(*this);
4593	if (bestSolution_)
4594	heuristic.setInputSolution(bestSolution_, bestObjective_);
4595	// Magic number
4596	heuristic.setFractionSmall(0.8);
4597	// `pumpTune' to stand-alone solver for explanations.
4598	heuristic.setFeasibilityPumpOptions(1008013);
4599	// Use numberNodes to say how many are original rows
4600	heuristic.setNumberNodes(continuousSolver_->getNumRows());
4601	#ifdef COIN_DEVELOP
4602	if (continuousSolver_->getNumRows() <
4603	solver_->getNumRows())
4604	printf("%d rows added ZZZZZ\n",
4605	solver_->getNumRows() - continuousSolver_->getNumRows());
4606	#endif
4607	int returnCode = heuristic.smallBranchAndBound(saveSolver,
4608	-1, newSolution,
4609	objectiveValue,
4610	checkCutoffForRestart, "Reduce");
4611	if (returnCode < 0) {
4612	#ifdef COIN_DEVELOP
4613	printf("Restart - not small enough to do search after fixing\n");
4614	#endif
4615	delete [] newSolution;
4616	} else {
4617	// 1 for sol'n, 2 for finished, 3 for both
4618	if ((returnCode&1) != 0) {
4619	// increment number of solutions so other heuristics can test
4620	numberSolutions_++;
4621	numberHeuristicSolutions_++;
4622	lastHeuristic_ = NULL;
4623	setBestSolution(CBC_ROUNDING, objectiveValue, newSolution) ;
4624	}
4625	delete [] newSolution;
4626	#ifdef CBC_THREAD
4627	if (master_) {
4628	lockThread();
4629	if (parallelMode() > 0) {
4630	while (master_->waitForThreadsInTree(0)) {
4631	lockThread();
4632	double dummyBest;
4633	tree_->cleanTree(this, -COIN_DBL_MAX, dummyBest) ;
4634	//unlockThread();
4635	}
4636	} else {
4637	double dummyBest;
4638	tree_->cleanTree(this, -COIN_DBL_MAX, dummyBest) ;
4639	}
4640	master_->waitForThreadsInTree(2);
4641	delete master_;
4642	master_ = NULL;
4643	masterThread_ = NULL;
4644	}
4645	#endif
4646	if (tree_->size()) {
4647	double dummyBest;
4648	tree_->cleanTree(this, -COIN_DBL_MAX, dummyBest) ;
4649	}
4650	break;
4651	}
4652	}
4653	delete saveSolver;
4654	saveSolver = NULL;
4655	}
4656	}
4657	/*
4658	Check for abort on limits: node count, solution count, time, integrality gap.
4659	*/
4660	if (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
4661	numberSolutions_ < intParam_[CbcMaxNumSol] &&
4662	!maximumSecondsReached() &&
4663	!stoppedOnGap_ && !eventHappened_ && (maximumNumberIterations_ < 0 \|\|
4664	numberIterations_ < maximumNumberIterations_))) {
4665	// out of loop
4666	break;
4667	}
4668	#ifdef BONMIN
4669	assert(!solverCharacteristics_->solutionAddsCuts() \|\| solverCharacteristics_->mipFeasible());
4670	#endif
4671	// Sets percentage of time when we try diving. Diving requires a bit of heap reorganisation, because
4672	// we need to replace the comparison function to dive, and that requires reordering to retain the
4673	// heap property.
4674	#define DIVE_WHEN 1000
4675	#define DIVE_STOP 2000
4676	int kNode = numberNodes_ % 4000;
4677	if (numberNodes_<100000 && kNode>DIVE_WHEN && kNode <= DIVE_STOP) {
4678	if (!parallelMode()) {
4679	if (kNode == DIVE_WHEN + 1 \|\| numberConsecutiveInfeasible > 1) {
4680	CbcCompareDefault * compare = dynamic_cast<CbcCompareDefault *>
4681	(nodeCompare_);
4682	// Don't interfere if user has replaced the compare function.
4683	if (compare) {
4684	//printf("Redoing tree\n");
4685	compare->startDive(this);
4686	numberConsecutiveInfeasible = 0;
4687	}
4688	}
4689	}
4690	}
4691	// replace current cutoff?
4692	if (cutoff > getCutoff()) {
4693	double newCutoff = getCutoff();
4694	if (analyzeResults_) {
4695	// see if we could fix any (more)
4696	int n = 0;
4697	double * newLower = analyzeResults_;
4698	double * objLower = newLower + numberIntegers_;
4699	double * newUpper = objLower + numberIntegers_;
4700	double * objUpper = newUpper + numberIntegers_;
4701	for (int i = 0; i < numberIntegers_; i++) {
4702	if (objLower[i] > newCutoff) {
4703	n++;
4704	if (objUpper[i] > newCutoff) {
4705	newCutoff = -COIN_DBL_MAX;
4706	break;
4707	}
4708	// add as global cut
4709	objLower[i]=-COIN_DBL_MAX;
4710	OsiRowCut rc;
4711	rc.setLb(newLower[i]);
4712	rc.setUb(COIN_DBL_MAX);
4713	double one=1.0;
4714	rc.setRow(1,integerVariable_+i,&one,false);
4715	rc.setGloballyValidAsInteger(2);
4716	globalCuts_.addCutIfNotDuplicate(rc) ;
4717	} else if (objUpper[i] > newCutoff) {
4718	n++;
4719	// add as global cut
4720	objUpper[i]=-COIN_DBL_MAX;
4721	OsiRowCut rc;
4722	rc.setLb(-COIN_DBL_MAX);
4723	rc.setUb(newUpper[i]);
4724	double one=1.0;
4725	rc.setRow(1,integerVariable_+i,&one,false);
4726	rc.setGloballyValidAsInteger(2);
4727	globalCuts_.addCutIfNotDuplicate(rc) ;
4728	}
4729	}
4730	if (newCutoff == -COIN_DBL_MAX) {
4731	COIN_DETAIL_PRINT(printf("Root analysis says finished\n"));
4732	} else if (n > numberFixedNow_) {
4733	COIN_DETAIL_PRINT(printf("%d more fixed by analysis - now %d\n", n - numberFixedNow_, n));
4734	numberFixedNow_ = n;
4735	}
4736	}
4737	if (eventHandler) {
4738	if (!eventHandler->event(CbcEventHandler::solution)) {
4739	eventHappened_ = true; // exit
4740	}
4741	newCutoff = getCutoff();
4742	}
4743	lockThread();
4744	/*
4745	Clean the tree to reflect the new solution, then see if the
4746	node comparison predicate wants to make any changes. If so,
4747	call setComparison for the side effect of rebuilding the heap.
4748	*/
4749	tree_->cleanTree(this,newCutoff,bestPossibleObjective_) ;
4750	if (nodeCompare_->newSolution(this) \|\|
4751	nodeCompare_->newSolution(this,continuousObjective_,
4752	continuousInfeasibilities_)) {
4753	tree_->setComparison(*nodeCompare_) ;
4754	}
4755	if (tree_->empty()) {
4756	continue;
4757	}
4758	unlockThread();
4759	}
4760	cutoff = getCutoff() ;
4761	/*
4762	Periodic activities: Opportunities to
4763	+ tweak the nodeCompare criteria,
4764	+ check if we've closed the integrality gap enough to quit,
4765	+ print a summary line to let the user know we're working
4766	*/
4767	if (numberNodes_ >= lastEvery1000) {
4768	lockThread();
4769	#ifdef COIN_HAS_CLP
4770	// See if we want dantzig row choice
4771	goToDantzig(1000, savePivotMethod);
4772	#endif
4773	lastEvery1000 = numberNodes_ + 1000;
4774	bool redoTree = nodeCompare_->every1000Nodes(this, numberNodes_) ;
4775	#ifdef CHECK_CUT_SIZE
4776	verifyCutSize (tree_, *this);
4777	#endif
4778	// redo tree if requested
4779	if (redoTree)
4780	tree_->setComparison(*nodeCompare_) ;
4781	unlockThread();
4782	}
4783	// Had hotstart before, now switched off
4784	if (saveCompare && !hotstartSolution_) {
4785	// hotstart switched off
4786	delete nodeCompare_; // off depth first
4787	nodeCompare_ = saveCompare;
4788	saveCompare = NULL;
4789	// redo tree
4790	lockThread();
4791	tree_->setComparison(*nodeCompare_) ;
4792	unlockThread();
4793	}
4794	if (numberNodes_ >= lastPrintEvery) {
4795	lastPrintEvery = numberNodes_ + printFrequency_;
4796	lockThread();
4797	int nNodes = tree_->size() ;
4798
4799	//MODIF PIERRE
4800	bestPossibleObjective_ = tree_->getBestPossibleObjective();
4801	#ifdef CBC_THREAD
4802	if (parallelMode() > 0 && master_) {
4803	// need to adjust for ones not on tree
4804	int numberThreads = master_->numberThreads();
4805	for (int i=0;i<numberThreads;i++) {
4806	CbcThread * child = master_->child(i);
4807	if (child->node()) {
4808	// adjust
4809	double value = child->node()->objectiveValue();
4810	bestPossibleObjective_ = CoinMin(bestPossibleObjective_, value);
4811	}
4812	}
4813	}
4814	#endif
4815	unlockThread();
4816	#if CBC_USEFUL_PRINTING>1
4817	if (getCutoff() < 1.0e20) {
4818	if (fabs(getCutoff() - (bestObjective_ - getCutoffIncrement())) > 1.0e-6 &&
4819	!parentModel_)
4820	printf("model cutoff in status %g, best %g, increment %g\n",
4821	getCutoff(), bestObjective_, getCutoffIncrement());
4822	assert (getCutoff() < bestObjective_ - getCutoffIncrement() +
4823	1.0e-6 + 1.0e-10*fabs(bestObjective_));
4824	}
4825	#endif
4826	if (!intParam_[CbcPrinting]) {
4827	// Parallel may not have any nodes
4828	if (!nNodes)
4829	bestPossibleObjective_ = lastBestPossibleObjective;
4830	else
4831	lastBestPossibleObjective = bestPossibleObjective_;
4832	messageHandler()->message(CBC_STATUS, messages())
4833	<< numberNodes_ << CoinMax(nNodes,1) << bestObjective_ << bestPossibleObjective_
4834	<< getCurrentSeconds()
4835	<< CoinMessageEol ;
4836	} else if (intParam_[CbcPrinting] == 1) {
4837	messageHandler()->message(CBC_STATUS2, messages())
4838	<< numberNodes_ << nNodes << bestObjective_ << bestPossibleObjective_
4839	<< tree_->lastDepth() << tree_->lastUnsatisfied()
4840	<< tree_->lastObjective() << numberIterations_
4841	<< getCurrentSeconds()
4842	<< CoinMessageEol ;
4843	} else if (!numberExtraIterations_) {
4844	messageHandler()->message(CBC_STATUS2, messages())
4845	<< numberNodes_ << nNodes << bestObjective_ << bestPossibleObjective_
4846	<< tree_->lastDepth() << tree_->lastUnsatisfied() << numberIterations_
4847	<< getCurrentSeconds()
4848	<< CoinMessageEol ;
4849	} else {
4850	messageHandler()->message(CBC_STATUS3, messages())
4851	<< numberNodes_ << numberFathoms_ << numberExtraNodes_ << nNodes
4852	<< bestObjective_ << bestPossibleObjective_
4853	<< tree_->lastDepth() << tree_->lastUnsatisfied() << numberIterations_ << numberExtraIterations_
4854	<< getCurrentSeconds()
4855	<< CoinMessageEol ;
4856	}
4857	#ifdef COIN_HAS_NTY
4858	if (symmetryInfo_)
4859	symmetryInfo_->statsOrbits(this,1);
4860	#endif
4861	#if PRINT_CONFLICT==1
4862	if (numberConflictCuts>lastNumberConflictCuts) {
4863	double length = lengthConflictCuts/numberConflictCuts;
4864	printf("%d new conflict cuts - total %d - average length %g\n",
4865	numberConflictCuts-lastNumberConflictCuts,
4866	numberConflictCuts,length);
4867	lastNumberConflictCuts = numberConflictCuts;
4868	}
4869	#endif
4870	if (eventHandler && !eventHandler->event(CbcEventHandler::treeStatus)) {
4871	eventHappened_ = true; // exit
4872	}
4873	}
4874	// See if can stop on gap
4875	if(canStopOnGap()) {
4876	stoppedOnGap_ = true ;
4877	}
4878
4879	#ifdef CHECK_NODE_FULL
4880	verifyTreeNodes(tree_, *this) ;
4881	# endif
4882	# ifdef CHECK_CUT_COUNTS
4883	verifyCutCounts(tree_, *this) ;
4884	# endif
4885	/*
4886	Now we come to the meat of the loop. To create the active subproblem, we'll
4887	pop the most promising node in the live set, rebuild the subproblem it
4888	represents, and then execute the current arm of the branch to create the
4889	active subproblem.
4890	*/
4891	CbcNode * node = NULL;
4892	#ifdef CBC_THREAD
4893	if (!parallelMode() \|\| parallelMode() == -1) {
4894	#endif
4895	node = tree_->bestNode(cutoff) ;
4896	// Possible one on tree worse than cutoff
4897	// Weird comparison function can leave ineligible nodes on tree
4898	if (!node \|\| node->objectiveValue() > cutoff)
4899	continue;
4900	// Do main work of solving node here
4901	doOneNode(this, node, createdNode);
4902	#ifdef JJF_ZERO
4903	if (node) {
4904	if (createdNode) {
4905	printf("Node %d depth %d, created %d depth %d\n",
4906	node->nodeNumber(), node->depth(),
4907	createdNode->nodeNumber(), createdNode->depth());
4908	} else {
4909	printf("Node %d depth %d, no created node\n",
4910	node->nodeNumber(), node->depth());
4911	}
4912	} else if (createdNode) {
4913	printf("Node exhausted, created %d depth %d\n",
4914	createdNode->nodeNumber(), createdNode->depth());
4915	} else {
4916	printf("Node exhausted, no created node\n");
4917	numberConsecutiveInfeasible = 2;
4918	}
4919	#endif
4920	//if (createdNode)
4921	//numberConsecutiveInfeasible=0;
4922	//else
4923	//numberConsecutiveInfeasible++;
4924	#ifdef CBC_THREAD
4925	} else if (parallelMode() > 0) {
4926	//lockThread();
4927	//node = tree_->bestNode(cutoff) ;
4928	// Possible one on tree worse than cutoff
4929	if (true \|\| !node \|\| node->objectiveValue() > cutoff) {
4930	assert (master_);
4931	if (master_) {
4932	int anyLeft = master_->waitForThreadsInTree(1);
4933	// may need to go round again
4934	if (anyLeft) {
4935	continue;
4936	} else {
4937	master_->stopThreads(-1);
4938	}
4939	}
4940	}
4941	//unlockThread();
4942	} else {
4943	// Deterministic parallel
4944	if ((tree_->size() < CoinMax(numberThreads_, 8)\|\|
4945	hotstartSolution_) && !goneParallel) {
4946	node = tree_->bestNode(cutoff) ;
4947	// Possible one on tree worse than cutoff
4948	if (!node \|\| node->objectiveValue() > cutoff)
4949	continue;
4950	// Do main work of solving node here
4951	doOneNode(this, node, createdNode);
4952	assert (createdNode);
4953	if (!createdNode->active()) {
4954	delete createdNode;
4955	createdNode = NULL;
4956	} else {
4957	// Say one more pointing to this
4958	node->nodeInfo()->increment() ;
4959	tree_->push(createdNode) ;
4960	}
4961	if (node->active()) {
4962	assert (node->nodeInfo());
4963	if (node->nodeInfo()->numberBranchesLeft()) {
4964	tree_->push(node) ;
4965	} else {
4966	node->setActive(false);
4967	}
4968	} else {
4969	if (node->nodeInfo()) {
4970	if (!node->nodeInfo()->numberBranchesLeft())
4971	node->nodeInfo()->allBranchesGone(); // can clean up
4972	// So will delete underlying stuff
4973	node->setActive(true);
4974	}
4975	delNode[nDeleteNode++] = node;
4976	node = NULL;
4977	}
4978	if (nDeleteNode >= MAX_DEL_NODE) {
4979	for (int i = 0; i < nDeleteNode; i++) {
4980	//printf("trying to del %d %x\n",i,delNode[i]);
4981	delete delNode[i];
4982	//printf("done to del %d %x\n",i,delNode[i]);
4983	}
4984	nDeleteNode = 0;
4985	}
4986	} else {
4987	// Split and solve
4988	master_->deterministicParallel();
4989	goneParallel = true;
4990	}
4991	}
4992	#endif
4993	}
4994	if (nDeleteNode) {
4995	for (int i = 0; i < nDeleteNode; i++) {
4996	delete delNode[i];
4997	}
4998	nDeleteNode = 0;
4999	}
5000	#ifdef CBC_THREAD
5001	if (master_) {
5002	master_->stopThreads(-1);
5003	master_->waitForThreadsInTree(2);
5004	// adjust time to allow for children on some systems
5005	//dblParam_[CbcStartSeconds] -= CoinCpuTimeJustChildren();
5006	}
5007	#endif
5008	/*
5009	End of the non-abort actions. The next block of code is executed if we've
5010	aborted because we hit one of the limits. Clean up by deleting the live set
5011	and break out of the node processing loop. Note that on an abort, node may
5012	have been pushed back onto the tree for further processing, in which case
5013	it'll be deleted in cleanTree. We need to check.
5014	*/
5015	if (!(numberNodes_ < intParam_[CbcMaxNumNode] &&
5016	numberSolutions_ < intParam_[CbcMaxNumSol] &&
5017	!maximumSecondsReached() &&
5018	!stoppedOnGap_ &&
5019	!eventHappened_ &&
5020	(maximumNumberIterations_ < 0 \|\| numberIterations_ < maximumNumberIterations_))
5021	) {
5022	if (tree_->size()) {
5023	double dummyBest;
5024	tree_->cleanTree(this, -COIN_DBL_MAX, dummyBest) ;
5025	#if 0 // Does not seem to be needed def CBC_THREAD
5026	if (parallelMode() > 0 && master_) {
5027	// see if any dangling nodes
5028	int numberThreads = master_->numberThreads();
5029	for (int i=0;i<numberThreads;i++) {
5030	CbcThread * child = master_->child(i);
5031	//if (child->createdNode())
5032	//printf("CHILD_NODE %p\n",child->createdNode());
5033	delete child->createdNode();
5034	}
5035	}
5036	#endif
5037	}
5038	delete nextRowCut_;
5039	/* order is important here:
5040	* maximumSecondsReached() should be checked before eventHappened_ and
5041	* isNodeLimitReached() should be checked after eventHappened_
5042	* reason is, that at timelimit, eventHappened_ is set to true to make Cbc stop fast
5043	* and if Ctrl+C is hit, then the nodelimit is set to -1 to make Cbc stop
5044	*/
5045	if (stoppedOnGap_) {
5046	messageHandler()->message(CBC_GAP, messages())
5047	<< bestObjective_ - bestPossibleObjective_
5048	<< dblParam_[CbcAllowableGap]
5049	<< dblParam_[CbcAllowableFractionGap]*100.0
5050	<< CoinMessageEol ;
5051	secondaryStatus_ = 2;
5052	status_ = 0 ;
5053	} else if (maximumSecondsReached()) {
5054	handler_->message(CBC_MAXTIME, messages_) << CoinMessageEol ;
5055	secondaryStatus_ = 4;
5056	status_ = 1 ;
5057	} else if (numberSolutions_ >= intParam_[CbcMaxNumSol]) {
5058	handler_->message(CBC_MAXSOLS, messages_) << CoinMessageEol ;
5059	secondaryStatus_ = 6;
5060	status_ = 1 ;
5061	} else if (isNodeLimitReached()) {
5062	handler_->message(CBC_MAXNODES, messages_) << CoinMessageEol ;
5063	secondaryStatus_ = 3;
5064	status_ = 1 ;
5065	} else if (maximumNumberIterations_ >= 0 && numberIterations_ >= maximumNumberIterations_) {
5066	handler_->message(CBC_MAXITERS, messages_) << CoinMessageEol ;
5067	secondaryStatus_ = 8;
5068	status_ = 1 ;
5069	} else {
5070	handler_->message(CBC_EVENT, messages_) << CoinMessageEol ;
5071	secondaryStatus_ = 5;
5072	status_ = 5 ;
5073	}
5074	}
5075	#ifdef CBC_THREAD
5076	if (master_) {
5077	delete master_;
5078	master_ = NULL;
5079	masterThread_ = NULL;
5080	}
5081	#endif
5082	/*
5083	That's it, we've exhausted the search tree, or broken out of the loop because
5084	we hit some limit on evaluation.
5085
5086	We may have got an intelligent tree so give it one more chance
5087	*/
5088	// Tell solver we are not in Branch and Cut
5089	solver_->setHintParam(OsiDoInBranchAndCut, false, OsiHintDo, NULL) ;
5090	tree_->endSearch();
5091	// If we did any sub trees - did we give up on any?
5092	if ( numberStoppedSubTrees_)
5093	status_ = 1;
5094	numberNodes_ += numberExtraNodes_;
5095	numberIterations_ += numberExtraIterations_;
5096	if (eventHandler) {
5097	eventHandler->event(CbcEventHandler::endSearch);
5098	}
5099	if (!status_) {
5100	// Set best possible unless stopped on gap
5101	if (secondaryStatus_ != 2)
5102	bestPossibleObjective_ = bestObjective_;
5103	handler_->message(CBC_END_GOOD, messages_)
5104	<< bestObjective_ << numberIterations_ << numberNodes_ << getCurrentSeconds()
5105	<< CoinMessageEol ;
5106	} else {
5107	handler_->message(CBC_END, messages_)
5108	<< bestObjective_ << bestPossibleObjective_
5109	<< numberIterations_ << numberNodes_ << getCurrentSeconds()
5110	<< CoinMessageEol ;
5111	}
5112	if ((moreSpecialOptions_&4194304)!=0) {
5113	// Conflict cuts
5114	int numberCuts = globalCuts_.sizeRowCuts();
5115	int nConflict=0;
5116	double sizeConflict = 0.0;
5117	for (int i=0;i<numberCuts;i++) {
5118	OsiRowCut2 * cut = globalCuts_.cut(i);
5119	if (cut->whichRow()==1) {
5120	nConflict++;
5121	sizeConflict += cut->row().getNumElements();
5122	}
5123	}
5124	if (nConflict) {
5125	sizeConflict /= nConflict;
5126	char general[200];
5127	sprintf(general, "%d conflict cuts generated - average length %g",
5128	nConflict,sizeConflict);
5129	messageHandler()->message(CBC_GENERAL,
5130	messages())
5131	<< general << CoinMessageEol ;
5132	}
5133	}
5134	if (numberStrongIterations_)
5135	handler_->message(CBC_STRONG_STATS, messages_)
5136	<< strongInfo_[0] << numberStrongIterations_ << strongInfo_[2]
5137	<< strongInfo_[1] << CoinMessageEol ;
5138	if (!numberExtraNodes_)
5139	handler_->message(CBC_OTHER_STATS, messages_)
5140	<< maximumDepthActual_
5141	<< numberDJFixed_ << CoinMessageEol ;
5142	else
5143	handler_->message(CBC_OTHER_STATS2, messages_)
5144	<< maximumDepthActual_
5145	<< numberDJFixed_ << numberFathoms_ << numberExtraNodes_ << numberExtraIterations_
5146	<< CoinMessageEol ;
5147	#ifdef COIN_HAS_NTY
5148	if (symmetryInfo_)
5149	symmetryInfo_->statsOrbits(this,1);
5150	#endif
5151	if (doStatistics == 100) {
5152	for (int i = 0; i < numberObjects_; i++) {
5153	CbcSimpleIntegerDynamicPseudoCost * obj =
5154	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object_[i]) ;
5155	if (obj)
5156	obj->print();
5157	}
5158	}
5159	if (statistics_) {
5160	// report in some way
5161	int * lookup = new int[numberObjects_];
5162	int i;
5163	for (i = 0; i < numberObjects_; i++)
5164	lookup[i] = -1;
5165	bool goodIds = false; //true;
5166	for (i = 0; i < numberObjects_; i++) {
5167	int iColumn = object_[i]->columnNumber();
5168	if (iColumn >= 0 && iColumn < numberColumns) {
5169	if (lookup[i] == -1) {
5170	lookup[i] = iColumn;
5171	} else {
5172	goodIds = false;
5173	break;
5174	}
5175	} else {
5176	goodIds = false;
5177	break;
5178	}
5179	}
5180	if (!goodIds) {
5181	delete [] lookup;
5182	lookup = NULL;
5183	}
5184	if (doStatistics >= 3) {
5185	printf(" node parent depth column value obj inf\n");
5186	for ( i = 0; i < numberNodes2_; i++) {
5187	statistics_[i]->print(lookup);
5188	}
5189	}
5190	if (doStatistics > 1) {
5191	// Find last solution
5192	int k;
5193	for (k = numberNodes2_ - 1; k >= 0; k--) {
5194	if (statistics_[k]->endingObjective() != COIN_DBL_MAX &&
5195	!statistics_[k]->endingInfeasibility())
5196	break;
5197	}
5198	if (k >= 0) {
5199	int depth = statistics_[k]->depth();
5200	int * which = new int[depth+1];
5201	for (i = depth; i >= 0; i--) {
5202	which[i] = k;
5203	k = statistics_[k]->parentNode();
5204	}
5205	printf(" node parent depth column value obj inf\n");
5206	for (i = 0; i <= depth; i++) {
5207	statistics_[which[i]]->print(lookup);
5208	}
5209	delete [] which;
5210	}
5211	}
5212	// now summary
5213	int maxDepth = 0;
5214	double averageSolutionDepth = 0.0;
5215	int numberSolutions = 0;
5216	double averageCutoffDepth = 0.0;
5217	double averageSolvedDepth = 0.0;
5218	int numberCutoff =