Context Navigation

source: trunk/Cbc/src/CbcNode.cpp @ 1650

Last change on this file since 1650 was 1650, checked in by forrest, 9 years ago
messages for fathoming plus soem attempts at robustness
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 197.5 KB

Line
1	/* $Id: CbcNode.cpp 1650 2011-05-10 12:33:58Z 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 CBC_CHECK_BASIS
18	#include <cassert>
19	#include <cfloat>
20	#define CUTS
21	#include "OsiSolverInterface.hpp"
22	#include "OsiChooseVariable.hpp"
23	#include "OsiAuxInfo.hpp"
24	#include "OsiSolverBranch.hpp"
25	#include "CoinWarmStartBasis.hpp"
26	#include "CoinTime.hpp"
27	#include "CbcModel.hpp"
28	#include "CbcNode.hpp"
29	#include "CbcStatistics.hpp"
30	#include "CbcStrategy.hpp"
31	#include "CbcBranchActual.hpp"
32	#include "CbcBranchDynamic.hpp"
33	#include "OsiRowCut.hpp"
34	#include "OsiRowCutDebugger.hpp"
35	#include "OsiCuts.hpp"
36	#include "CbcCountRowCut.hpp"
37	#include "CbcFeasibilityBase.hpp"
38	#include "CbcMessage.hpp"
39	#ifdef COIN_HAS_CLP
40	#include "OsiClpSolverInterface.hpp"
41	#include "ClpSimplexOther.hpp"
42	#endif
43	using namespace std;
44	#include "CglCutGenerator.hpp"
45
46
47	CbcNode::CbcNode() :
48	nodeInfo_(NULL),
49	objectiveValue_(1.0e100),
50	guessedObjectiveValue_(1.0e100),
51	sumInfeasibilities_(0.0),
52	branch_(NULL),
53	depth_(-1),
54	numberUnsatisfied_(0),
55	nodeNumber_(-1),
56	state_(0)
57	{
58	#ifdef CHECK_NODE
59	printf("CbcNode %x Constructor\n", this);
60	#endif
61	}
62	// Print
63	void
64	CbcNode::print() const
65	{
66	printf("number %d obj %g depth %d sumun %g nunsat %d state %d\n",
67	nodeNumber_, objectiveValue_, depth_, sumInfeasibilities_, numberUnsatisfied_, state_);
68	}
69	CbcNode::CbcNode(CbcModel * model,
70	CbcNode * lastNode) :
71	nodeInfo_(NULL),
72	objectiveValue_(1.0e100),
73	guessedObjectiveValue_(1.0e100),
74	sumInfeasibilities_(0.0),
75	branch_(NULL),
76	depth_(-1),
77	numberUnsatisfied_(0),
78	nodeNumber_(-1),
79	state_(0)
80	{
81	#ifdef CHECK_NODE
82	printf("CbcNode %x Constructor from model\n", this);
83	#endif
84	model->setObjectiveValue(this, lastNode);
85
86	if (lastNode) {
87	if (lastNode->nodeInfo_) {
88	lastNode->nodeInfo_->increment();
89	}
90	}
91	nodeNumber_ = model->getNodeCount();
92	}
93
94	#define CBC_NEW_CREATEINFO
95	#ifdef CBC_NEW_CREATEINFO
96
97	/*
98	New createInfo, with basis manipulation hidden inside mergeBasis. Allows
99	solvers to override and carry over all information from one basis to
100	another.
101	*/
102
103	void
104	CbcNode::createInfo (CbcModel *model,
105	CbcNode *lastNode,
106	const CoinWarmStartBasis *lastws,
107	const double lastLower, const double lastUpper,
108	int numberOldActiveCuts, int numberNewCuts)
109
110	{
111	OsiSolverInterface *solver = model->solver();
112	CbcStrategy *strategy = model->strategy();
113	/*
114	The root --- no parent. Create full basis and bounds information.
115	*/
116	if (!lastNode) {
117	if (!strategy)
118	nodeInfo_ = new CbcFullNodeInfo(model, solver->getNumRows());
119	else
120	nodeInfo_ = strategy->fullNodeInfo(model, solver->getNumRows());
121	} else {
122	/*
123	Not the root. Create an edit from the parent's basis & bound information.
124	This is not quite as straightforward as it seems. We need to reintroduce
125	cuts we may have dropped out of the basis, in the correct position, because
126	this whole process is strictly positional. Start by grabbing the current
127	basis.
128	*/
129	bool mustDeleteBasis;
130	const CoinWarmStartBasis *ws =
131	dynamic_cast<const CoinWarmStartBasis*>(solver->getPointerToWarmStart(mustDeleteBasis));
132	assert(ws != NULL); // make sure not volume
133	//int numberArtificials = lastws->getNumArtificial();
134	int numberColumns = solver->getNumCols();
135	int numberRowsAtContinuous = model->numberRowsAtContinuous();
136	int currentNumberCuts = model->currentNumberCuts();
137	# ifdef CBC_CHECK_BASIS
138	std::cout
139	<< "Before expansion: orig " << numberRowsAtContinuous
140	<< ", old " << numberOldActiveCuts
141	<< ", new " << numberNewCuts
142	<< ", current " << currentNumberCuts << "." << std::endl ;
143	ws->print();
144	# endif
145	/*
146	Clone the basis and resize it to hold the structural constraints, plus
147	all the cuts: old cuts, both active and inactive (currentNumberCuts),
148	and new cuts (numberNewCuts). This will become the expanded basis.
149	*/
150	CoinWarmStartBasis *expanded =
151	dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ;
152	int iCompact = numberRowsAtContinuous + numberOldActiveCuts + numberNewCuts ;
153	// int nPartial = numberRowsAtContinuous+currentNumberCuts;
154	int iFull = numberRowsAtContinuous + currentNumberCuts + numberNewCuts;
155	// int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4);
156	// printf("l %d full %d\n",maxBasisLength,iFull);
157	expanded->resize(iFull, numberColumns);
158	# ifdef CBC_CHECK_BASIS
159	std::cout
160	<< "\tFull basis " << iFull << " rows, "
161	<< numberColumns << " columns; compact "
162	<< iCompact << " rows." << std::endl ;
163	# endif
164	/*
165	Now flesh out the expanded basis. The clone already has the
166	correct status information for the variables and for the structural
167	(numberRowsAtContinuous) constraints. Any indices beyond nPartial must be
168	cuts created while processing this node --- they can be copied en bloc
169	into the correct position in the expanded basis. The space reserved for
170	xferRows is a gross overestimate.
171	*/
172	CoinWarmStartBasis::XferVec xferRows ;
173	xferRows.reserve(iFull - numberRowsAtContinuous + 1) ;
174	if (numberNewCuts) {
175	xferRows.push_back(
176	CoinWarmStartBasis::XferEntry(iCompact - numberNewCuts,
177	iFull - numberNewCuts, numberNewCuts)) ;
178	}
179	/*
180	From nPartial down, record the entries we want to copy from the current
181	basis (the entries for the active cuts; non-zero in the list returned
182	by addedCuts). Fill the expanded basis with entries showing a status of
183	basic for the deactivated (loose) cuts.
184	*/
185	CbcCountRowCut **cut = model->addedCuts();
186	iFull -= (numberNewCuts + 1) ;
187	iCompact -= (numberNewCuts + 1) ;
188	int runLen = 0 ;
189	CoinWarmStartBasis::XferEntry entry(-1, -1, -1) ;
190	while (iFull >= numberRowsAtContinuous) {
191	for ( ; iFull >= numberRowsAtContinuous &&
192	cut[iFull-numberRowsAtContinuous] ; iFull--)
193	runLen++ ;
194	if (runLen) {
195	iCompact -= runLen ;
196	entry.first = iCompact + 1 ;
197	entry.second = iFull + 1 ;
198	entry.third = runLen ;
199	runLen = 0 ;
200	xferRows.push_back(entry) ;
201	}
202	for ( ; iFull >= numberRowsAtContinuous &&
203	!cut[iFull-numberRowsAtContinuous] ; iFull--)
204	expanded->setArtifStatus(iFull, CoinWarmStartBasis::basic);
205	}
206	/*
207	Finally, call mergeBasis to copy over entries from the current basis to
208	the expanded basis. Since we cloned the expanded basis from the active basis
209	and haven't changed the number of variables, only row status entries need
210	to be copied.
211	*/
212	expanded->mergeBasis(ws, &xferRows, 0) ;
213
214	#ifdef CBC_CHECK_BASIS
215	std::cout << "Expanded basis:" << std::endl ;
216	expanded->print() ;
217	std::cout << "Diffing against:" << std::endl ;
218	lastws->print() ;
219	#endif
220	assert (expanded->getNumArtificial() >= lastws->getNumArtificial());
221	#ifdef CLP_INVESTIGATE
222	if (!expanded->fullBasis()) {
223	int iFull = numberRowsAtContinuous + currentNumberCuts + numberNewCuts;
224	printf("cont %d old %d new %d current %d full inc %d full %d\n",
225	numberRowsAtContinuous, numberOldActiveCuts, numberNewCuts,
226	currentNumberCuts, iFull, iFull - numberNewCuts);
227	}
228	#endif
229
230	/*
231	Now that we have two bases in proper positional correspondence, creating
232	the actual diff is dead easy.
233
234	Note that we're going to compare the expanded basis here to the stripped
235	basis (lastws) produced by addCuts. It doesn't affect the correctness (the
236	diff process has no knowledge of the meaning of an entry) but it does
237	mean that we'll always generate a whack of diff entries because the expanded
238	basis is considerably larger than the stripped basis.
239	*/
240	CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ;
241
242	/*
243	Diff the bound vectors. It's assumed the number of structural variables
244	is not changing. For branching objects that change bounds on integer
245	variables, we should see at least one bound change as a consequence
246	of applying the branch that generated this subproblem from its parent.
247	This need not hold for other types of branching objects (hyperplane
248	branches, for example).
249	*/
250	const double * lower = solver->getColLower();
251	const double * upper = solver->getColUpper();
252
253	double boundChanges = new double [2numberColumns] ;
254	int variables = new int [2numberColumns] ;
255	int numberChangedBounds = 0;
256
257	int i;
258	for (i = 0; i < numberColumns; i++) {
259	if (lower[i] != lastLower[i]) {
260	variables[numberChangedBounds] = i;
261	boundChanges[numberChangedBounds++] = lower[i];
262	}
263	if (upper[i] != lastUpper[i]) {
264	variables[numberChangedBounds] = i \| 0x80000000;
265	boundChanges[numberChangedBounds++] = upper[i];
266	}
267	#ifdef CBC_DEBUG
268	if (lower[i] != lastLower[i]) {
269	std::cout
270	<< "lower on " << i << " changed from "
271	<< lastLower[i] << " to " << lower[i] << std::endl ;
272	}
273	if (upper[i] != lastUpper[i]) {
274	std::cout
275	<< "upper on " << i << " changed from "
276	<< lastUpper[i] << " to " << upper[i] << std::endl ;
277	}
278	#endif
279	}
280	#ifdef CBC_DEBUG
281	std::cout << numberChangedBounds << " changed bounds." << std::endl ;
282	#endif
283	//if (lastNode->branchingObject()->boundBranch())
284	//assert (numberChangedBounds);
285	/*
286	Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and
287	return.
288	*/
289	if (!strategy)
290	nodeInfo_ =
291	new CbcPartialNodeInfo(lastNode->nodeInfo_, this, numberChangedBounds,
292	variables, boundChanges, basisDiff) ;
293	else
294	nodeInfo_ =
295	strategy->partialNodeInfo(model, lastNode->nodeInfo_, this,
296	numberChangedBounds, variables, boundChanges,
297	basisDiff) ;
298	delete basisDiff ;
299	delete [] boundChanges;
300	delete [] variables;
301	delete expanded ;
302	if (mustDeleteBasis)
303	delete ws;
304	}
305	// Set node number
306	nodeInfo_->setNodeNumber(model->getNodeCount2());
307	state_ \|= 2; // say active
308	}
309
310	#else // CBC_NEW_CREATEINFO
311
312	/*
313	Original createInfo, with bare manipulation of basis vectors. Fails if solver
314	maintains additional information in basis.
315	*/
316
317	void
318	CbcNode::createInfo (CbcModel *model,
319	CbcNode *lastNode,
320	const CoinWarmStartBasis *lastws,
321	const double lastLower, const double lastUpper,
322	int numberOldActiveCuts, int numberNewCuts)
323	{
324	OsiSolverInterface * solver = model->solver();
325	CbcStrategy * strategy = model->strategy();
326	/*
327	The root --- no parent. Create full basis and bounds information.
328	*/
329	if (!lastNode) {
330	if (!strategy)
331	nodeInfo_ = new CbcFullNodeInfo(model, solver->getNumRows());
332	else
333	nodeInfo_ = strategy->fullNodeInfo(model, solver->getNumRows());
334	}
335	/*
336	Not the root. Create an edit from the parent's basis & bound information.
337	This is not quite as straightforward as it seems. We need to reintroduce
338	cuts we may have dropped out of the basis, in the correct position, because
339	this whole process is strictly positional. Start by grabbing the current
340	basis.
341	*/
342	else {
343	bool mustDeleteBasis;
344	const CoinWarmStartBasis* ws =
345	dynamic_cast<const CoinWarmStartBasis*>(solver->getPointerToWarmStart(mustDeleteBasis));
346	assert(ws != NULL); // make sure not volume
347	//int numberArtificials = lastws->getNumArtificial();
348	int numberColumns = solver->getNumCols();
349
350	const double * lower = solver->getColLower();
351	const double * upper = solver->getColUpper();
352
353	int i;
354	/*
355	Create a clone and resize it to hold all the structural constraints, plus
356	all the cuts: old cuts, both active and inactive (currentNumberCuts), and
357	new cuts (numberNewCuts).
358
359	TODO: You'd think that the set of constraints (logicals) in the expanded
360	basis should match the set represented in lastws. At least, that's
361	what I thought. But at the point I first looked hard at this bit of
362	code, it turned out that lastws was the stripped basis produced at
363	the end of addCuts(), rather than the raw basis handed back by
364	addCuts1(). The expanded basis here is equivalent to the raw basis of
365	addCuts1(). I said ``whoa, that's not good, I must have introduced a
366	bug'' and went back to John's code to see where I'd gone wrong.
367	And discovered the same `error' in his code.
368
369	After a bit of thought, my conclusion is that correctness is not
370	affected by whether lastws is the stripped or raw basis. The diffs
371	have no semantics --- just a set of changes that need to be made
372	to convert lastws into expanded. I think the only effect is that we
373	store a lot more diffs (everything in expanded that's not covered by
374	the stripped basis). But I need to give this more thought. There
375	may well be some subtle error cases.
376
377	In the mean time, I've twiddled addCuts() to set lastws to the raw
378	basis. Makes me (Lou) less nervous to compare apples to apples.
379	*/
380	CoinWarmStartBasis *expanded =
381	dynamic_cast<CoinWarmStartBasis *>(ws->clone()) ;
382	int numberRowsAtContinuous = model->numberRowsAtContinuous();
383	int iFull = numberRowsAtContinuous + model->currentNumberCuts() +
384	numberNewCuts;
385	//int numberArtificialsNow = iFull;
386	//int maxBasisLength = ((iFull+15)>>4)+((numberColumns+15)>>4);
387	//printf("l %d full %d\n",maxBasisLength,iFull);
388	if (expanded)
389	expanded->resize(iFull, numberColumns);
390	#ifdef CBC_CHECK_BASIS
391	printf("Before expansion: orig %d, old %d, new %d, current %d\n",
392	numberRowsAtContinuous, numberOldActiveCuts, numberNewCuts,
393	model->currentNumberCuts()) ;
394	ws->print();
395	#endif
396	/*
397	Now fill in the expanded basis. Any indices beyond nPartial must
398	be cuts created while processing this node --- they can be copied directly
399	into the expanded basis. From nPartial down, pull the status of active cuts
400	from ws, interleaving with a B entry for the deactivated (loose) cuts.
401	*/
402	int numberDropped = model->currentNumberCuts() - numberOldActiveCuts;
403	int iCompact = iFull - numberDropped;
404	CbcCountRowCut ** cut = model->addedCuts();
405	int nPartial = model->currentNumberCuts() + numberRowsAtContinuous;
406	iFull--;
407	for (; iFull >= nPartial; iFull--) {
408	CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact);
409	//assert (status != CoinWarmStartBasis::basic); // may be permanent cut
410	expanded->setArtifStatus(iFull, status);
411	}
412	for (; iFull >= numberRowsAtContinuous; iFull--) {
413	if (cut[iFull-numberRowsAtContinuous]) {
414	CoinWarmStartBasis::Status status = ws->getArtifStatus(--iCompact);
415	// If no cut generator being used then we may have basic variables
416	//if (model->getMaximumCutPasses()&&
417	// status == CoinWarmStartBasis::basic)
418	//printf("cut basic\n");
419	expanded->setArtifStatus(iFull, status);
420	} else {
421	expanded->setArtifStatus(iFull, CoinWarmStartBasis::basic);
422	}
423	}
424	#ifdef CBC_CHECK_BASIS
425	printf("Expanded basis\n");
426	expanded->print() ;
427	printf("Diffing against\n") ;
428	lastws->print() ;
429	#endif
430	/*
431	Now that we have two bases in proper positional correspondence, creating
432	the actual diff is dead easy.
433	*/
434
435	CoinWarmStartDiff *basisDiff = expanded->generateDiff(lastws) ;
436	/*
437	Diff the bound vectors. It's assumed the number of structural variables is
438	not changing. Assuming that branching objects all involve integer variables,
439	we should see at least one bound change as a consequence of processing this
440	subproblem. Different types of branching objects could break this assertion.
441	Not true at all - we have not applied current branch - JJF.
442	*/
443	double boundChanges = new double [2numberColumns] ;
444	int variables = new int [2numberColumns] ;
445	int numberChangedBounds = 0;
446	for (i = 0; i < numberColumns; i++) {
447	if (lower[i] != lastLower[i]) {
448	variables[numberChangedBounds] = i;
449	boundChanges[numberChangedBounds++] = lower[i];
450	}
451	if (upper[i] != lastUpper[i]) {
452	variables[numberChangedBounds] = i \| 0x80000000;
453	boundChanges[numberChangedBounds++] = upper[i];
454	}
455	#ifdef CBC_DEBUG
456	if (lower[i] != lastLower[i])
457	printf("lower on %d changed from %g to %g\n",
458	i, lastLower[i], lower[i]);
459	if (upper[i] != lastUpper[i])
460	printf("upper on %d changed from %g to %g\n",
461	i, lastUpper[i], upper[i]);
462	#endif
463	}
464	#ifdef CBC_DEBUG
465	printf("%d changed bounds\n", numberChangedBounds) ;
466	#endif
467	//if (lastNode->branchingObject()->boundBranch())
468	//assert (numberChangedBounds);
469	/*
470	Hand the lot over to the CbcPartialNodeInfo constructor, then clean up and
471	return.
472	*/
473	if (!strategy)
474	nodeInfo_ =
475	new CbcPartialNodeInfo(lastNode->nodeInfo_, this, numberChangedBounds,
476	variables, boundChanges, basisDiff) ;
477	else
478	nodeInfo_ = strategy->partialNodeInfo(model, lastNode->nodeInfo_, this, numberChangedBounds,
479	variables, boundChanges, basisDiff) ;
480	delete basisDiff ;
481	delete [] boundChanges;
482	delete [] variables;
483	delete expanded ;
484	if (mustDeleteBasis)
485	delete ws;
486	}
487	// Set node number
488	nodeInfo_->setNodeNumber(model->getNodeCount2());
489	state_ \|= 2; // say active
490	}
491
492	#endif // CBC_NEW_CREATEINFO
493	/*
494	The routine scans through the object list of the model looking for objects
495	that indicate infeasibility. It tests each object using strong branching
496	and selects the one with the least objective degradation. A corresponding
497	branching object is left attached to lastNode.
498
499	If strong branching is disabled, a candidate object is chosen essentially
500	at random (whatever object ends up in pos'n 0 of the candidate array).
501
502	If a branching candidate is found to be monotone, bounds are set to fix the
503	variable and the routine immediately returns (the caller is expected to
504	reoptimize).
505
506	If a branching candidate is found to result in infeasibility in both
507	directions, the routine immediately returns an indication of infeasibility.
508
509	Returns: 0 both branch directions are feasible
510	-1 branching variable is monotone
511	-2 infeasible
512
513	Original comments:
514	Here could go cuts etc etc
515	For now just fix on objective from strong branching.
516	*/
517
518	int CbcNode::chooseBranch (CbcModel model, CbcNode lastNode, int numberPassesLeft)
519
520	{
521	if (lastNode)
522	depth_ = lastNode->depth_ + 1;
523	else
524	depth_ = 0;
525	delete branch_;
526	branch_ = NULL;
527	OsiSolverInterface * solver = model->solver();
528	# ifdef COIN_HAS_CLP
529	OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
530	int saveClpOptions = 0;
531	if (osiclp) {
532	// for faster hot start
533	saveClpOptions = osiclp->specialOptions();
534	osiclp->setSpecialOptions(saveClpOptions \| 8192);
535	}
536	# else
537	OsiSolverInterface *osiclp = NULL ;
538	# endif
539	double saveObjectiveValue = solver->getObjValue();
540	double objectiveValue = CoinMax(solver->getObjSense() * saveObjectiveValue, objectiveValue_);
541	const double * lower = solver->getColLower();
542	const double * upper = solver->getColUpper();
543	// See what user thinks
544	int anyAction = model->problemFeasibility()->feasible(model, 0);
545	if (anyAction) {
546	// will return -2 if infeasible , 0 if treat as integer
547	return anyAction - 1;
548	}
549	double integerTolerance =
550	model->getDblParam(CbcModel::CbcIntegerTolerance);
551	// point to useful information
552	OsiBranchingInformation usefulInfo = model->usefulInformation();
553	// and modify
554	usefulInfo.depth_ = depth_;
555	int i;
556	bool beforeSolution = model->getSolutionCount() == 0;
557	int numberStrong = model->numberStrong();
558	// switch off strong if hotstart
559	const double * hotstartSolution = model->hotstartSolution();
560	const int * hotstartPriorities = model->hotstartPriorities();
561	int numberObjects = model->numberObjects();
562	int numberColumns = model->getNumCols();
563	double * saveUpper = new double[numberColumns];
564	double * saveLower = new double[numberColumns];
565	for (i = 0; i < numberColumns; i++) {
566	saveLower[i] = lower[i];
567	saveUpper[i] = upper[i];
568	}
569
570	// Save solution in case heuristics need good solution later
571
572	double * saveSolution = new double[numberColumns];
573	memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
574	model->reserveCurrentSolution(saveSolution);
575	if (hotstartSolution) {
576	numberStrong = 0;
577	if ((model->moreSpecialOptions()&1024) != 0) {
578	int nBad = 0;
579	int nUnsat = 0;
580	int nDiff = 0;
581	for (int i = 0; i < numberObjects; i++) {
582	OsiObject * object = model->modifiableObject(i);
583	const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
584	if (thisOne) {
585	int iColumn = thisOne->columnNumber();
586	double targetValue = hotstartSolution[iColumn];
587	double value = saveSolution[iColumn];
588	if (fabs(value - floor(value + 0.5)) > 1.0e-6) {
589	nUnsat++;
590	#ifdef CLP_INVESTIGATE
591	printf("H %d is %g target %g\n", iColumn, value, targetValue);
592	#endif
593	} else if (fabs(targetValue - value) > 1.0e-6) {
594	nDiff++;
595	}
596	if (targetValue < saveLower[iColumn] \|\|
597	targetValue > saveUpper[iColumn]) {
598	#ifdef CLP_INVESTIGATE
599	printf("%d has target %g and current bounds %g and %g\n",
600	iColumn, targetValue, saveLower[iColumn], saveUpper[iColumn]);
601	#endif
602	nBad++;
603	}
604	}
605	}
606	#ifdef CLP_INVESTIGATE
607	printf("Hot %d unsatisfied, %d outside limits, %d different\n",
608	nUnsat, nBad, nDiff);
609	#endif
610	if (nBad) {
611	// switch off as not possible
612	hotstartSolution = NULL;
613	model->setHotstartSolution(NULL, NULL);
614	usefulInfo.hotstartSolution_ = NULL;
615	}
616	}
617	}
618	int numberStrongDone = 0;
619	int numberUnfinished = 0;
620	int numberStrongInfeasible = 0;
621	int numberStrongIterations = 0;
622	int saveNumberStrong = numberStrong;
623	bool checkFeasibility = numberObjects > model->numberIntegers();
624	int maximumStrong = CoinMax(CoinMin(numberStrong, numberObjects), 1);
625	/*
626	Get a branching decision object. Use the default decision criteria unless
627	the user has loaded a decision method into the model.
628	*/
629	CbcBranchDecision *decision = model->branchingMethod();
630	CbcDynamicPseudoCostBranchingObject * dynamicBranchingObject =
631	dynamic_cast<CbcDynamicPseudoCostBranchingObject *>(decision);
632	if (!decision \|\| dynamicBranchingObject)
633	decision = new CbcBranchDefaultDecision();
634	decision->initialize(model);
635	CbcStrongInfo * choice = new CbcStrongInfo[maximumStrong];
636	// May go round twice if strong branching fixes all local candidates
637	bool finished = false;
638	double estimatedDegradation = 0.0;
639	while (!finished) {
640	finished = true;
641	// Some objects may compute an estimate of best solution from here
642	estimatedDegradation = 0.0;
643	//int numberIntegerInfeasibilities=0; // without odd ones
644	numberStrongDone = 0;
645	numberUnfinished = 0;
646	numberStrongInfeasible = 0;
647	numberStrongIterations = 0;
648
649	// We may go round this loop twice (only if we think we have solution)
650	for (int iPass = 0; iPass < 2; iPass++) {
651
652	// compute current state
653	//int numberObjectInfeasibilities; // just odd ones
654	//model->feasibleSolution(
655	// numberIntegerInfeasibilities,
656	// numberObjectInfeasibilities);
657	// Some objects may compute an estimate of best solution from here
658	estimatedDegradation = 0.0;
659	numberUnsatisfied_ = 0;
660	// initialize sum of "infeasibilities"
661	sumInfeasibilities_ = 0.0;
662	int bestPriority = COIN_INT_MAX;
663	/*
664	Scan for branching objects that indicate infeasibility. Choose the best
665	maximumStrong candidates, using priority as the first criteria, then
666	integer infeasibility.
667
668	The algorithm is to fill the choice array with a set of good candidates (by
669	infeasibility) with priority bestPriority. Finding a candidate with
670	priority better (less) than bestPriority flushes the choice array. (This
671	serves as initialization when the first candidate is found.)
672
673	A new candidate is added to choices only if its infeasibility exceeds the
674	current max infeasibility (mostAway). When a candidate is added, it
675	replaces the candidate with the smallest infeasibility (tracked by
676	iSmallest).
677	*/
678	int iSmallest = 0;
679	double mostAway = 1.0e-100;
680	for (i = 0 ; i < maximumStrong ; i++)
681	choice[i].possibleBranch = NULL ;
682	numberStrong = 0;
683	bool canDoOneHot = false;
684	for (i = 0; i < numberObjects; i++) {
685	OsiObject * object = model->modifiableObject(i);
686	int preferredWay;
687	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
688	int priorityLevel = object->priority();
689	if (hotstartSolution) {
690	// we are doing hot start
691	const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
692	if (thisOne) {
693	int iColumn = thisOne->columnNumber();
694	bool canDoThisHot = true;
695	double targetValue = hotstartSolution[iColumn];
696	if (saveUpper[iColumn] > saveLower[iColumn]) {
697	double value = saveSolution[iColumn];
698	if (hotstartPriorities)
699	priorityLevel = hotstartPriorities[iColumn];
700	//double originalLower = thisOne->originalLower();
701	//double originalUpper = thisOne->originalUpper();
702	// switch off if not possible
703	if (targetValue >= saveLower[iColumn] && targetValue <= saveUpper[iColumn]) {
704	/* priority outranks rest always if negative
705	otherwise can be downgraded if at correct level.
706	Infeasibility may be increased to choose 1.0 values first.
707	choose one near wanted value
708	*/
709	if (fabs(value - targetValue) > integerTolerance) {
710	//if (infeasibility>0.01)
711	//infeasibility = fabs(1.0e6-fabs(value-targetValue));
712	//else
713	infeasibility = fabs(value - targetValue);
714	//if (targetValue==1.0)
715	//infeasibility += 1.0;
716	if (value > targetValue) {
717	preferredWay = -1;
718	} else {
719	preferredWay = 1;
720	}
721	priorityLevel = CoinAbs(priorityLevel);
722	} else if (priorityLevel < 0) {
723	priorityLevel = CoinAbs(priorityLevel);
724	if (targetValue == saveLower[iColumn]) {
725	infeasibility = integerTolerance + 1.0e-12;
726	preferredWay = -1;
727	} else if (targetValue == saveUpper[iColumn]) {
728	infeasibility = integerTolerance + 1.0e-12;
729	preferredWay = 1;
730	} else {
731	// can't
732	priorityLevel += 10000000;
733	canDoThisHot = false;
734	}
735	} else {
736	priorityLevel += 10000000;
737	canDoThisHot = false;
738	}
739	} else {
740	// switch off if not possible
741	canDoThisHot = false;
742	}
743	if (canDoThisHot)
744	canDoOneHot = true;
745	} else if (targetValue < saveLower[iColumn] \|\| targetValue > saveUpper[iColumn]) {
746	}
747	} else {
748	priorityLevel += 10000000;
749	}
750	}
751	if (infeasibility) {
752	// Increase estimated degradation to solution
753	estimatedDegradation += CoinMin(object->upEstimate(), object->downEstimate());
754	numberUnsatisfied_++;
755	sumInfeasibilities_ += infeasibility;
756	// Better priority? Flush choices.
757	if (priorityLevel < bestPriority) {
758	int j;
759	iSmallest = 0;
760	for (j = 0; j < maximumStrong; j++) {
761	choice[j].upMovement = 0.0;
762	delete choice[j].possibleBranch;
763	choice[j].possibleBranch = NULL;
764	}
765	bestPriority = priorityLevel;
766	mostAway = 1.0e-100;
767	numberStrong = 0;
768	} else if (priorityLevel > bestPriority) {
769	continue;
770	}
771	// Check for suitability based on infeasibility.
772	if (infeasibility > mostAway) {
773	//add to list
774	choice[iSmallest].upMovement = infeasibility;
775	delete choice[iSmallest].possibleBranch;
776	CbcObject * obj =
777	dynamic_cast <CbcObject *>(object) ;
778	assert (obj);
779	choice[iSmallest].possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
780	numberStrong = CoinMax(numberStrong, iSmallest + 1);
781	// Save which object it was
782	choice[iSmallest].objectNumber = i;
783	int j;
784	iSmallest = -1;
785	mostAway = 1.0e50;
786	for (j = 0; j < maximumStrong; j++) {
787	if (choice[j].upMovement < mostAway) {
788	mostAway = choice[j].upMovement;
789	iSmallest = j;
790	}
791	}
792	}
793	}
794	}
795	if (!canDoOneHot && hotstartSolution) {
796	// switch off as not possible
797	hotstartSolution = NULL;
798	model->setHotstartSolution(NULL, NULL);
799	usefulInfo.hotstartSolution_ = NULL;
800	}
801	if (numberUnsatisfied_) {
802	// some infeasibilities - go to next steps
803	#ifdef CLP_INVESTIGATE
804	if (hotstartSolution) {
805	int k = choice[0].objectNumber;
806	OsiObject * object = model->modifiableObject(k);
807	const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
808	assert (thisOne);
809	int iColumn = thisOne->columnNumber();
810	double targetValue = hotstartSolution[iColumn];
811	double value = saveSolution[iColumn];
812	printf("Branch on %d has target %g (value %g) and current bounds %g and %g\n",
813	iColumn, targetValue, value, saveLower[iColumn], saveUpper[iColumn]);
814	}
815	#endif
816	break;
817	} else if (!iPass) {
818	// looks like a solution - get paranoid
819	bool roundAgain = false;
820	// get basis
821	CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
822	if (!ws)
823	break;
824	for (i = 0; i < numberColumns; i++) {
825	double value = saveSolution[i];
826	if (value < lower[i]) {
827	saveSolution[i] = lower[i];
828	roundAgain = true;
829	ws->setStructStatus(i, CoinWarmStartBasis::atLowerBound);
830	} else if (value > upper[i]) {
831	saveSolution[i] = upper[i];
832	roundAgain = true;
833	ws->setStructStatus(i, CoinWarmStartBasis::atUpperBound);
834	}
835	}
836	if (roundAgain && saveNumberStrong) {
837	// restore basis
838	solver->setWarmStart(ws);
839	delete ws;
840	solver->resolve();
841	memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
842	model->reserveCurrentSolution(saveSolution);
843	if (!solver->isProvenOptimal()) {
844	// infeasible
845	anyAction = -2;
846	break;
847	}
848	} else {
849	delete ws;
850	break;
851	}
852	}
853	}
854	/* Some solvers can do the strong branching calculations faster if
855	they do them all at once. At present only Clp does for ordinary
856	integers but I think this coding would be easy to modify
857	*/
858	bool allNormal = true; // to say if we can do fast strong branching
859	// Say which one will be best
860	int bestChoice = 0;
861	double worstInfeasibility = 0.0;
862	for (i = 0; i < numberStrong; i++) {
863	choice[i].numIntInfeasUp = numberUnsatisfied_;
864	choice[i].numIntInfeasDown = numberUnsatisfied_;
865	choice[i].fix = 0; // say not fixed
866	if (!dynamic_cast <const CbcSimpleInteger *> (model->object(choice[i].objectNumber)))
867	allNormal = false; // Something odd so lets skip clever fast branching
868	if ( !model->object(choice[i].objectNumber)->boundBranch())
869	numberStrong = 0; // switch off
870	if ( choice[i].possibleBranch->numberBranches() > 2)
871	numberStrong = 0; // switch off
872	// Do best choice in case switched off
873	if (choice[i].upMovement > worstInfeasibility) {
874	worstInfeasibility = choice[i].upMovement;
875	bestChoice = i;
876	}
877	}
878	// If we have hit max time don't do strong branching
879	bool hitMaxTime = (model->getCurrentSeconds() >
880	model->getDblParam(CbcModel::CbcMaximumSeconds));
881	// also give up if we are looping round too much
882	if (hitMaxTime \|\| numberPassesLeft <= 0)
883	numberStrong = 0;
884	/*
885	Is strong branching enabled? If so, set up and do it. Otherwise, we'll
886	fall through to simple branching.
887
888	Setup for strong branching involves saving the current basis (for restoration
889	afterwards) and setting up for hot starts.
890	*/
891	if (numberStrong && saveNumberStrong) {
892
893	bool solveAll = false; // set true to say look at all even if some fixed (experiment)
894	solveAll = true;
895	// worth trying if too many times
896	// Save basis
897	CoinWarmStart * ws = solver->getWarmStart();
898	// save limit
899	int saveLimit;
900	solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit);
901	if (beforeSolution && saveLimit < 100)
902	solver->setIntParam(OsiMaxNumIterationHotStart, 100); // go to end
903	# ifdef COIN_HAS_CLP
904	/* If we are doing all strong branching in one go then we create new arrays
905	to store information. If clp NULL then doing old way.
906	Going down -
907	outputSolution[2*i] is final solution.
908	outputStuff[2*i] is status (0 - finished, 1 infeas, other unknown
909	outputStuff[2*i+numberStrong] is number iterations
910	On entry newUpper[i] is new upper bound, on exit obj change
911	Going up -
912	outputSolution[2*i+1] is final solution.
913	outputStuff[2*i+1] is status (0 - finished, 1 infeas, other unknown
914	outputStuff[2*i+1+numberStrong] is number iterations
915	On entry newLower[i] is new lower bound, on exit obj change
916	*/
917	ClpSimplex * clp = NULL;
918	double * newLower = NULL;
919	double * newUpper = NULL;
920	double ** outputSolution = NULL;
921	int * outputStuff = NULL;
922	// Go back to normal way if user wants it
923	if (osiclp && (osiclp->specialOptions()&16) != 0 && osiclp->specialOptions() > 0)
924	allNormal = false;
925	if (osiclp && !allNormal) {
926	// say do fast
927	int easy = 1;
928	osiclp->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
929	}
930	if (osiclp && allNormal) {
931	clp = osiclp->getModelPtr();
932	// Clp - do a different way
933	newLower = new double[numberStrong];
934	newUpper = new double[numberStrong];
935	outputSolution = new double * [2*numberStrong];
936	outputStuff = new int [4*numberStrong];
937	int * which = new int[numberStrong];
938	int startFinishOptions;
939	int specialOptions = osiclp->specialOptions();
940	int clpOptions = clp->specialOptions();
941	int returnCode = 0;
942	#define CRUNCH
943	#ifdef CRUNCH
944	// Crunch down problem
945	int numberRows = clp->numberRows();
946	// Use dual region
947	double * rhs = clp->dualRowSolution();
948	int * whichRow = new int[3*numberRows];
949	int * whichColumn = new int[2*numberColumns];
950	int nBound;
951	ClpSimplex * small = static_cast<ClpSimplexOther *> (clp)->crunch(rhs, whichRow, whichColumn, nBound, true);
952	if (!small) {
953	anyAction = -2;
954	//printf("XXXX Inf by inspection\n");
955	delete [] whichColumn;
956	whichColumn = NULL;
957	delete [] whichRow;
958	whichRow = NULL;
959	break;
960	} else {
961	clp = small;
962	}
963	#else
964	int saveLogLevel = clp->logLevel();
965	int saveMaxIts = clp->maximumIterations();
966	#endif
967	clp->setLogLevel(0);
968	if ((specialOptions&1) == 0) {
969	startFinishOptions = 0;
970	clp->setSpecialOptions(clpOptions \| (64 \| 1024));
971	} else {
972	startFinishOptions = 1 + 2 + 4;
973	//startFinishOptions=1+4; // for moment re-factorize
974	if ((specialOptions&4) == 0)
975	clp->setSpecialOptions(clpOptions \| (64 \| 128 \| 512 \| 1024 \| 4096));
976	else
977	clp->setSpecialOptions(clpOptions \| (64 \| 128 \| 512 \| 1024 \| 2048 \| 4096));
978	}
979	// User may want to clean up before strong branching
980	if ((clp->specialOptions()&32) != 0) {
981	clp->primal(1);
982	if (clp->numberIterations())
983	model->messageHandler()->message(CBC_ITERATE_STRONG, *model->messagesPointer())
984	<< clp->numberIterations()
985	<< CoinMessageEol;
986	}
987	clp->setMaximumIterations(saveLimit);
988	#ifdef CRUNCH
989	int * backColumn = whichColumn + numberColumns;
990	#endif
991	for (i = 0; i < numberStrong; i++) {
992	int iObject = choice[i].objectNumber;
993	const OsiObject * object = model->object(iObject);
994	const CbcSimpleInteger * simple = static_cast <const CbcSimpleInteger *> (object);
995	int iSequence = simple->columnNumber();
996	newLower[i] = ceil(saveSolution[iSequence]);
997	newUpper[i] = floor(saveSolution[iSequence]);
998	#ifdef CRUNCH
999	iSequence = backColumn[iSequence];
1000	assert (iSequence >= 0);
1001	#endif
1002	which[i] = iSequence;
1003	outputSolution[2*i] = new double [numberColumns];
1004	outputSolution[2*i+1] = new double [numberColumns];
1005	}
1006	//clp->writeMps("bad");
1007	returnCode = clp->strongBranching(numberStrong, which,
1008	newLower, newUpper, outputSolution,
1009	outputStuff, outputStuff + 2 * numberStrong, !solveAll, false,
1010	startFinishOptions);
1011	#ifndef CRUNCH
1012	clp->setSpecialOptions(clpOptions); // restore
1013	clp->setMaximumIterations(saveMaxIts);
1014	clp->setLogLevel(saveLogLevel);
1015	#endif
1016	if (returnCode == -2) {
1017	// bad factorization!!!
1018	// Doing normal way
1019	// Mark hot start
1020	solver->markHotStart();
1021	clp = NULL;
1022	} else {
1023	#ifdef CRUNCH
1024	// extract solution
1025	//bool checkSol=true;
1026	for (i = 0; i < numberStrong; i++) {
1027	int iObject = choice[i].objectNumber;
1028	const OsiObject * object = model->object(iObject);
1029	const CbcSimpleInteger * simple = static_cast <const CbcSimpleInteger *> (object);
1030	int iSequence = simple->columnNumber();
1031	which[i] = iSequence;
1032	double * sol = outputSolution[2*i];
1033	double * sol2 = outputSolution[2*i+1];
1034	//bool x=true;
1035	//bool x2=true;
1036	for (int iColumn = numberColumns - 1; iColumn >= 0; iColumn--) {
1037	int jColumn = backColumn[iColumn];
1038	if (jColumn >= 0) {
1039	sol[iColumn] = sol[jColumn];
1040	sol2[iColumn] = sol2[jColumn];
1041	} else {
1042	sol[iColumn] = saveSolution[iColumn];
1043	sol2[iColumn] = saveSolution[iColumn];
1044	}
1045	}
1046	}
1047	#endif
1048	}
1049	#ifdef CRUNCH
1050	delete [] whichColumn;
1051	delete [] whichRow;
1052	delete small;
1053	#endif
1054	delete [] which;
1055	} else {
1056	// Doing normal way
1057	// Mark hot start
1058	solver->markHotStart();
1059	}
1060	# else /* COIN_HAS_CLP */
1061
1062	OsiSolverInterface *clp = NULL ;
1063	double **outputSolution = NULL ;
1064	int *outputStuff = NULL ;
1065	double * newLower = NULL ;
1066	double * newUpper = NULL ;
1067
1068	solver->markHotStart();
1069
1070	# endif /* COIN_HAS_CLP */
1071	/*
1072	Open a loop to do the strong branching LPs. For each candidate variable,
1073	solve an LP with the variable forced down, then up. If a direction turns
1074	out to be infeasible or monotonic (i.e., over the dual objective cutoff),
1075	force the objective change to be big (1.0e100). If we determine the problem
1076	is infeasible, or find a monotone variable, escape the loop.
1077
1078	TODO: The `restore bounds' part might be better encapsulated as an
1079	unbranch() method. Branching objects more exotic than simple integers
1080	or cliques might not restrict themselves to variable bounds.
1081
1082	TODO: Virtuous solvers invalidate the current solution (or give bogus
1083	results :-) when the bounds are changed out from under them. So we
1084	need to do all the work associated with finding a new solution before
1085	restoring the bounds.
1086	*/
1087	for (i = 0 ; i < numberStrong ; i++) {
1088	double objectiveChange ;
1089	double newObjectiveValue = 1.0e100;
1090	// status is 0 finished, 1 infeasible and other
1091	int iStatus;
1092	/*
1093	Try the down direction first. (Specify the initial branching alternative as
1094	down with a call to way(-1). Each subsequent call to branch() performs the
1095	specified branch and advances the branch object state to the next branch
1096	alternative.)
1097	*/
1098	if (!clp) {
1099	choice[i].possibleBranch->way(-1) ;
1100	choice[i].possibleBranch->branch() ;
1101	bool feasible = true;
1102	if (checkFeasibility) {
1103	// check branching did not make infeasible
1104	int iColumn;
1105	int numberColumns = solver->getNumCols();
1106	const double * columnLower = solver->getColLower();
1107	const double * columnUpper = solver->getColUpper();
1108	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1109	if (columnLower[iColumn] > columnUpper[iColumn] + 1.0e-5)
1110	feasible = false;
1111	}
1112	}
1113	if (feasible) {
1114	solver->solveFromHotStart() ;
1115	numberStrongDone++;
1116	numberStrongIterations += solver->getIterationCount();
1117	/*
1118	We now have an estimate of objective degradation that we can use for strong
1119	branching. If we're over the cutoff, the variable is monotone up.
1120	If we actually made it to optimality, check for a solution, and if we have
1121	a good one, call setBestSolution to process it. Note that this may reduce the
1122	cutoff, so we check again to see if we can declare this variable monotone.
1123	*/
1124	if (solver->isProvenOptimal())
1125	iStatus = 0; // optimal
1126	else if (solver->isIterationLimitReached()
1127	&& !solver->isDualObjectiveLimitReached())
1128	iStatus = 2; // unknown
1129	else
1130	iStatus = 1; // infeasible
1131	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1132	choice[i].numItersDown = solver->getIterationCount();
1133	} else {
1134	iStatus = 1; // infeasible
1135	newObjectiveValue = 1.0e100;
1136	choice[i].numItersDown = 0;
1137	}
1138	} else {
1139	iStatus = outputStuff[2*i];
1140	choice[i].numItersDown = outputStuff[2numberStrong+2i];
1141	numberStrongDone++;
1142	numberStrongIterations += choice[i].numItersDown;
1143	newObjectiveValue = objectiveValue + newUpper[i];
1144	solver->setColSolution(outputSolution[2*i]);
1145	}
1146	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0);
1147	if (!iStatus) {
1148	choice[i].finishedDown = true ;
1149	if (newObjectiveValue >= model->getCutoff()) {
1150	objectiveChange = 1.0e100; // say infeasible
1151	numberStrongInfeasible++;
1152	} else {
1153	// See if integer solution
1154	if (model->feasibleSolution(choice[i].numIntInfeasDown,
1155	choice[i].numObjInfeasDown)
1156	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
1157	model->setBestSolution(CBC_STRONGSOL,
1158	newObjectiveValue,
1159	solver->getColSolution()) ;
1160	// only needed for odd solvers
1161	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1162	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0) ;
1163	model->setLastHeuristic(NULL);
1164	model->incrementUsed(solver->getColSolution());
1165	if (newObjectiveValue >= model->getCutoff()) { // new cutoff
1166	objectiveChange = 1.0e100 ;
1167	numberStrongInfeasible++;
1168	}
1169	}
1170	}
1171	} else if (iStatus == 1) {
1172	objectiveChange = 1.0e100 ;
1173	numberStrongInfeasible++;
1174	} else {
1175	// Can't say much as we did not finish
1176	choice[i].finishedDown = false ;
1177	numberUnfinished++;
1178	}
1179	choice[i].downMovement = objectiveChange ;
1180
1181	// restore bounds
1182	if (!clp) {
1183	for (int j = 0; j < numberColumns; j++) {
1184	if (saveLower[j] != lower[j])
1185	solver->setColLower(j, saveLower[j]);
1186	if (saveUpper[j] != upper[j])
1187	solver->setColUpper(j, saveUpper[j]);
1188	}
1189	}
1190	//printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
1191	// choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersDown,
1192	// choice[i].downMovement,choice[i].finishedDown,choice[i].numIntInfeasDown,
1193	// choice[i].numObjInfeasDown);
1194
1195	// repeat the whole exercise, forcing the variable up
1196	if (!clp) {
1197	bool feasible = true;
1198	// If odd branching then maybe just one possibility
1199	if (choice[i].possibleBranch->numberBranchesLeft() > 0) {
1200	choice[i].possibleBranch->branch();
1201	if (checkFeasibility) {
1202	// check branching did not make infeasible
1203	int iColumn;
1204	int numberColumns = solver->getNumCols();
1205	const double * columnLower = solver->getColLower();
1206	const double * columnUpper = solver->getColUpper();
1207	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1208	if (columnLower[iColumn] > columnUpper[iColumn] + 1.0e-5)
1209	feasible = false;
1210	}
1211	}
1212	} else {
1213	// second branch infeasible
1214	feasible = false;
1215	}
1216	if (feasible) {
1217	solver->solveFromHotStart() ;
1218	numberStrongDone++;
1219	numberStrongIterations += solver->getIterationCount();
1220	/*
1221	We now have an estimate of objective degradation that we can use for strong
1222	branching. If we're over the cutoff, the variable is monotone up.
1223	If we actually made it to optimality, check for a solution, and if we have
1224	a good one, call setBestSolution to process it. Note that this may reduce the
1225	cutoff, so we check again to see if we can declare this variable monotone.
1226	*/
1227	if (solver->isProvenOptimal())
1228	iStatus = 0; // optimal
1229	else if (solver->isIterationLimitReached()
1230	&& !solver->isDualObjectiveLimitReached())
1231	iStatus = 2; // unknown
1232	else
1233	iStatus = 1; // infeasible
1234	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1235	choice[i].numItersUp = solver->getIterationCount();
1236	} else {
1237	iStatus = 1; // infeasible
1238	newObjectiveValue = 1.0e100;
1239	choice[i].numItersDown = 0;
1240	}
1241	} else {
1242	iStatus = outputStuff[2*i+1];
1243	choice[i].numItersUp = outputStuff[2numberStrong+2i+1];
1244	numberStrongDone++;
1245	numberStrongIterations += choice[i].numItersUp;
1246	newObjectiveValue = objectiveValue + newLower[i];
1247	solver->setColSolution(outputSolution[2*i+1]);
1248	}
1249	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0);
1250	if (!iStatus) {
1251	choice[i].finishedUp = true ;
1252	if (newObjectiveValue >= model->getCutoff()) {
1253	objectiveChange = 1.0e100; // say infeasible
1254	numberStrongInfeasible++;
1255	} else {
1256	// See if integer solution
1257	if (model->feasibleSolution(choice[i].numIntInfeasUp,
1258	choice[i].numObjInfeasUp)
1259	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
1260	model->setBestSolution(CBC_STRONGSOL,
1261	newObjectiveValue,
1262	solver->getColSolution()) ;
1263	// only needed for odd solvers
1264	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
1265	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0) ;
1266	model->setLastHeuristic(NULL);
1267	model->incrementUsed(solver->getColSolution());
1268	if (newObjectiveValue >= model->getCutoff()) { // new cutoff
1269	objectiveChange = 1.0e100 ;
1270	numberStrongInfeasible++;
1271	}
1272	}
1273	}
1274	} else if (iStatus == 1) {
1275	objectiveChange = 1.0e100 ;
1276	numberStrongInfeasible++;
1277	} else {
1278	// Can't say much as we did not finish
1279	choice[i].finishedUp = false ;
1280	numberUnfinished++;
1281	}
1282	choice[i].upMovement = objectiveChange ;
1283
1284	// restore bounds
1285	if (!clp) {
1286	for (int j = 0; j < numberColumns; j++) {
1287	if (saveLower[j] != lower[j])
1288	solver->setColLower(j, saveLower[j]);
1289	if (saveUpper[j] != upper[j])
1290	solver->setColUpper(j, saveUpper[j]);
1291	}
1292	}
1293
1294	//printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
1295	// choice[i].objectNumber,iStatus,newObjectiveValue,choice[i].numItersUp,
1296	// choice[i].upMovement,choice[i].finishedUp,choice[i].numIntInfeasUp,
1297	// choice[i].numObjInfeasUp);
1298
1299	/*
1300	End of evaluation for this candidate variable. Possibilities are:
1301	* Both sides below cutoff; this variable is a candidate for branching.
1302	* Both sides infeasible or above the objective cutoff: no further action
1303	here. Break from the evaluation loop and assume the node will be purged
1304	by the caller.
1305	* One side below cutoff: Install the branch (i.e., fix the variable). Break
1306	from the evaluation loop and assume the node will be reoptimised by the
1307	caller.
1308	*/
1309	// reset
1310	choice[i].possibleBranch->resetNumberBranchesLeft();
1311	if (choice[i].upMovement < 1.0e100) {
1312	if (choice[i].downMovement < 1.0e100) {
1313	// feasible - no action
1314	} else {
1315	// up feasible, down infeasible
1316	anyAction = -1;
1317	//printf("Down infeasible for choice %d sequence %d\n",i,
1318	// model->object(choice[i].objectNumber)->columnNumber());
1319	if (!solveAll) {
1320	choice[i].possibleBranch->way(1);
1321	choice[i].possibleBranch->branch();
1322	break;
1323	} else {
1324	choice[i].fix = 1;
1325	}
1326	}
1327	} else {
1328	if (choice[i].downMovement < 1.0e100) {
1329	// down feasible, up infeasible
1330	anyAction = -1;
1331	//printf("Up infeasible for choice %d sequence %d\n",i,
1332	// model->object(choice[i].objectNumber)->columnNumber());
1333	if (!solveAll) {
1334	choice[i].possibleBranch->way(-1);
1335	choice[i].possibleBranch->branch();
1336	break;
1337	} else {
1338	choice[i].fix = -1;
1339	}
1340	} else {
1341	// neither side feasible
1342	anyAction = -2;
1343	//printf("Both infeasible for choice %d sequence %d\n",i,
1344	// model->object(choice[i].objectNumber)->columnNumber());
1345	break;
1346	}
1347	}
1348	bool hitMaxTime = (model->getCurrentSeconds() >
1349	model->getDblParam(CbcModel::CbcMaximumSeconds));
1350	if (hitMaxTime) {
1351	numberStrong = i + 1;
1352	break;
1353	}
1354	}
1355	if (!clp) {
1356	// Delete the snapshot
1357	solver->unmarkHotStart();
1358	} else {
1359	delete [] newLower;
1360	delete [] newUpper;
1361	delete [] outputStuff;
1362	int i;
1363	for (i = 0; i < 2*numberStrong; i++)
1364	delete [] outputSolution[i];
1365	delete [] outputSolution;
1366	}
1367	solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit);
1368	// restore basis
1369	solver->setWarmStart(ws);
1370	// Unless infeasible we will carry on
1371	// But we could fix anyway
1372	if (anyAction == -1 && solveAll) {
1373	// apply and take off
1374	for (i = 0 ; i < numberStrong ; i++) {
1375	if (choice[i].fix) {
1376	choice[i].possibleBranch->way(choice[i].fix) ;
1377	choice[i].possibleBranch->branch() ;
1378	}
1379	}
1380	bool feasible = true;
1381	if (checkFeasibility) {
1382	// check branching did not make infeasible
1383	int iColumn;
1384	int numberColumns = solver->getNumCols();
1385	const double * columnLower = solver->getColLower();
1386	const double * columnUpper = solver->getColUpper();
1387	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
1388	if (columnLower[iColumn] > columnUpper[iColumn] + 1.0e-5)
1389	feasible = false;
1390	}
1391	}
1392	if (feasible) {
1393	// can do quick optimality check
1394	int easy = 2;
1395	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
1396	solver->resolve() ;
1397	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
1398	feasible = solver->isProvenOptimal();
1399	}
1400	if (feasible) {
1401	memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
1402	model->reserveCurrentSolution(saveSolution);
1403	memcpy(saveLower, solver->getColLower(), numberColumns*sizeof(double));
1404	memcpy(saveUpper, solver->getColUpper(), numberColumns*sizeof(double));
1405	// Clean up all candidates whih are fixed
1406	int numberLeft = 0;
1407	for (i = 0 ; i < numberStrong ; i++) {
1408	CbcStrongInfo thisChoice = choice[i];
1409	choice[i].possibleBranch = NULL;
1410	const OsiObject * object = model->object(thisChoice.objectNumber);
1411	int preferredWay;
1412	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
1413	if (!infeasibility) {
1414	// take out
1415	delete thisChoice.possibleBranch;
1416	} else {
1417	choice[numberLeft++] = thisChoice;
1418	}
1419	}
1420	numberStrong = numberLeft;
1421	for (; i < maximumStrong; i++) {
1422	delete choice[i].possibleBranch;
1423	choice[i].possibleBranch = NULL;
1424	}
1425	// If all fixed then round again
1426	if (!numberLeft) {
1427	finished = false;
1428	numberStrong = 0;
1429	saveNumberStrong = 0;
1430	maximumStrong = 1;
1431	} else {
1432	anyAction = 0;
1433	}
1434	// If these two uncommented then different action
1435	anyAction = -1;
1436	finished = true;
1437	//printf("some fixed but continuing %d left\n",numberLeft);
1438	} else {
1439	anyAction = -2; // say infeasible
1440	}
1441	}
1442	delete ws;
1443	int numberNodes = model->getNodeCount();
1444	// update number of strong iterations etc
1445	model->incrementStrongInfo(numberStrongDone, numberStrongIterations,
1446	anyAction == -2 ? 0 : numberStrongInfeasible, anyAction == -2);
1447
1448	/*
1449	anyAction >= 0 indicates that strong branching didn't produce any monotone
1450	variables. Sift through the candidates for the best one.
1451
1452	QUERY: Setting numberNodes looks to be a distributed noop. numberNodes is
1453	local to this code block. Perhaps should be numberNodes_ from model?
1454	Unclear what this calculation is doing.
1455	*/
1456	if (anyAction >= 0) {
1457
1458	// get average cost per iteration and assume stopped ones
1459	// would stop after 50% more iterations at average cost??? !!! ???
1460	double averageCostPerIteration = 0.0;
1461	double totalNumberIterations = 1.0;
1462	int smallestNumberInfeasibilities = COIN_INT_MAX;
1463	for (i = 0; i < numberStrong; i++) {
1464	totalNumberIterations += choice[i].numItersDown +
1465	choice[i].numItersUp ;
1466	averageCostPerIteration += choice[i].downMovement +
1467	choice[i].upMovement;
1468	smallestNumberInfeasibilities =
1469	CoinMin(CoinMin(choice[i].numIntInfeasDown ,
1470	choice[i].numIntInfeasUp ),
1471	smallestNumberInfeasibilities);
1472	}
1473	//if (smallestNumberInfeasibilities>=numberIntegerInfeasibilities)
1474	//numberNodes=1000000; // switch off search for better solution
1475	numberNodes = 1000000; // switch off anyway
1476	averageCostPerIteration /= totalNumberIterations;
1477	// all feasible - choose best bet
1478
1479	// New method does all at once so it can be more sophisticated
1480	// in deciding how to balance actions.
1481	// But it does need arrays
1482	double * changeUp = new double [numberStrong];
1483	int * numberInfeasibilitiesUp = new int [numberStrong];
1484	double * changeDown = new double [numberStrong];
1485	int * numberInfeasibilitiesDown = new int [numberStrong];
1486	CbcBranchingObject ** objects = new CbcBranchingObject * [ numberStrong];
1487	for (i = 0 ; i < numberStrong ; i++) {
1488	int iColumn = choice[i].possibleBranch->variable() ;
1489	model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
1490	<< i << iColumn
1491	<< choice[i].downMovement << choice[i].numIntInfeasDown
1492	<< choice[i].upMovement << choice[i].numIntInfeasUp
1493	<< choice[i].possibleBranch->value()
1494	<< CoinMessageEol;
1495	changeUp[i] = choice[i].upMovement;
1496	numberInfeasibilitiesUp[i] = choice[i].numIntInfeasUp;
1497	changeDown[i] = choice[i].downMovement;
1498	numberInfeasibilitiesDown[i] = choice[i].numIntInfeasDown;
1499	objects[i] = choice[i].possibleBranch;
1500	}
1501	int whichObject = decision->bestBranch(objects, numberStrong, numberUnsatisfied_,
1502	changeUp, numberInfeasibilitiesUp,
1503	changeDown, numberInfeasibilitiesDown,
1504	objectiveValue_);
1505	// move branching object and make sure it will not be deleted
1506	if (whichObject >= 0) {
1507	branch_ = objects[whichObject];
1508	if (model->messageHandler()->logLevel() > 3)
1509	printf("Choosing column %d\n", choice[whichObject].possibleBranch->variable()) ;
1510	choice[whichObject].possibleBranch = NULL;
1511	}
1512	delete [] changeUp;
1513	delete [] numberInfeasibilitiesUp;
1514	delete [] changeDown;
1515	delete [] numberInfeasibilitiesDown;
1516	delete [] objects;
1517	}
1518	# ifdef COIN_HAS_CLP
1519	if (osiclp && !allNormal) {
1520	// back to normal
1521	osiclp->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
1522	}
1523	# endif
1524	}
1525	/*
1526	Simple branching. Probably just one, but we may have got here
1527	because of an odd branch e.g. a cut
1528	*/
1529	else {
1530	// not strong
1531	// C) create branching object
1532	branch_ = choice[bestChoice].possibleBranch;
1533	choice[bestChoice].possibleBranch = NULL;
1534	}
1535	}
1536	// Set guessed solution value
1537	guessedObjectiveValue_ = objectiveValue_ + estimatedDegradation;
1538	/*
1539	Cleanup, then we're outta here.
1540	*/
1541	if (!model->branchingMethod() \|\| dynamicBranchingObject)
1542	delete decision;
1543
1544	for (i = 0; i < maximumStrong; i++)
1545	delete choice[i].possibleBranch;
1546	delete [] choice;
1547	delete [] saveLower;
1548	delete [] saveUpper;
1549
1550	// restore solution
1551	solver->setColSolution(saveSolution);
1552	delete [] saveSolution;
1553	# ifdef COIN_HAS_CLP
1554	if (osiclp)
1555	osiclp->setSpecialOptions(saveClpOptions);
1556	# endif
1557	return anyAction;
1558	}
1559
1560	/*
1561	Version for dynamic pseudo costs.
1562
1563	**** For now just return if anything odd
1564	later allow even if odd
1565
1566	The routine scans through the object list of the model looking for objects
1567	that indicate infeasibility. It tests each object using strong branching
1568	and selects the one with the least objective degradation. A corresponding
1569	branching object is left attached to lastNode.
1570	This version gives preference in evaluation to variables which
1571	have not been evaluated many times. It also uses numberStrong
1572	to say give up if last few tries have not changed incumbent.
1573	See Achterberg, Koch and Martin.
1574
1575	If strong branching is disabled, a candidate object is chosen essentially
1576	at random (whatever object ends up in pos'n 0 of the candidate array).
1577
1578	If a branching candidate is found to be monotone, bounds are set to fix the
1579	variable and the routine immediately returns (the caller is expected to
1580	reoptimize).
1581
1582	If a branching candidate is found to result in infeasibility in both
1583	directions, the routine immediately returns an indication of infeasibility.
1584
1585	Returns: 0 both branch directions are feasible
1586	-1 branching variable is monotone
1587	-2 infeasible
1588	-3 Use another method
1589
1590	For now just fix on objective from strong branching.
1591	*/
1592
1593	int CbcNode::chooseDynamicBranch (CbcModel model, CbcNode lastNode,
1594	OsiSolverBranch * & /branches/,
1595	int numberPassesLeft)
1596
1597	{
1598	if (lastNode)
1599	depth_ = lastNode->depth_ + 1;
1600	else
1601	depth_ = 0;
1602	// Go to other choose if hot start
1603	if (model->hotstartSolution() &&
1604	(((model->moreSpecialOptions()&1024) == 0) \|\| false))
1605	return -3;
1606	delete branch_;
1607	branch_ = NULL;
1608	OsiSolverInterface * solver = model->solver();
1609	// get information on solver type
1610	const OsiAuxInfo * auxInfo = solver->getAuxiliaryInfo();
1611	const OsiBabSolver * auxiliaryInfo = dynamic_cast<const OsiBabSolver *> (auxInfo);
1612	if (!auxiliaryInfo) {
1613	// use one from CbcModel
1614	auxiliaryInfo = model->solverCharacteristics();
1615	}
1616	int numberObjects = model->numberObjects();
1617	// If very odd set of objects then use older chooseBranch
1618	bool useOldWay = false;
1619	// point to useful information
1620	OsiBranchingInformation usefulInfo = model->usefulInformation();
1621	if (numberObjects > model->numberIntegers()) {
1622	for (int i = model->numberIntegers(); i < numberObjects; i++) {
1623	OsiObject * object = model->modifiableObject(i);
1624	CbcObject * obj = dynamic_cast <CbcObject *>(object) ;
1625	if (!obj \|\| !obj->optionalObject()) {
1626	int preferredWay;
1627	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
1628	if (infeasibility) {
1629	useOldWay = true;
1630	break;
1631	}
1632	}
1633	}
1634	}
1635	if ((model->specialOptions()&128) != 0)
1636	useOldWay = false; // allow
1637	// For now return if not simple
1638	if (useOldWay)
1639	return -3;
1640	// Modify useful info
1641	usefulInfo.depth_ = depth_;
1642	if ((model->specialOptions()&128) != 0) {
1643	// SOS - shadow prices
1644	int numberRows = solver->getNumRows();
1645	const double * pi = usefulInfo.pi_;
1646	double sumPi = 0.0;
1647	for (int i = 0; i < numberRows; i++)
1648	sumPi += fabs(pi[i]);
1649	sumPi /= static_cast<double> (numberRows);
1650	// and scale back
1651	sumPi *= 0.01;
1652	usefulInfo.defaultDual_ = sumPi; // switch on
1653	int numberColumns = solver->getNumCols();
1654	int size = CoinMax(numberColumns, 2 * numberRows);
1655	usefulInfo.usefulRegion_ = new double [size];
1656	CoinZeroN(usefulInfo.usefulRegion_, size);
1657	usefulInfo.indexRegion_ = new int [size];
1658	// pi may change
1659	usefulInfo.pi_ = CoinCopyOfArray(usefulInfo.pi_, numberRows);
1660	}
1661	assert (auxiliaryInfo);
1662	double cutoff = model->getCutoff();
1663	const double * lower = solver->getColLower();
1664	const double * upper = solver->getColUpper();
1665	// See if user thinks infeasible
1666	int anyAction = model->problemFeasibility()->feasible(model, 0);
1667	if (anyAction) {
1668	// will return -2 if infeasible , 0 if treat as integer
1669	return anyAction - 1;
1670	}
1671	int i;
1672	int saveStateOfSearch = model->stateOfSearch() % 10;
1673	int numberStrong = model->numberStrong();
1674	/* Ranging is switched off.
1675	The idea is that you can find out the effect of one iteration
1676	on each unsatisfied variable cheaply. Then use this
1677	if you have not got much else to go on.
1678	*/
1679	//#define RANGING
1680	#ifdef RANGING
1681	// must have clp
1682	#ifndef COIN_HAS_CLP
1683	# warning("Ranging switched off as not Clp");
1684	#undef RANGING
1685	#endif
1686	// Pass number
1687	int kPass = 0;
1688	int numberRows = solver->getNumRows();
1689	#endif
1690	int numberColumns = model->getNumCols();
1691	double * saveUpper = new double[numberColumns];
1692	double * saveLower = new double[numberColumns];
1693	for (i = 0; i < numberColumns; i++) {
1694	saveLower[i] = lower[i];
1695	saveUpper[i] = upper[i];
1696	}
1697
1698	// Save solution in case heuristics need good solution later
1699
1700	double * saveSolution = new double[numberColumns];
1701	memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
1702	model->reserveCurrentSolution(saveSolution);
1703	const double * hotstartSolution = model->hotstartSolution();
1704	const int * hotstartPriorities = model->hotstartPriorities();
1705	double integerTolerance =
1706	model->getDblParam(CbcModel::CbcIntegerTolerance);
1707	if (hotstartSolution) {
1708	if ((model->moreSpecialOptions()&1024) != 0) {
1709	int nBad = 0;
1710	int nUnsat = 0;
1711	int nDiff = 0;
1712	for (int i = 0; i < numberObjects; i++) {
1713	OsiObject * object = model->modifiableObject(i);
1714	const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
1715	if (thisOne) {
1716	int iColumn = thisOne->columnNumber();
1717	double targetValue = hotstartSolution[iColumn];
1718	double value = saveSolution[iColumn];
1719	if (fabs(value - floor(value + 0.5)) > 1.0e-6) {
1720	nUnsat++;
1721	#ifdef CLP_INVESTIGATE
1722	printf("H %d is %g target %g\n", iColumn, value, targetValue);
1723	#endif
1724	} else if (fabs(targetValue - value) > 1.0e-6) {
1725	nDiff++;
1726	}
1727	if (targetValue < saveLower[iColumn] \|\|
1728	targetValue > saveUpper[iColumn]) {
1729	#ifdef CLP_INVESTIGATE
1730	printf("%d has target %g and current bounds %g and %g\n",
1731	iColumn, targetValue, saveLower[iColumn], saveUpper[iColumn]);
1732	#endif
1733	nBad++;
1734	}
1735	}
1736	}
1737	#ifdef CLP_INVESTIGATE
1738	printf("Hot %d unsatisfied, %d outside limits, %d different\n",
1739	nUnsat, nBad, nDiff);
1740	#endif
1741	if (nBad) {
1742	// switch off as not possible
1743	hotstartSolution = NULL;
1744	model->setHotstartSolution(NULL, NULL);
1745	usefulInfo.hotstartSolution_ = NULL;
1746	}
1747	}
1748	}
1749	/*
1750	Get a branching decision object. Use the default dynamic decision criteria unless
1751	the user has loaded a decision method into the model.
1752	*/
1753	CbcBranchDecision *decision = model->branchingMethod();
1754	if (!decision)
1755	decision = new CbcBranchDynamicDecision();
1756	int xMark = 0;
1757	// Get arrays to sort
1758	double * sort = new double[numberObjects];
1759	int * whichObject = new int[numberObjects];
1760	#ifdef RANGING
1761	int xPen = 0;
1762	int * objectMark = new int[2*numberObjects+1];
1763	#endif
1764	// Arrays with movements
1765	double * upEstimate = new double[numberObjects];
1766	double * downEstimate = new double[numberObjects];
1767	double estimatedDegradation = 0.0;
1768	int numberNodes = model->getNodeCount();
1769	int saveLogLevel = model->logLevel();
1770	#ifdef JJF_ZERO
1771	if ((numberNodes % 500) == 0) {
1772	model->setLogLevel(6);
1773	// Get average up and down costs
1774	double averageUp = 0.0;
1775	double averageDown = 0.0;
1776	int numberUp = 0;
1777	int numberDown = 0;
1778	int i;
1779	for ( i = 0; i < numberObjects; i++) {
1780	OsiObject * object = model->modifiableObject(i);
1781	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
1782	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
1783	assert(dynamicObject);
1784	int numberUp2 = 0;
1785	int numberDown2 = 0;
1786	double up = 0.0;
1787	double down = 0.0;
1788	if (dynamicObject->numberTimesUp()) {
1789	numberUp++;
1790	averageUp += dynamicObject->upDynamicPseudoCost();
1791	numberUp2 += dynamicObject->numberTimesUp();
1792	up = dynamicObject->upDynamicPseudoCost();
1793	}
1794	if (dynamicObject->numberTimesDown()) {
1795	numberDown++;
1796	averageDown += dynamicObject->downDynamicPseudoCost();
1797	numberDown2 += dynamicObject->numberTimesDown();
1798	down = dynamicObject->downDynamicPseudoCost();
1799	}
1800	if (numberUp2 \|\| numberDown2)
1801	printf("col %d - up %d times cost %g, - down %d times cost %g\n",
1802	dynamicObject->columnNumber(), numberUp2, up, numberDown2, down);
1803	}
1804	if (numberUp)
1805	averageUp /= static_cast<double> (numberUp);
1806	else
1807	averageUp = 1.0;
1808	if (numberDown)
1809	averageDown /= static_cast<double> (numberDown);
1810	else
1811	averageDown = 1.0;
1812	printf("total - up %d vars average %g, - down %d vars average %g\n",
1813	numberUp, averageUp, numberDown, averageDown);
1814	}
1815	#endif
1816	int numberBeforeTrust = model->numberBeforeTrust();
1817	// May go round twice if strong branching fixes all local candidates
1818	bool finished = false;
1819	int numberToFix = 0;
1820	# ifdef COIN_HAS_CLP
1821	OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
1822	int saveClpOptions = 0;
1823	if (osiclp) {
1824	// for faster hot start
1825	saveClpOptions = osiclp->specialOptions();
1826	osiclp->setSpecialOptions(saveClpOptions \| 8192);
1827	}
1828	# else
1829	OsiSolverInterface *osiclp = NULL ;
1830	# endif
1831	//const CglTreeProbingInfo * probingInfo = NULL; //model->probingInfo();
1832	// Old code left in with DEPRECATED_STRATEGY
1833	assert (model->searchStrategy() == -1 \|\|
1834	model->searchStrategy() == 1 \|\|
1835	model->searchStrategy() == 2);
1836	#ifdef DEPRECATED_STRATEGY
1837	int saveSearchStrategy2 = model->searchStrategy();
1838	#endif
1839	// Get average up and down costs
1840	{
1841	double averageUp = 0.0;
1842	double averageDown = 0.0;
1843	int numberUp = 0;
1844	int numberDown = 0;
1845	int i;
1846	for ( i = 0; i < numberObjects; i++) {
1847	OsiObject * object = model->modifiableObject(i);
1848	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
1849	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
1850	if (dynamicObject) {
1851	if (dynamicObject->numberTimesUp()) {
1852	numberUp++;
1853	averageUp += dynamicObject->upDynamicPseudoCost();
1854	}
1855	if (dynamicObject->numberTimesDown()) {
1856	numberDown++;
1857	averageDown += dynamicObject->downDynamicPseudoCost();
1858	}
1859	}
1860	}
1861	if (numberUp)
1862	averageUp /= static_cast<double> (numberUp);
1863	else
1864	averageUp = 1.0;
1865	if (numberDown)
1866	averageDown /= static_cast<double> (numberDown);
1867	else
1868	averageDown = 1.0;
1869	for ( i = 0; i < numberObjects; i++) {
1870	OsiObject * object = model->modifiableObject(i);
1871	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
1872	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
1873	if (dynamicObject) {
1874	if (!dynamicObject->numberTimesUp())
1875	dynamicObject->setUpDynamicPseudoCost(averageUp);
1876	if (!dynamicObject->numberTimesDown())
1877	dynamicObject->setDownDynamicPseudoCost(averageDown);
1878	}
1879	}
1880	}
1881	/*
1882	1 strong
1883	2 no strong
1884	3 strong just before solution
1885	4 no strong just before solution
1886	5 strong first time or before solution
1887	6 strong first time
1888	*/
1889	int useShadow = model->moreSpecialOptions() & 7;
1890	if (useShadow > 2) {
1891	if (model->getSolutionCount()) {
1892	if (numberNodes \|\| useShadow < 5) {
1893	useShadow = 0;
1894	// zap pseudo shadow prices
1895	model->pseudoShadow(-1);
1896	// and switch off
1897	model->setMoreSpecialOptions(model->moreSpecialOptions()&(~1023));
1898	} else {
1899	useShadow = 1;
1900	}
1901	} else if (useShadow < 5) {
1902	useShadow -= 2;
1903	} else {
1904	useShadow = 1;
1905	}
1906	}
1907	if (useShadow) {
1908	// pseudo shadow prices
1909	model->pseudoShadow((model->moreSpecialOptions() >> 3)&63);
1910	}
1911	#ifdef DEPRECATED_STRATEGY
1912	{ // in for tabbing
1913	} else if (saveSearchStrategy2 < 1999) {
1914	// pseudo shadow prices
1915	model->pseudoShadow(NULL, NULL);
1916	} else if (saveSearchStrategy2 < 2999) {
1917	// leave old ones
1918	} else if (saveSearchStrategy2 < 3999) {
1919	// pseudo shadow prices at root
1920	if (!numberNodes)
1921	model->pseudoShadow(NULL, NULL);
1922	} else {
1923	abort();
1924	}
1925	if (saveSearchStrategy2 >= 0)
1926	saveSearchStrategy2 = saveSearchStrategy2 % 1000;
1927	if (saveSearchStrategy2 == 999)
1928	saveSearchStrategy2 = -1;
1929	int saveSearchStrategy = saveSearchStrategy2 < 99 ? saveSearchStrategy2 : saveSearchStrategy2 - 100;
1930	#endif //DEPRECATED_STRATEGY
1931	int numberNotTrusted = 0;
1932	int numberStrongDone = 0;
1933	int numberUnfinished = 0;
1934	int numberStrongInfeasible = 0;
1935	int numberStrongIterations = 0;
1936	// so we can save lots of stuff
1937	CbcStrongInfo choice;
1938	CbcDynamicPseudoCostBranchingObject * choiceObject = NULL;
1939	if (model->allDynamic()) {
1940	CbcSimpleIntegerDynamicPseudoCost * object = NULL;
1941	choiceObject = new CbcDynamicPseudoCostBranchingObject(model, 0, -1, 0.5, object);
1942	}
1943	choice.possibleBranch = choiceObject;
1944	numberPassesLeft = CoinMax(numberPassesLeft, 2);
1945	while (!finished) {
1946	numberPassesLeft--;
1947	finished = true;
1948	decision->initialize(model);
1949	// Some objects may compute an estimate of best solution from here
1950	estimatedDegradation = 0.0;
1951	numberToFix = 0;
1952	int numberToDo = 0;
1953	int iBestNot = -1;
1954	int iBestGot = -1;
1955	double best = 0.0;
1956	numberNotTrusted = 0;
1957	numberStrongDone = 0;
1958	numberUnfinished = 0;
1959	numberStrongInfeasible = 0;
1960	numberStrongIterations = 0;
1961	#ifdef RANGING
1962	int * which = objectMark + numberObjects + 1;
1963	int neededPenalties;
1964	int optionalPenalties;
1965	#endif
1966	// We may go round this loop three times (only if we think we have solution)
1967	for (int iPass = 0; iPass < 3; iPass++) {
1968
1969	// Some objects may compute an estimate of best solution from here
1970	estimatedDegradation = 0.0;
1971	numberUnsatisfied_ = 0;
1972	// initialize sum of "infeasibilities"
1973	sumInfeasibilities_ = 0.0;
1974	int bestPriority = COIN_INT_MAX;
1975	#ifdef JJF_ZERO
1976	int number01 = 0;
1977	const cliqueEntry * entry = NULL;
1978	const int * toZero = NULL;
1979	const int * toOne = NULL;
1980	const int * backward = NULL;
1981	int numberUnsatisProbed = 0;
1982	int numberUnsatisNotProbed = 0; // 0-1
1983	if (probingInfo) {
1984	number01 = probingInfo->numberIntegers();
1985	entry = probingInfo->fixEntries();
1986	toZero = probingInfo->toZero();
1987	toOne = probingInfo->toOne();
1988	backward = probingInfo->backward();
1989	if (!toZero[number01] \|\| number01 < numberObjects \|\| true) {
1990	// no info
1991	probingInfo = NULL;
1992	}
1993	}
1994	#endif
1995	/*
1996	Scan for branching objects that indicate infeasibility. Choose candidates
1997	using priority as the first criteria, then integer infeasibility.
1998
1999	The algorithm is to fill the array with a set of good candidates (by
2000	infeasibility) with priority bestPriority. Finding a candidate with
2001	priority better (less) than bestPriority flushes the choice array. (This
2002	serves as initialization when the first candidate is found.)
2003
2004	*/
2005	numberToDo = 0;
2006	#ifdef RANGING
2007	neededPenalties = 0;
2008	optionalPenalties = numberObjects;
2009	#endif
2010	iBestNot = -1;
2011	double bestNot = 0.0;
2012	iBestGot = -1;
2013	best = 0.0;
2014	/* Problem type as set by user or found by analysis. This will be extended
2015	0 - not known
2016	1 - Set partitioning <=
2017	2 - Set partitioning ==
2018	3 - Set covering
2019	4 - all +- 1 or all +1 and odd
2020	*/
2021	int problemType = model->problemType();
2022	bool canDoOneHot = false;
2023	for (i = 0; i < numberObjects; i++) {
2024	OsiObject * object = model->modifiableObject(i);
2025	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2026	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2027	int preferredWay;
2028	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
2029	int priorityLevel = object->priority();
2030	if (hotstartSolution) {
2031	// we are doing hot start
2032	const CbcSimpleInteger * thisOne = dynamic_cast <const CbcSimpleInteger *> (object);
2033	if (thisOne) {
2034	int iColumn = thisOne->columnNumber();
2035	bool canDoThisHot = true;
2036	double targetValue = hotstartSolution[iColumn];
2037	if (saveUpper[iColumn] > saveLower[iColumn]) {
2038	double value = saveSolution[iColumn];
2039	if (hotstartPriorities)
2040	priorityLevel = hotstartPriorities[iColumn];
2041	//double originalLower = thisOne->originalLower();
2042	//double originalUpper = thisOne->originalUpper();
2043	// switch off if not possible
2044	if (targetValue >= saveLower[iColumn] && targetValue <= saveUpper[iColumn]) {
2045	/* priority outranks rest always if negative
2046	otherwise can be downgraded if at correct level.
2047	Infeasibility may be increased to choose 1.0 values first.
2048	choose one near wanted value
2049	*/
2050	if (fabs(value - targetValue) > integerTolerance) {
2051	//if (infeasibility>0.01)
2052	//infeasibility = fabs(1.0e6-fabs(value-targetValue));
2053	//else
2054	infeasibility = fabs(value - targetValue);
2055	//if (targetValue==1.0)
2056	//infeasibility += 1.0;
2057	if (value > targetValue) {
2058	preferredWay = -1;
2059	} else {
2060	preferredWay = 1;
2061	}
2062	priorityLevel = CoinAbs(priorityLevel);
2063	} else if (priorityLevel < 0) {
2064	priorityLevel = CoinAbs(priorityLevel);
2065	if (targetValue == saveLower[iColumn]) {
2066	infeasibility = integerTolerance + 1.0e-12;
2067	preferredWay = -1;
2068	} else if (targetValue == saveUpper[iColumn]) {
2069	infeasibility = integerTolerance + 1.0e-12;
2070	preferredWay = 1;
2071	} else {
2072	// can't
2073	priorityLevel += 10000000;
2074	canDoThisHot = false;
2075	}
2076	} else {
2077	priorityLevel += 10000000;
2078	canDoThisHot = false;
2079	}
2080	} else {
2081	// switch off if not possible
2082	canDoThisHot = false;
2083	}
2084	if (canDoThisHot)
2085	canDoOneHot = true;
2086	} else if (targetValue < saveLower[iColumn] \|\| targetValue > saveUpper[iColumn]) {
2087	}
2088	} else {
2089	priorityLevel += 10000000;
2090	}
2091	}
2092	#define ZERO_ONE 0
2093	#define ZERO_FAKE 1.0e20;
2094	#if ZERO_ONE==1
2095	// branch on 0-1 first (temp)
2096	if (fabs(saveSolution[dynamicObject->columnNumber()]) < 1.0)
2097	priorityLevel--;
2098	#endif
2099	#if ZERO_ONE==2
2100	if (fabs(saveSolution[dynamicObject->columnNumber()]) < 1.0)
2101	infeasibility *= ZERO_FAKE;
2102	#endif
2103	if (infeasibility) {
2104	int iColumn = numberColumns + i;
2105	bool gotDown = false;
2106	int numberThisDown = 0;
2107	bool gotUp = false;
2108	int numberThisUp = 0;
2109	double downGuess = object->downEstimate();
2110	double upGuess = object->upEstimate();
2111	if (dynamicObject) {
2112	// Use this object's numberBeforeTrust
2113	int numberBeforeTrust = dynamicObject->numberBeforeTrust();
2114	iColumn = dynamicObject->columnNumber();
2115	gotDown = false;
2116	numberThisDown = dynamicObject->numberTimesDown();
2117	if (numberThisDown >= numberBeforeTrust)
2118	gotDown = true;
2119	gotUp = false;
2120	numberThisUp = dynamicObject->numberTimesUp();
2121	if (numberThisUp >= numberBeforeTrust)
2122	gotUp = true;
2123	if (!depth_ && false) {
2124	// try closest to 0.5
2125	double part = saveSolution[iColumn] - floor(saveSolution[iColumn]);
2126	infeasibility = fabs(0.5 - part);
2127	}
2128	if (problemType > 0 && problemType < 4 && false) {
2129	// try closest to 0.5
2130	double part = saveSolution[iColumn] - floor(saveSolution[iColumn]);
2131	infeasibility = 0.5 - fabs(0.5 - part);
2132	}
2133	#ifdef JJF_ZERO
2134	if (probingInfo) {
2135	int iSeq = backward[iColumn];
2136	assert (iSeq >= 0);
2137	infeasibility = 1.0 + (toZero[iSeq+1] - toZero[iSeq]) +
2138	5.0 * CoinMin(toOne[iSeq] - toZero[iSeq], toZero[iSeq+1] - toOne[iSeq]);
2139	if (toZero[iSeq+1] > toZero[iSeq]) {
2140	numberUnsatisProbed++;
2141	} else {
2142	numberUnsatisNotProbed++;
2143	}
2144	}
2145	#endif
2146	} else {
2147	// see if SOS
2148	CbcSOS * sosObject =
2149	dynamic_cast <CbcSOS *>(object) ;
2150	if (sosObject) {
2151	gotDown = false;
2152	numberThisDown = sosObject->numberTimesDown();
2153	if (numberThisDown >= numberBeforeTrust)
2154	gotDown = true;
2155	gotUp = false;
2156	numberThisUp = sosObject->numberTimesUp();
2157	if (numberThisUp >= numberBeforeTrust)
2158	gotUp = true;
2159	} else {
2160	gotDown = true;
2161	numberThisDown = 999999;
2162	downGuess = 1.0e20;
2163	gotUp = true;
2164	numberThisUp = 999999;
2165	upGuess = 1.0e20;
2166	numberPassesLeft = 0;
2167	}
2168	}
2169	// Increase estimated degradation to solution
2170	estimatedDegradation += CoinMin(downGuess, upGuess);
2171	downEstimate[i] = downGuess;
2172	upEstimate[i] = upGuess;
2173	numberUnsatisfied_++;
2174	sumInfeasibilities_ += infeasibility;
2175	// Better priority? Flush choices.
2176	if (priorityLevel < bestPriority) {
2177	numberToDo = 0;
2178	bestPriority = priorityLevel;
2179	iBestGot = -1;
2180	best = 0.0;
2181	numberNotTrusted = 0;
2182	#ifdef RANGING
2183	neededPenalties = 0;
2184	optionalPenalties = numberObjects;
2185	#endif
2186	} else if (priorityLevel > bestPriority) {
2187	continue;
2188	}
2189	if (!gotUp \|\| !gotDown)
2190	numberNotTrusted++;
2191	// Check for suitability based on infeasibility.
2192	if ((gotDown && gotUp) && numberStrong > 0) {
2193	sort[numberToDo] = -infeasibility;
2194	if (infeasibility > best) {
2195	best = infeasibility;
2196	iBestGot = numberToDo;
2197	}
2198	#ifdef RANGING
2199	if (dynamicObject) {
2200	objectMark[--optionalPenalties] = numberToDo;
2201	which[optionalPenalties] = iColumn;
2202	}
2203	#endif
2204	} else {
2205	#ifdef RANGING
2206	if (dynamicObject) {
2207	objectMark[neededPenalties] = numberToDo;
2208	which[neededPenalties++] = iColumn;
2209	}
2210	#endif
2211	sort[numberToDo] = -10.0 * infeasibility;
2212	if (!(numberThisUp + numberThisDown))
2213	sort[numberToDo] *= 100.0; // make even more likely
2214	if (iColumn < numberColumns) {
2215	double part = saveSolution[iColumn] - floor(saveSolution[iColumn]);
2216	if (1.0 - fabs(part - 0.5) > bestNot) {
2217	iBestNot = numberToDo;
2218	bestNot = 1.0 - fabs(part - 0.5);
2219	}
2220	} else {
2221	// SOS
2222	if (-sort[numberToDo] > bestNot) {
2223	iBestNot = numberToDo;
2224	bestNot = -sort[numberToDo];
2225	}
2226	}
2227	}
2228	if (model->messageHandler()->logLevel() > 3) {
2229	printf("%d (%d) down %d %g up %d %g - infeas %g - sort %g solution %g\n",
2230	i, iColumn, numberThisDown, object->downEstimate(), numberThisUp, object->upEstimate(),
2231	infeasibility, sort[numberToDo], saveSolution[iColumn]);
2232	}
2233	whichObject[numberToDo++] = i;
2234	} else {
2235	// for debug
2236	downEstimate[i] = -1.0;
2237	upEstimate[i] = -1.0;
2238	}
2239	}
2240	if (numberUnsatisfied_) {
2241	//if (probingInfo&&false)
2242	//printf("nunsat %d, %d probed, %d other 0-1\n",numberUnsatisfied_,
2243	// numberUnsatisProbed,numberUnsatisNotProbed);
2244	// some infeasibilities - go to next steps
2245	if (!canDoOneHot && hotstartSolution) {
2246	// switch off as not possible
2247	hotstartSolution = NULL;
2248	model->setHotstartSolution(NULL, NULL);
2249	usefulInfo.hotstartSolution_ = NULL;
2250	}
2251	break;
2252	} else if (!iPass) {
2253	// may just need resolve
2254	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2255	//double newObjValue = solver->getObjSense()*solver->getObjValue();
2256	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2257	if (!solver->isProvenOptimal()) {
2258	// infeasible
2259	anyAction = -2;
2260	break;
2261	}
2262	// Double check looks OK - just look at rows with all integers
2263	if (model->allDynamic()) {
2264	double * solution = CoinCopyOfArray(saveSolution, numberColumns);
2265	for (int i = 0; i < numberColumns; i++) {
2266	if (model->isInteger(i))
2267	solution[i] = floor(solution[i] + 0.5);
2268	}
2269	int numberRows = solver->getNumRows();
2270	double * rowActivity = new double [numberRows];
2271	CoinZeroN(rowActivity, numberRows);
2272	solver->getMatrixByCol()->times(solution, rowActivity);
2273	//const double * element = model->solver()->getMatrixByCol()->getElements();
2274	const int * row = model->solver()->getMatrixByCol()->getIndices();
2275	const CoinBigIndex * columnStart = model->solver()->getMatrixByCol()->getVectorStarts();
2276	const int * columnLength = model->solver()->getMatrixByCol()->getVectorLengths();
2277	int nFree = 0;
2278	int nFreeNon = 0;
2279	int nFixedNon = 0;
2280	double mostAway = 0.0;
2281	int whichAway = -1;
2282	const double * columnLower = solver->getColLower();
2283	const double * columnUpper = solver->getColUpper();
2284	for (int i = 0; i < numberColumns; i++) {
2285	if (!model->isInteger(i)) {
2286	// mark rows as flexible
2287	CoinBigIndex start = columnStart[i];
2288	CoinBigIndex end = start + columnLength[i];
2289	for (CoinBigIndex j = start; j < end; j++) {
2290	int iRow = row[j];
2291	rowActivity[iRow] = COIN_DBL_MAX;
2292	}
2293	} else if (columnLower[i] < columnUpper[i]) {
2294	if (solution[i] != saveSolution[i]) {
2295	nFreeNon++;
2296	if (fabs(solution[i] - saveSolution[i]) > mostAway) {
2297	mostAway = fabs(solution[i] - saveSolution[i]);
2298	whichAway = i;
2299	}
2300	} else {
2301	nFree++;
2302	}
2303	} else if (solution[i] != saveSolution[i]) {
2304	nFixedNon++;
2305	}
2306	}
2307	const double * lower = solver->getRowLower();
2308	const double * upper = solver->getRowUpper();
2309	bool satisfied = true;
2310	for (int i = 0; i < numberRows; i++) {
2311	double value = rowActivity[i];
2312	if (value != COIN_DBL_MAX) {
2313	if (value > upper[i] + 1.0e-5 \|\| value < lower[i] - 1.0e-5) {
2314	satisfied = false;
2315	}
2316	}
2317	}
2318	delete [] rowActivity;
2319	delete [] solution;
2320	if (!satisfied) {
2321	#ifdef CLP_INVESTIGATE
2322	printf("%d free ok %d free off target %d fixed off target\n",
2323	nFree, nFreeNon, nFixedNon);
2324	#endif
2325	if (nFreeNon) {
2326	// try branching on these
2327	delete branch_;
2328	for (int i = 0; i < numberObjects; i++) {
2329	OsiObject * object = model->modifiableObject(i);
2330	CbcSimpleIntegerDynamicPseudoCost * obj =
2331	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2332	assert (obj);
2333	int iColumn = obj->columnNumber();
2334	if (iColumn == whichAway) {
2335	int preferredWay = (saveSolution[iColumn] > solution[iColumn])
2336	? -1 : +1;
2337	usefulInfo.integerTolerance_ = 0.0;
2338	branch_ = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
2339	break;
2340	}
2341	}
2342	anyAction = 0;
2343	break;
2344	}
2345	}
2346	}
2347	} else if (iPass == 1) {
2348	// looks like a solution - get paranoid
2349	bool roundAgain = false;
2350	// get basis
2351	CoinWarmStartBasis * ws = dynamic_cast<CoinWarmStartBasis*>(solver->getWarmStart());
2352	if (!ws)
2353	break;
2354	double tolerance;
2355	solver->getDblParam(OsiPrimalTolerance, tolerance);
2356	for (i = 0; i < numberColumns; i++) {
2357	double value = saveSolution[i];
2358	if (value < lower[i] - tolerance) {
2359	saveSolution[i] = lower[i];
2360	roundAgain = true;
2361	ws->setStructStatus(i, CoinWarmStartBasis::atLowerBound);
2362	} else if (value > upper[i] + tolerance) {
2363	saveSolution[i] = upper[i];
2364	roundAgain = true;
2365	ws->setStructStatus(i, CoinWarmStartBasis::atUpperBound);
2366	}
2367	}
2368	if (roundAgain) {
2369	// restore basis
2370	solver->setWarmStart(ws);
2371	solver->setColSolution(saveSolution);
2372	delete ws;
2373	bool takeHint;
2374	OsiHintStrength strength;
2375	solver->getHintParam(OsiDoDualInResolve, takeHint, strength);
2376	solver->setHintParam(OsiDoDualInResolve, false, OsiHintDo) ;
2377	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2378	//double newObjValue = solver->getObjSense()*solver->getObjValue();
2379	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2380	solver->setHintParam(OsiDoDualInResolve, takeHint, strength) ;
2381	if (!solver->isProvenOptimal()) {
2382	// infeasible
2383	anyAction = -2;
2384	break;
2385	}
2386	} else {
2387	delete ws;
2388	break;
2389	}
2390	}
2391	}
2392	if (anyAction == -2) {
2393	break;
2394	}
2395	// skip if solution
2396	if (!numberUnsatisfied_)
2397	break;
2398	int skipAll = (numberNotTrusted == 0 \|\| numberToDo == 1) ? 1 : 0;
2399	bool doneHotStart = false;
2400	//DEPRECATED_STRATEGYint searchStrategy = saveSearchStrategy>=0 ? (saveSearchStrategy%10) : -1;
2401	int searchStrategy = model->searchStrategy();
2402	// But adjust depending on ratio of iterations
2403	if (searchStrategy > 0) {
2404	if (numberBeforeTrust >= 5 && numberBeforeTrust <= 10) {
2405	if (searchStrategy != 2) {
2406	assert (searchStrategy == 1);
2407	if (depth_ > 5) {
2408	int numberIterations = model->getIterationCount();
2409	int numberStrongIterations = model->numberStrongIterations();
2410	if (numberStrongIterations > numberIterations + 10000) {
2411	searchStrategy = 2;
2412	skipAll = 1;
2413	} else if (numberStrongIterations*4 + 1000 < numberIterations) {
2414	searchStrategy = 3;
2415	skipAll = 0;
2416	}
2417	} else {
2418	searchStrategy = 3;
2419	skipAll = 0;
2420	}
2421	}
2422	}
2423	}
2424	// worth trying if too many times
2425	// Save basis
2426	CoinWarmStart * ws = NULL;
2427	// save limit
2428	int saveLimit = 0;
2429	solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit);
2430	if (!numberPassesLeft)
2431	skipAll = 1;
2432	if (!skipAll) {
2433	ws = solver->getWarmStart();
2434	int limit = 100;
2435	if (!saveStateOfSearch && saveLimit < limit && saveLimit == 100)
2436	solver->setIntParam(OsiMaxNumIterationHotStart, limit);
2437	}
2438	// Say which one will be best
2439	int whichChoice = 0;
2440	int bestChoice;
2441	if (iBestGot >= 0)
2442	bestChoice = iBestGot;
2443	else
2444	bestChoice = iBestNot;
2445	assert (bestChoice >= 0);
2446	// If we have hit max time don't do strong branching
2447	bool hitMaxTime = (model->getCurrentSeconds() >
2448	model->getDblParam(CbcModel::CbcMaximumSeconds));
2449	// also give up if we are looping round too much
2450	if (hitMaxTime \|\| numberPassesLeft <= 0 \|\| useShadow == 2) {
2451	int iObject = whichObject[bestChoice];
2452	OsiObject * object = model->modifiableObject(iObject);
2453	int preferredWay;
2454	object->infeasibility(&usefulInfo, preferredWay);
2455	CbcObject * obj =
2456	dynamic_cast <CbcObject *>(object) ;
2457	assert (obj);
2458	branch_ = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
2459	{
2460	CbcBranchingObject * branchObj =
2461	dynamic_cast <CbcBranchingObject *>(branch_) ;
2462	assert (branchObj);
2463	branchObj->way(preferredWay);
2464	}
2465	delete ws;
2466	ws = NULL;
2467	break;
2468	} else {
2469	// say do fast
2470	int easy = 1;
2471	if (!skipAll)
2472	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
2473	int iDo;
2474	#ifdef RANGING
2475	bool useRanging = model->allDynamic() && !skipAll;
2476	if (useRanging) {
2477	double currentObjective = solver->getObjValue() * solver->getObjSense();
2478	double gap = cutoff - currentObjective;
2479	// relax a bit
2480	gap *= 1.0000001;
2481	gap = CoinMax(1.0e-5, gap);
2482	// off penalties if too much
2483	double needed = neededPenalties;
2484	needed *= numberRows;
2485	if (numberNodes) {
2486	if (needed > 1.0e6) {
2487	neededPenalties = 0;
2488	} else if (gap < 1.0e5) {
2489	// maybe allow some not needed
2490	int extra = static_cast<int> ((1.0e6 - needed) / numberRows);
2491	int nStored = numberObjects - optionalPenalties;
2492	extra = CoinMin(extra, nStored);
2493	for (int i = 0; i < extra; i++) {
2494	objectMark[neededPenalties] = objectMark[optionalPenalties+i];
2495	which[neededPenalties++] = which[optionalPenalties+i];;
2496	}
2497	}
2498	}
2499	if (osiclp && neededPenalties) {
2500	assert (!doneHotStart);
2501	xPen += neededPenalties;
2502	which--;
2503	which[0] = neededPenalties;
2504	osiclp->passInRanges(which);
2505	// Mark hot start and get ranges
2506	if (kPass) {
2507	// until can work out why solution can go funny
2508	int save = osiclp->specialOptions();
2509	osiclp->setSpecialOptions(save \| 256);
2510	solver->markHotStart();
2511	osiclp->setSpecialOptions(save);
2512	} else {
2513	solver->markHotStart();
2514	}
2515	doneHotStart = true;
2516	xMark++;
2517	kPass++;
2518	osiclp->passInRanges(NULL);
2519	const double * downCost = osiclp->upRange();
2520	const double * upCost = osiclp->downRange();
2521	bool problemFeasible = true;
2522	int numberFixed = 0;
2523	for (int i = 0; i < neededPenalties; i++) {
2524	int j = objectMark[i];
2525	int iObject = whichObject[j];
2526	OsiObject * object = model->modifiableObject(iObject);
2527	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2528	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2529	// Use this object's numberBeforeTrust
2530	int numberBeforeTrust = dynamicObject->numberBeforeTrust();
2531	int iSequence = dynamicObject->columnNumber();
2532	double value = saveSolution[iSequence];
2533	value -= floor(value);
2534	double upPenalty = CoinMin(upCost[i], 1.0e110) * (1.0 - value);
2535	double downPenalty = CoinMin(downCost[i], 1.0e110) * value;
2536	int numberThisDown = dynamicObject->numberTimesDown();
2537	int numberThisUp = dynamicObject->numberTimesUp();
2538	if (!numberBeforeTrust) {
2539	// override
2540	downEstimate[iObject] = downPenalty;
2541	upEstimate[iObject] = upPenalty;
2542	double min1 = CoinMin(downEstimate[iObject],
2543	upEstimate[iObject]);
2544	double max1 = CoinMax(downEstimate[iObject],
2545	upEstimate[iObject]);
2546	min1 = 0.8 * min1 + 0.2 * max1;
2547	sort[j] = - min1;
2548	} else if (numberThisDown < numberBeforeTrust \|\|
2549	numberThisUp < numberBeforeTrust) {
2550	double invTrust = 1.0 / static_cast<double> (numberBeforeTrust);
2551	if (numberThisDown < numberBeforeTrust) {
2552	double fraction = numberThisDown * invTrust;
2553	downEstimate[iObject] = fraction * downEstimate[iObject] + (1.0 - fraction) * downPenalty;
2554	}
2555	if (numberThisUp < numberBeforeTrust) {
2556	double fraction = numberThisUp * invTrust;
2557	upEstimate[iObject] = fraction * upEstimate[iObject] + (1.0 - fraction) * upPenalty;
2558	}
2559	double min1 = CoinMin(downEstimate[iObject],
2560	upEstimate[iObject]);
2561	double max1 = CoinMax(downEstimate[iObject],
2562	upEstimate[iObject]);
2563	min1 = 0.8 * min1 + 0.2 * max1;
2564	min1 *= 10.0;
2565	if (!(numberThisDown + numberThisUp))
2566	min1 *= 100.0;
2567	sort[j] = - min1;
2568	}
2569	if (CoinMax(downPenalty, upPenalty) > gap) {
2570	COIN_DETAIL_PRINT(printf("gap %g object %d has down range %g, up %g\n",
2571	gap, i, downPenalty, upPenalty));
2572	//sort[j] -= 1.0e50; // make more likely to be chosen
2573	int number;
2574	if (downPenalty > gap) {
2575	number = dynamicObject->numberTimesDown();
2576	if (upPenalty > gap)
2577	problemFeasible = false;
2578	CbcBranchingObject * branch = dynamicObject->createCbcBranch(solver, &usefulInfo, 1);
2579	//branch->fix(solver,saveLower,saveUpper,1);
2580	delete branch;
2581	} else {
2582	number = dynamicObject->numberTimesUp();
2583	CbcBranchingObject * branch = dynamicObject->createCbcBranch(solver, &usefulInfo, 1);
2584	//branch->fix(solver,saveLower,saveUpper,-1);
2585	delete branch;
2586	}
2587	if (number >= numberBeforeTrust)
2588	dynamicObject->setNumberBeforeTrust(number + 1);
2589	numberFixed++;
2590	}
2591	if (!numberNodes)
2592	COIN_DETAIL_PRINT(printf("%d pen down ps %g -> %g up ps %g -> %g\n",
2593	iObject, downPenalty, downPenalty, upPenalty, upPenalty));
2594	}
2595	if (numberFixed && problemFeasible) {
2596	assert(doneHotStart);
2597	solver->unmarkHotStart();
2598	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2599	//double newObjValue = solver->getObjSense()*solver->getObjValue();
2600	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2601	solver->markHotStart();
2602	problemFeasible = solver->isProvenOptimal();
2603	}
2604	if (!problemFeasible) {
2605	COIN_DETAIL_PRINT(fprintf(stdout, "both ways infeas on ranging - code needed\n"));
2606	anyAction = -2;
2607	if (!choiceObject) {
2608	delete choice.possibleBranch;
2609	choice.possibleBranch = NULL;
2610	}
2611	//printf("Both infeasible for choice %d sequence %d\n",i,
2612	// model->object(choice.objectNumber)->columnNumber());
2613	// Delete the snapshot
2614	solver->unmarkHotStart();
2615	// back to normal
2616	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
2617	// restore basis
2618	solver->setWarmStart(ws);
2619	doneHotStart = false;
2620	delete ws;
2621	ws = NULL;
2622	break;
2623	}
2624	}
2625	}
2626	#endif /* RANGING */
2627	{
2628	int numberIterations = model->getIterationCount();
2629	//numberIterations += (model->numberExtraIterations()>>2);
2630	const int * strongInfo = model->strongInfo();
2631	//int numberDone = strongInfo[0]-strongInfo[3];
2632	int numberFixed = strongInfo[1] - strongInfo[4];
2633	int numberInfeasible = strongInfo[2] - strongInfo[5];
2634	assert (!strongInfo[3]);
2635	assert (!strongInfo[4]);
2636	assert (!strongInfo[5]);
2637	int numberStrongIterations = model->numberStrongIterations();
2638	int numberRows = solver->getNumRows();
2639	if (numberStrongIterations > numberIterations + CoinMin(100, 10*numberRows) && depth_ >= 4 && numberNodes > 100) {
2640	if (20numberInfeasible + 4numberFixed < numberNodes) {
2641	// Say never do
2642	skipAll = -1;
2643	}
2644	}
2645	}
2646	// make sure best will be first
2647	if (iBestGot >= 0)
2648	sort[iBestGot] = -COIN_DBL_MAX;
2649	// Actions 0 - exit for repeat, 1 resolve and try old choice,2 exit for continue
2650	if (anyAction)
2651	numberToDo = 0; // skip as we will be trying again
2652	// Sort
2653	CoinSort_2(sort, sort + numberToDo, whichObject);
2654	// Change in objective opposite infeasible
2655	double worstFeasible = 0.0;
2656	// Just first if strong off
2657	if (!numberStrong)
2658	numberToDo = CoinMin(numberToDo, 1);
2659	if (searchStrategy == 2)
2660	numberToDo = CoinMin(numberToDo, 20);
2661	iDo = 0;
2662	int saveLimit2;
2663	solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit2);
2664	bool doQuickly = false; // numberToDo>2*numberStrong;
2665	if (searchStrategy == 2)
2666	doQuickly = true;
2667	int numberTest = numberNotTrusted > 0 ? numberStrong : (numberStrong + 1) / 2;
2668	if (searchStrategy == 3) {
2669	// Previously decided we need strong
2670	doQuickly = false;
2671	numberTest = numberStrong;
2672	}
2673	// Try nearly always off
2674	if (skipAll >= 0) {
2675	if (searchStrategy < 2) {
2676	//if ((numberNodes%20)!=0) {
2677	if ((model->specialOptions()&8) == 0) {
2678	numberTest = 0;
2679	doQuickly = true;
2680	}
2681	//} else {
2682	//doQuickly=false;
2683	//numberTest=2*numberStrong;
2684	//skipAll=0;
2685	//}
2686	}
2687	} else {
2688	// Just take first
2689	doQuickly = true;
2690	numberTest = 1;
2691	}
2692	int testDepth = (skipAll >= 0) ? 8 : 4;
2693	if (depth_ < testDepth && numberStrong) {
2694	if (searchStrategy != 2) {
2695	doQuickly = false;
2696	int numberRows = solver->getNumRows();
2697	// whether to do this or not is important - think
2698	if (numberRows < 300 \|\| numberRows + numberColumns < 2500) {
2699	if (depth_ < 7)
2700	numberStrong = CoinMin(3 * numberStrong, numberToDo);
2701	if (!depth_)
2702	numberStrong = CoinMin(6 * numberStrong, numberToDo);
2703	}
2704	numberTest = numberStrong;
2705	skipAll = 0;
2706	}
2707	}
2708	// Do at least 5 strong
2709	if (numberColumns < 1000 && (depth_ < 15 \|\| numberNodes < 1000000))
2710	numberTest = CoinMax(numberTest, 5);
2711	if ((model->specialOptions()&8) == 0) {
2712	if (skipAll) {
2713	numberTest = 0;
2714	doQuickly = true;
2715	}
2716	} else {
2717	// do 5 as strong is fixing
2718	numberTest = CoinMax(numberTest, 5);
2719	}
2720	// see if switched off
2721	if (skipAll < 0) {
2722	numberTest = 0;
2723	doQuickly = true;
2724	}
2725	int realMaxHotIterations = 999999;
2726	if (skipAll < 0)
2727	numberToDo = 1;
2728	#ifdef DO_ALL_AT_ROOT
2729	if (!numberNodes) {
2730	printf("DOX %d test %d unsat %d - nobj %d\n",
2731	numberToDo, numberTest, numberUnsatisfied_,
2732	numberObjects);
2733	numberTest = numberToDo;
2734	}
2735	#endif
2736	for ( iDo = 0; iDo < numberToDo; iDo++) {
2737	int iObject = whichObject[iDo];
2738	OsiObject * object = model->modifiableObject(iObject);
2739	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
2740	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
2741	int iColumn = dynamicObject ? dynamicObject->columnNumber() : numberColumns + iObject;
2742	int preferredWay;
2743	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
2744	// may have become feasible
2745	if (!infeasibility)
2746	continue;
2747	#ifndef NDEBUG
2748	if (iColumn < numberColumns) {
2749	const double * solution = model->testSolution();
2750	assert (saveSolution[iColumn] == solution[iColumn]);
2751	}
2752	#endif
2753	CbcSimpleInteger * obj =
2754	dynamic_cast <CbcSimpleInteger *>(object) ;
2755	if (obj) {
2756	if (choiceObject) {
2757	obj->fillCreateBranch(choiceObject, &usefulInfo, preferredWay);
2758	choiceObject->setObject(dynamicObject);
2759	} else {
2760	choice.possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
2761	}
2762	} else {
2763	CbcObject * obj =
2764	dynamic_cast <CbcObject *>(object) ;
2765	assert (obj);
2766	choice.possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
2767	}
2768	// Save which object it was
2769	choice.objectNumber = iObject;
2770	choice.numIntInfeasUp = numberUnsatisfied_;
2771	choice.numIntInfeasDown = numberUnsatisfied_;
2772	choice.upMovement = upEstimate[iObject];
2773	choice.downMovement = downEstimate[iObject];
2774	assert (choice.upMovement >= 0.0);
2775	assert (choice.downMovement >= 0.0);
2776	choice.fix = 0; // say not fixed
2777	// see if can skip strong branching
2778	int canSkip = choice.possibleBranch->fillStrongInfo(choice);
2779	if ((numberTest <= 0 \|\| skipAll)) {
2780	if (iDo > 20) {
2781	#ifdef DO_ALL_AT_ROOT
2782	if (!numberNodes)
2783	printf("DOY test %d - iDo %d\n",
2784	numberTest, iDo);
2785	#endif
2786	if (!choiceObject) {
2787	delete choice.possibleBranch;
2788	choice.possibleBranch = NULL;
2789	}
2790	break; // give up anyway
2791	}
2792	}
2793	if (model->messageHandler()->logLevel() > 3 && numberBeforeTrust && dynamicObject)
2794	dynamicObject->print(1, choice.possibleBranch->value());
2795	if (skipAll < 0)
2796	canSkip = true;
2797	if (!canSkip) {
2798	if (!doneHotStart) {
2799	// Mark hot start
2800	doneHotStart = true;
2801	solver->markHotStart();
2802	if (!solver->isProvenOptimal()) {
2803	skipAll=-2;
2804	canSkip = true;
2805	}
2806	xMark++;
2807	}
2808	}
2809	if (!canSkip) {
2810	numberTest--;
2811	// just do a few
2812	if (searchStrategy == 2)
2813	solver->setIntParam(OsiMaxNumIterationHotStart, 10);
2814	double objectiveChange ;
2815	double newObjectiveValue = 1.0e100;
2816	int j;
2817	// status is 0 finished, 1 infeasible and other
2818	int iStatus;
2819	/*
2820	Try the down direction first. (Specify the initial branching alternative as
2821	down with a call to way(-1). Each subsequent call to branch() performs the
2822	specified branch and advances the branch object state to the next branch
2823	alternative.)
2824	*/
2825	choice.possibleBranch->way(-1) ;
2826	choice.possibleBranch->branch() ;
2827	solver->solveFromHotStart() ;
2828	bool needHotStartUpdate = false;
2829	numberStrongDone++;
2830	numberStrongIterations += solver->getIterationCount();
2831	/*
2832	We now have an estimate of objective degradation that we can use for strong
2833	branching. If we're over the cutoff, the variable is monotone up.
2834	If we actually made it to optimality, check for a solution, and if we have
2835	a good one, call setBestSolution to process it. Note that this may reduce the
2836	cutoff, so we check again to see if we can declare this variable monotone.
2837	*/
2838	if (solver->isProvenOptimal())
2839	iStatus = 0; // optimal
2840	else if (solver->isIterationLimitReached()
2841	&& !solver->isDualObjectiveLimitReached()) {
2842	iStatus = 2; // unknown
2843	} else {
2844	iStatus = 1; // infeasible
2845	}
2846	if (iStatus != 2 && solver->getIterationCount() >
2847	realMaxHotIterations)
2848	numberUnfinished++;
2849	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
2850	choice.numItersDown = solver->getIterationCount();
2851	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0);
2852	// Update branching information if wanted
2853	CbcBranchingObject * cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
2854	if (cbcobj) {
2855	CbcObject * object = cbcobj->object();
2856	assert (object) ;
2857	CbcObjectUpdateData update = object->createUpdateInformation(solver, this, cbcobj);
2858	update.objectNumber_ = choice.objectNumber;
2859	model->addUpdateInformation(update);
2860	} else {
2861	decision->updateInformation( solver, this);
2862	}
2863	if (!iStatus) {
2864	choice.finishedDown = true ;
2865	if (newObjectiveValue >= cutoff) {
2866	objectiveChange = 1.0e100; // say infeasible
2867	numberStrongInfeasible++;
2868	} else {
2869	#define CBCNODE_TIGHTEN_BOUNDS
2870	#ifdef CBCNODE_TIGHTEN_BOUNDS
2871	// Can we tighten bounds?
2872	if (iColumn < numberColumns && cutoff < 1.0e20
2873	&& objectiveChange > 1.0e-5) {
2874	double value = saveSolution[iColumn];
2875	double down = value - floor(value);
2876	double changePer = objectiveChange / (down + 1.0e-7);
2877	double distance = (cutoff - objectiveValue_) / changePer;
2878	distance += 1.0e-3;
2879	if (distance < 5.0) {
2880	double newLower = ceil(value - distance);
2881	if (newLower > saveLower[iColumn]) {
2882	//printf("Could increase lower bound on %d from %g to %g\n",
2883	// iColumn,saveLower[iColumn],newLower);
2884	saveLower[iColumn] = newLower;
2885	solver->setColLower(iColumn, newLower);
2886	}
2887	}
2888	}
2889	#endif
2890	// See if integer solution
2891	if (model->feasibleSolution(choice.numIntInfeasDown,
2892	choice.numObjInfeasDown)
2893	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
2894	if (auxiliaryInfo->solutionAddsCuts()) {
2895	needHotStartUpdate = true;
2896	solver->unmarkHotStart();
2897	}
2898	model->setLogLevel(saveLogLevel);
2899	model->setBestSolution(CBC_STRONGSOL,
2900	newObjectiveValue,
2901	solver->getColSolution()) ;
2902	if (needHotStartUpdate) {
2903	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2904	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
2905	//objectiveValue_ = CoinMax(objectiveValue_,newObjectiveValue);
2906	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0);
2907	model->feasibleSolution(choice.numIntInfeasDown,
2908	choice.numObjInfeasDown);
2909	}
2910	model->setLastHeuristic(NULL);
2911	model->incrementUsed(solver->getColSolution());
2912	cutoff = model->getCutoff();
2913	if (newObjectiveValue >= cutoff) { // new cutoff
2914	objectiveChange = 1.0e100 ;
2915	numberStrongInfeasible++;
2916	}
2917	}
2918	}
2919	} else if (iStatus == 1) {
2920	objectiveChange = 1.0e100 ;
2921	numberStrongInfeasible++;
2922	} else {
2923	// Can't say much as we did not finish
2924	choice.finishedDown = false ;
2925	numberUnfinished++;
2926	}
2927	choice.downMovement = objectiveChange ;
2928
2929	// restore bounds
2930	for ( j = 0; j < numberColumns; j++) {
2931	if (saveLower[j] != lower[j])
2932	solver->setColLower(j, saveLower[j]);
2933	if (saveUpper[j] != upper[j])
2934	solver->setColUpper(j, saveUpper[j]);
2935	}
2936	if (needHotStartUpdate) {
2937	needHotStartUpdate = false;
2938	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2939	//double newObjValue = solver->getObjSense()*solver->getObjValue();
2940	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2941	//we may again have an integer feasible solution
2942	int numberIntegerInfeasibilities;
2943	int numberObjectInfeasibilities;
2944	if (model->feasibleSolution(
2945	numberIntegerInfeasibilities,
2946	numberObjectInfeasibilities)) {
2947	#ifdef BONMIN
2948	//In this case node has become integer feasible, let us exit the loop
2949	std::cout << "Node has become integer feasible" << std::endl;
2950	numberUnsatisfied_ = 0;
2951	break;
2952	#endif
2953	double objValue = solver->getObjValue();
2954	model->setLogLevel(saveLogLevel);
2955	model->setBestSolution(CBC_STRONGSOL,
2956	objValue,
2957	solver->getColSolution()) ;
2958	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
2959	//double newObjValue = solver->getObjSense()*solver->getObjValue();
2960	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
2961	cutoff = model->getCutoff();
2962	}
2963	solver->markHotStart();
2964	}
2965	#ifdef DO_ALL_AT_ROOT
2966	if (!numberNodes)
2967	printf("Down on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
2968	choice.objectNumber, iStatus, newObjectiveValue, choice.numItersDown,
2969	choice.downMovement, choice.finishedDown, choice.numIntInfeasDown,
2970	choice.numObjInfeasDown);
2971	#endif
2972
2973	// repeat the whole exercise, forcing the variable up
2974	choice.possibleBranch->branch();
2975	solver->solveFromHotStart() ;
2976	numberStrongDone++;
2977	numberStrongIterations += solver->getIterationCount();
2978	/*
2979	We now have an estimate of objective degradation that we can use for strong
2980	branching. If we're over the cutoff, the variable is monotone up.
2981	If we actually made it to optimality, check for a solution, and if we have
2982	a good one, call setBestSolution to process it. Note that this may reduce the
2983	cutoff, so we check again to see if we can declare this variable monotone.
2984	*/
2985	if (solver->isProvenOptimal())
2986	iStatus = 0; // optimal
2987	else if (solver->isIterationLimitReached()
2988	&& !solver->isDualObjectiveLimitReached()) {
2989	iStatus = 2; // unknown
2990	} else {
2991	iStatus = 1; // infeasible
2992	}
2993	if (iStatus != 2 && solver->getIterationCount() >
2994	realMaxHotIterations)
2995	numberUnfinished++;
2996	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
2997	choice.numItersUp = solver->getIterationCount();
2998	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0);
2999	// Update branching information if wanted
3000	cbcobj = dynamic_cast<CbcBranchingObject *> (choice.possibleBranch);
3001	if (cbcobj) {
3002	CbcObject * object = cbcobj->object();
3003	assert (object) ;
3004	CbcObjectUpdateData update = object->createUpdateInformation(solver, this, cbcobj);
3005	update.objectNumber_ = choice.objectNumber;
3006	model->addUpdateInformation(update);
3007	} else {
3008	decision->updateInformation( solver, this);
3009	}
3010	if (!iStatus) {
3011	choice.finishedUp = true ;
3012	if (newObjectiveValue >= cutoff) {
3013	objectiveChange = 1.0e100; // say infeasible
3014	numberStrongInfeasible++;
3015	} else {
3016	#ifdef CBCNODE_TIGHTEN_BOUNDS
3017	// Can we tighten bounds?
3018	if (iColumn < numberColumns && cutoff < 1.0e20
3019	&& objectiveChange > 1.0e-5) {
3020	double value = saveSolution[iColumn];
3021	double up = ceil(value) - value;
3022	double changePer = objectiveChange / (up + 1.0e-7);
3023	double distance = (cutoff - objectiveValue_) / changePer;
3024	distance += 1.0e-3;
3025	if (distance < 5.0) {
3026	double newUpper = floor(value + distance);
3027	if (newUpper < saveUpper[iColumn]) {
3028	//printf("Could decrease upper bound on %d from %g to %g\n",
3029	// iColumn,saveUpper[iColumn],newUpper);
3030	saveUpper[iColumn] = newUpper;
3031	solver->setColUpper(iColumn, newUpper);
3032	}
3033	}
3034	}
3035	#endif
3036	// See if integer solution
3037	if (model->feasibleSolution(choice.numIntInfeasUp,
3038	choice.numObjInfeasUp)
3039	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
3040	#ifdef BONMIN
3041	std::cout << "Node has become integer feasible" << std::endl;
3042	numberUnsatisfied_ = 0;
3043	break;
3044	#endif
3045	if (auxiliaryInfo->solutionAddsCuts()) {
3046	needHotStartUpdate = true;
3047	solver->unmarkHotStart();
3048	}
3049	model->setLogLevel(saveLogLevel);
3050	model->setBestSolution(CBC_STRONGSOL,
3051	newObjectiveValue,
3052	solver->getColSolution()) ;
3053	if (choice.finishedDown) {
3054	double cutoff = model->getCutoff();
3055	double downObj = objectiveValue_
3056	+ choice.downMovement ;
3057	if (downObj >= cutoff) {
3058	choice.downMovement = 1.0e100 ;
3059	numberStrongInfeasible++;
3060	}
3061
3062	}
3063	if (needHotStartUpdate) {
3064	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3065	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3066	//objectiveValue_ = CoinMax(objectiveValue_,newObjectiveValue);
3067	objectiveChange = CoinMax(newObjectiveValue - objectiveValue_, 0.0);
3068	model->feasibleSolution(choice.numIntInfeasDown,
3069	choice.numObjInfeasDown);
3070	}
3071	model->setLastHeuristic(NULL);
3072	model->incrementUsed(solver->getColSolution());
3073	cutoff = model->getCutoff();
3074	if (newObjectiveValue >= cutoff) { // new cutoff
3075	objectiveChange = 1.0e100 ;
3076	numberStrongInfeasible++;
3077	}
3078	}
3079	}
3080	} else if (iStatus == 1) {
3081	objectiveChange = 1.0e100 ;
3082	numberStrongInfeasible++;
3083	} else {
3084	// Can't say much as we did not finish
3085	choice.finishedUp = false ;
3086	numberUnfinished++;
3087	}
3088	choice.upMovement = objectiveChange ;
3089
3090	// restore bounds
3091	for ( j = 0; j < numberColumns; j++) {
3092	if (saveLower[j] != lower[j])
3093	solver->setColLower(j, saveLower[j]);
3094	if (saveUpper[j] != upper[j])
3095	solver->setColUpper(j, saveUpper[j]);
3096	}
3097	if (needHotStartUpdate) {
3098	needHotStartUpdate = false;
3099	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3100	//double newObjValue = solver->getObjSense()*solver->getObjValue();
3101	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3102	//we may again have an integer feasible solution
3103	int numberIntegerInfeasibilities;
3104	int numberObjectInfeasibilities;
3105	if (model->feasibleSolution(
3106	numberIntegerInfeasibilities,
3107	numberObjectInfeasibilities)) {
3108	double objValue = solver->getObjValue();
3109	model->setLogLevel(saveLogLevel);
3110	model->setBestSolution(CBC_STRONGSOL,
3111	objValue,
3112	solver->getColSolution()) ;
3113	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3114	//double newObjValue = solver->getObjSense()*solver->getObjValue();
3115	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3116	cutoff = model->getCutoff();
3117	}
3118	solver->markHotStart();
3119	}
3120
3121	#ifdef DO_ALL_AT_ROOT
3122	if (!numberNodes)
3123	printf("Up on %d, status is %d, obj %g its %d cost %g finished %d inf %d infobj %d\n",
3124	choice.objectNumber, iStatus, newObjectiveValue, choice.numItersUp,
3125	choice.upMovement, choice.finishedUp, choice.numIntInfeasUp,
3126	choice.numObjInfeasUp);
3127	#endif
3128	}
3129
3130	solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit2);
3131	/*
3132	End of evaluation for this candidate variable. Possibilities are:
3133	* Both sides below cutoff; this variable is a candidate for branching.
3134	* Both sides infeasible or above the objective cutoff: no further action
3135	here. Break from the evaluation loop and assume the node will be purged
3136	by the caller.
3137	* One side below cutoff: Install the branch (i.e., fix the variable). Break
3138	from the evaluation loop and assume the node will be reoptimised by the
3139	caller.
3140	*/
3141	// reset
3142	choice.possibleBranch->resetNumberBranchesLeft();
3143	if (choice.upMovement < 1.0e100) {
3144	if (choice.downMovement < 1.0e100) {
3145	// In case solution coming in was odd
3146	choice.upMovement = CoinMax(0.0, choice.upMovement);
3147	choice.downMovement = CoinMax(0.0, choice.downMovement);
3148	#if ZERO_ONE==2
3149	// branch on 0-1 first (temp)
3150	if (fabs(choice.possibleBranch->value()) < 1.0) {
3151	choice.upMovement *= ZERO_FAKE;
3152	choice.downMovement *= ZERO_FAKE;
3153	}
3154	#endif
3155	// feasible - see which best
3156	if (!canSkip) {
3157	if (model->messageHandler()->logLevel() > 3)
3158	printf("sort %g downest %g upest %g ", sort[iDo], downEstimate[iObject],
3159	upEstimate[iObject]);
3160	model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3161	<< iObject << iColumn
3162	<< choice.downMovement << choice.numIntInfeasDown
3163	<< choice.upMovement << choice.numIntInfeasUp
3164	<< choice.possibleBranch->value()
3165	<< CoinMessageEol;
3166	}
3167	int betterWay;
3168	{
3169	CbcBranchingObject * branchObj =
3170	dynamic_cast <CbcBranchingObject *>(branch_) ;
3171	if (branch_)
3172	assert (branchObj);
3173	betterWay = decision->betterBranch(choice.possibleBranch,
3174	branchObj,
3175	choice.upMovement,
3176	choice.numIntInfeasUp ,
3177	choice.downMovement,
3178	choice.numIntInfeasDown );
3179	}
3180	if (betterWay) {
3181	// C) create branching object
3182	if (choiceObject) {
3183	delete branch_;
3184	branch_ = choice.possibleBranch->clone();
3185	} else {
3186	delete branch_;
3187	branch_ = choice.possibleBranch;
3188	choice.possibleBranch = NULL;
3189	}
3190	{
3191	CbcBranchingObject * branchObj =
3192	dynamic_cast <CbcBranchingObject *>(branch_) ;
3193	assert (branchObj);
3194	//branchObj->way(preferredWay);
3195	branchObj->way(betterWay);
3196	}
3197	bestChoice = choice.objectNumber;
3198	whichChoice = iDo;
3199	if (numberStrong <= 1) {
3200	delete ws;
3201	ws = NULL;
3202	break;
3203	}
3204	} else {
3205	if (!choiceObject) {
3206	delete choice.possibleBranch;
3207	choice.possibleBranch = NULL;
3208	}
3209	if (iDo >= 2*numberStrong) {
3210	delete ws;
3211	ws = NULL;
3212	break;
3213	}
3214	if (!dynamicObject \|\| dynamicObject->numberTimesUp() > 1) {
3215	if (iDo - whichChoice >= numberStrong) {
3216	if (!choiceObject) {
3217	delete choice.possibleBranch;
3218	choice.possibleBranch = NULL;
3219	}
3220	break; // give up
3221	}
3222	} else {
3223	if (iDo - whichChoice >= 2*numberStrong) {
3224	delete ws;
3225	ws = NULL;
3226	if (!choiceObject) {
3227	delete choice.possibleBranch;
3228	choice.possibleBranch = NULL;
3229	}
3230	break; // give up
3231	}
3232	}
3233	}
3234	} else {
3235	// up feasible, down infeasible
3236	anyAction = -1;
3237	worstFeasible = CoinMax(worstFeasible, choice.upMovement);
3238	model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3239	<< iObject << iColumn
3240	<< choice.downMovement << choice.numIntInfeasDown
3241	<< choice.upMovement << choice.numIntInfeasUp
3242	<< choice.possibleBranch->value()
3243	<< CoinMessageEol;
3244	//printf("Down infeasible for choice %d sequence %d\n",i,
3245	// model->object(choice.objectNumber)->columnNumber());
3246	choice.fix = 1;
3247	numberToFix++;
3248	choice.possibleBranch->fix(solver, saveLower, saveUpper, 1);
3249	if (!choiceObject) {
3250	delete choice.possibleBranch;
3251	choice.possibleBranch = NULL;
3252	} else {
3253	//choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
3254	choice.possibleBranch = choiceObject;
3255	}
3256	assert(doneHotStart);
3257	solver->unmarkHotStart();
3258	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3259	//double newObjValue = solver->getObjSense()*solver->getObjValue();
3260	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3261	solver->markHotStart();
3262	// may be infeasible (if other way stopped on iterations)
3263	if (!solver->isProvenOptimal()\|\|solver->isDualObjectiveLimitReached()) {
3264	// neither side feasible
3265	anyAction = -2;
3266	if (!choiceObject) {
3267	delete choice.possibleBranch;
3268	choice.possibleBranch = NULL;
3269	}
3270	//printf("Both infeasible for choice %d sequence %d\n",i,
3271	// model->object(choice.objectNumber)->columnNumber());
3272	delete ws;
3273	ws = NULL;
3274	break;
3275	}
3276	}
3277	} else {
3278	if (choice.downMovement < 1.0e100) {
3279	// down feasible, up infeasible
3280	anyAction = -1;
3281	worstFeasible = CoinMax(worstFeasible, choice.downMovement);
3282	model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3283	<< iObject << iColumn
3284	<< choice.downMovement << choice.numIntInfeasDown
3285	<< choice.upMovement << choice.numIntInfeasUp
3286	<< choice.possibleBranch->value()
3287	<< CoinMessageEol;
3288	choice.fix = -1;
3289	numberToFix++;
3290	choice.possibleBranch->fix(solver, saveLower, saveUpper, -1);
3291	if (!choiceObject) {
3292	delete choice.possibleBranch;
3293	choice.possibleBranch = NULL;
3294	} else {
3295	//choiceObject = new CbcDynamicPseudoCostBranchingObject(*choiceObject);
3296	choice.possibleBranch = choiceObject;
3297	}
3298	assert(doneHotStart);
3299	solver->unmarkHotStart();
3300	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper);
3301	//double newObjValue = solver->getObjSense()*solver->getObjValue();
3302	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3303	solver->markHotStart();
3304	// may be infeasible (if other way stopped on iterations)
3305	if (!solver->isProvenOptimal()\|\|solver->isDualObjectiveLimitReached()) {
3306	// neither side feasible
3307	anyAction = -2;
3308	if (!choiceObject) {
3309	delete choice.possibleBranch;
3310	choice.possibleBranch = NULL;
3311	}
3312	delete ws;
3313	ws = NULL;
3314	break;
3315	}
3316	} else {
3317	// neither side feasible
3318	anyAction = -2;
3319	if (!choiceObject) {
3320	delete choice.possibleBranch;
3321	choice.possibleBranch = NULL;
3322	}
3323	delete ws;
3324	ws = NULL;
3325	break;
3326	}
3327	}
3328	// Check max time
3329	hitMaxTime = (model->getCurrentSeconds() >
3330	model->getDblParam(CbcModel::CbcMaximumSeconds));
3331	if (hitMaxTime) {
3332	// make sure rest are fast
3333	doQuickly = true;
3334	for ( int jDo = iDo + 1; jDo < numberToDo; jDo++) {
3335	int iObject = whichObject[iDo];
3336	OsiObject * object = model->modifiableObject(iObject);
3337	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3338	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3339	if (dynamicObject)
3340	dynamicObject->setNumberBeforeTrust(0);
3341	}
3342	numberTest = 0;
3343	}
3344	if (!choiceObject) {
3345	delete choice.possibleBranch;
3346	}
3347	}
3348	if (model->messageHandler()->logLevel() > 3) {
3349	if (anyAction == -2) {
3350	printf("infeasible\n");
3351	} else if (anyAction == -1) {
3352	printf("%d fixed AND choosing %d iDo %d iChosenWhen %d numberToDo %d\n", numberToFix, bestChoice,
3353	iDo, whichChoice, numberToDo);
3354	} else {
3355	int iObject = whichObject[whichChoice];
3356	OsiObject * object = model->modifiableObject(iObject);
3357	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3358	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3359	if (dynamicObject) {
3360	int iColumn = dynamicObject->columnNumber();
3361	printf("choosing %d (column %d) iChosenWhen %d numberToDo %d\n", bestChoice,
3362	iColumn, whichChoice, numberToDo);
3363	}
3364	}
3365	}
3366	if (doneHotStart) {
3367	// Delete the snapshot
3368	solver->unmarkHotStart();
3369	// back to normal
3370	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
3371	// restore basis
3372	solver->setWarmStart(ws);
3373	}
3374	solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit);
3375	// Unless infeasible we will carry on
3376	// But we could fix anyway
3377	if (numberToFix && !hitMaxTime) {
3378	if (anyAction != -2) {
3379	// apply and take off
3380	bool feasible = true;
3381	// can do quick optimality check
3382	int easy = 2;
3383	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, &easy) ;
3384	model->resolve(NULL, 11, saveSolution, saveLower, saveUpper) ;
3385	//double newObjValue = solver->getObjSense()*solver->getObjValue();
3386	//objectiveValue_ = CoinMax(objectiveValue_,newObjValue);
3387	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
3388	feasible = solver->isProvenOptimal();
3389	if (feasible) {
3390	anyAction = 0;
3391	} else {
3392	anyAction = -2;
3393	finished = true;
3394	}
3395	}
3396	}
3397	// If fixed then round again
3398	// See if candidate still possible
3399	if (branch_) {
3400	const OsiObject * object = model->object(bestChoice);
3401	int preferredWay;
3402	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
3403	if (!infeasibility) {
3404	// take out
3405	delete branch_;
3406	branch_ = NULL;
3407	} else {
3408	CbcBranchingObject * branchObj =
3409	dynamic_cast <CbcBranchingObject *>(branch_) ;
3410	assert (branchObj);
3411	branchObj->way(preferredWay);
3412	#ifdef CBCNODE_TIGHTEN_BOUNDS
3413	bool fixed = branchObj->tighten(solver);
3414	if (fixed) {
3415	//printf("Variable now fixed!\n");
3416	// take out
3417	delete branch_;
3418	branch_ = NULL;
3419	}
3420	#endif
3421	}
3422	}
3423	if (!branch_ && anyAction != -2 && !hitMaxTime) {
3424	finished = false;
3425	}
3426	delete ws;
3427	}
3428	}
3429	// update number of strong iterations etc
3430	model->incrementStrongInfo(numberStrongDone, numberStrongIterations,
3431	anyAction == -2 ? 0 : numberToFix, anyAction == -2);
3432	if (model->searchStrategy() == -1) {
3433	#ifndef COIN_DEVELOP
3434	if (solver->messageHandler()->logLevel() > 1)
3435	#endif
3436	printf("%d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
3437	numberStrongDone, numberStrongIterations, numberStrongInfeasible, numberUnfinished,
3438	numberNotTrusted);
3439	// decide what to do
3440	int strategy = 1;
3441	if (((numberUnfinished*4 > numberStrongDone &&
3442	numberStrongInfeasible*40 < numberStrongDone) \|\|
3443	numberStrongInfeasible < 0) && model->numberStrong() < 10 && model->numberBeforeTrust() <= 20 && model->numberObjects() > CoinMax(1000, solver->getNumRows())) {
3444	strategy = 2;
3445	#ifdef COIN_DEVELOP
3446	//if (model->logLevel()>1)
3447	printf("going to strategy 2\n");
3448	#endif
3449	// Weaken
3450	model->setNumberStrong(2);
3451	model->setNumberBeforeTrust(1);
3452	model->synchronizeNumberBeforeTrust();
3453	}
3454	if (numberNodes)
3455	strategy = 1; // should only happen after hot start
3456	model->setSearchStrategy(strategy);
3457	} else if (numberStrongDone) {
3458	//printf("%d strongB, %d iters, %d inf, %d not finished, %d not trusted\n",
3459	// numberStrongDone,numberStrongIterations,numberStrongInfeasible,numberUnfinished,
3460	// numberNotTrusted);
3461	}
3462	if (model->searchStrategy() == 1 && numberNodes > 500 && numberNodes < -510) {
3463	#ifndef COIN_DEVELOP
3464	if (solver->messageHandler()->logLevel() > 1)
3465	#endif
3466	printf("after %d nodes - %d strong, %d iters, %d inf, %d not finished, %d not trusted\n",
3467	numberNodes, numberStrongDone, numberStrongIterations, numberStrongInfeasible, numberUnfinished,
3468	numberNotTrusted);
3469	// decide what to do
3470	if (numberUnfinished*10 > numberStrongDone + 1 \|\|
3471	!numberStrongInfeasible) {
3472	COIN_DETAIL_PRINT(printf("going to strategy 2\n"));
3473	// Weaken
3474	model->setNumberStrong(2);
3475	model->setNumberBeforeTrust(1);
3476	model->synchronizeNumberBeforeTrust();
3477	model->setSearchStrategy(2);
3478	}
3479	}
3480	if (numberUnfinished*10 < numberStrongDone &&
3481	numberStrongIterations*20 < model->getIterationCount()) {
3482	//printf("increasing trust\n");
3483	model->synchronizeNumberBeforeTrust(2);
3484	}
3485
3486	// Set guessed solution value
3487	guessedObjectiveValue_ = objectiveValue_ + estimatedDegradation;
3488
3489	/*
3490	Cleanup, then we're finished
3491	*/
3492	if (!model->branchingMethod())
3493	delete decision;
3494
3495	delete choiceObject;
3496	delete [] sort;
3497	delete [] whichObject;
3498	#ifdef RANGING
3499	delete [] objectMark;
3500	#endif
3501	delete [] saveLower;
3502	delete [] saveUpper;
3503	delete [] upEstimate;
3504	delete [] downEstimate;
3505	# ifdef COIN_HAS_CLP
3506	if (osiclp) {
3507	osiclp->setSpecialOptions(saveClpOptions);
3508	}
3509	# endif
3510	// restore solution
3511	solver->setColSolution(saveSolution);
3512	model->reserveCurrentSolution(saveSolution);
3513	delete [] saveSolution;
3514	model->setStateOfSearch(saveStateOfSearch);
3515	model->setLogLevel(saveLogLevel);
3516	// delete extra regions
3517	if (usefulInfo.usefulRegion_) {
3518	delete [] usefulInfo.usefulRegion_;
3519	delete [] usefulInfo.indexRegion_;
3520	delete [] usefulInfo.pi_;
3521	usefulInfo.usefulRegion_ = NULL;
3522	usefulInfo.indexRegion_ = NULL;
3523	usefulInfo.pi_ = NULL;
3524	}
3525	useShadow = model->moreSpecialOptions() & 7;
3526	if ((useShadow == 5 && model->getSolutionCount()) \|\| useShadow == 6) {
3527	// zap pseudo shadow prices
3528	model->pseudoShadow(-1);
3529	// and switch off
3530	model->setMoreSpecialOptions(model->moreSpecialOptions()&(~1023));
3531	}
3532	return anyAction;
3533	}
3534	int CbcNode::analyze (CbcModel model, double results)
3535	{
3536	int i;
3537	int numberIterationsAllowed = model->numberAnalyzeIterations();
3538	OsiSolverInterface * solver = model->solver();
3539	objectiveValue_ = solver->getObjSense() * solver->getObjValue();
3540	double cutoff = model->getCutoff();
3541	const double * lower = solver->getColLower();
3542	const double * upper = solver->getColUpper();
3543	const double * dj = solver->getReducedCost();
3544	int numberObjects = model->numberObjects();
3545	int numberColumns = model->getNumCols();
3546	// Initialize arrays
3547	int numberIntegers = model->numberIntegers();
3548	int * back = new int[numberColumns];
3549	const int * integerVariable = model->integerVariable();
3550	for (i = 0; i < numberColumns; i++)
3551	back[i] = -1;
3552	// What results is
3553	double * newLower = results;
3554	double * objLower = newLower + numberIntegers;
3555	double * newUpper = objLower + numberIntegers;
3556	double * objUpper = newUpper + numberIntegers;
3557	for (i = 0; i < numberIntegers; i++) {
3558	int iColumn = integerVariable[i];
3559	back[iColumn] = i;
3560	newLower[i] = 0.0;
3561	objLower[i] = -COIN_DBL_MAX;
3562	newUpper[i] = 0.0;
3563	objUpper[i] = -COIN_DBL_MAX;
3564	}
3565	double * saveUpper = new double[numberColumns];
3566	double * saveLower = new double[numberColumns];
3567	int anyAction = 0;
3568	// Save solution in case heuristics need good solution later
3569
3570	double * saveSolution = new double[numberColumns];
3571	memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
3572	model->reserveCurrentSolution(saveSolution);
3573	for (i = 0; i < numberColumns; i++) {
3574	saveLower[i] = lower[i];
3575	saveUpper[i] = upper[i];
3576	}
3577	// Get arrays to sort
3578	double * sort = new double[numberObjects];
3579	int * whichObject = new int[numberObjects];
3580	int numberToFix = 0;
3581	int numberToDo = 0;
3582	double integerTolerance =
3583	model->getDblParam(CbcModel::CbcIntegerTolerance);
3584	// point to useful information
3585	OsiBranchingInformation usefulInfo = model->usefulInformation();
3586	// and modify
3587	usefulInfo.depth_ = depth_;
3588
3589	// compute current state
3590	int numberObjectInfeasibilities; // just odd ones
3591	int numberIntegerInfeasibilities;
3592	model->feasibleSolution(
3593	numberIntegerInfeasibilities,
3594	numberObjectInfeasibilities);
3595	# ifdef COIN_HAS_CLP
3596	OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (solver);
3597	int saveClpOptions = 0;
3598	bool fastIterations = (model->specialOptions() & 8) != 0;
3599	if (osiclp && fastIterations) {
3600	// for faster hot start
3601	saveClpOptions = osiclp->specialOptions();
3602	osiclp->setSpecialOptions(saveClpOptions \| 8192);
3603	}
3604	# else
3605	bool fastIterations = false ;
3606	# endif
3607	/*
3608	Scan for branching objects that indicate infeasibility. Choose candidates
3609	using priority as the first criteria, then integer infeasibility.
3610
3611	The algorithm is to fill the array with a set of good candidates (by
3612	infeasibility) with priority bestPriority. Finding a candidate with
3613	priority better (less) than bestPriority flushes the choice array. (This
3614	serves as initialization when the first candidate is found.)
3615
3616	*/
3617	numberToDo = 0;
3618	for (i = 0; i < numberObjects; i++) {
3619	OsiObject * object = model->modifiableObject(i);
3620	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3621	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3622	if (!dynamicObject)
3623	continue;
3624	int preferredWay;
3625	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
3626	int iColumn = dynamicObject->columnNumber();
3627	if (saveUpper[iColumn] == saveLower[iColumn])
3628	continue;
3629	if (infeasibility)
3630	sort[numberToDo] = -1.0e10 - infeasibility;
3631	else
3632	sort[numberToDo] = -fabs(dj[iColumn]);
3633	whichObject[numberToDo++] = i;
3634	}
3635	// Save basis
3636	CoinWarmStart * ws = solver->getWarmStart();
3637	int saveLimit;
3638	solver->getIntParam(OsiMaxNumIterationHotStart, saveLimit);
3639	int targetIterations = CoinMax(500, numberIterationsAllowed / numberObjects);
3640	if (saveLimit < targetIterations)
3641	solver->setIntParam(OsiMaxNumIterationHotStart, targetIterations);
3642	// Mark hot start
3643	solver->markHotStart();
3644	// Sort
3645	CoinSort_2(sort, sort + numberToDo, whichObject);
3646	double * currentSolution = model->currentSolution();
3647	double objMin = 1.0e50;
3648	double objMax = -1.0e50;
3649	bool needResolve = false;
3650	/*
3651	Now calculate the cost forcing the variable up and down.
3652	*/
3653	int iDo;
3654	for (iDo = 0; iDo < numberToDo; iDo++) {
3655	CbcStrongInfo choice;
3656	int iObject = whichObject[iDo];
3657	OsiObject * object = model->modifiableObject(iObject);
3658	CbcSimpleIntegerDynamicPseudoCost * dynamicObject =
3659	dynamic_cast <CbcSimpleIntegerDynamicPseudoCost *>(object) ;
3660	if (!dynamicObject)
3661	continue;
3662	int iColumn = dynamicObject->columnNumber();
3663	int preferredWay;
3664	/*
3665	Update the information held in the object.
3666	*/
3667	object->infeasibility(&usefulInfo, preferredWay);
3668	double value = currentSolution[iColumn];
3669	double nearest = floor(value + 0.5);
3670	double lowerValue = floor(value);
3671	bool satisfied = false;
3672	if (fabs(value - nearest) <= integerTolerance \|\| value < saveLower[iColumn] \|\| value > saveUpper[iColumn]) {
3673	satisfied = true;
3674	double newValue;
3675	if (nearest < saveUpper[iColumn]) {
3676	newValue = nearest + 1.0001 * integerTolerance;
3677	lowerValue = nearest;
3678	} else {
3679	newValue = nearest - 1.0001 * integerTolerance;
3680	lowerValue = nearest - 1;
3681	}
3682	currentSolution[iColumn] = newValue;
3683	}
3684	double upperValue = lowerValue + 1.0;
3685	//CbcSimpleInteger * obj =
3686	//dynamic_cast <CbcSimpleInteger *>(object) ;
3687	//if (obj) {
3688	//choice.possibleBranch=obj->createCbcBranch(solver,&usefulInfo,preferredWay);
3689	//} else {
3690	CbcObject * obj =
3691	dynamic_cast <CbcObject *>(object) ;
3692	assert (obj);
3693	choice.possibleBranch = obj->createCbcBranch(solver, &usefulInfo, preferredWay);
3694	//}
3695	currentSolution[iColumn] = value;
3696	// Save which object it was
3697	choice.objectNumber = iObject;
3698	choice.numIntInfeasUp = numberUnsatisfied_;
3699	choice.numIntInfeasDown = numberUnsatisfied_;
3700	choice.downMovement = 0.0;
3701	choice.upMovement = 0.0;
3702	choice.numItersDown = 0;
3703	choice.numItersUp = 0;
3704	choice.fix = 0; // say not fixed
3705	double objectiveChange ;
3706	double newObjectiveValue = 1.0e100;
3707	int j;
3708	// status is 0 finished, 1 infeasible and other
3709	int iStatus;
3710	/*
3711	Try the down direction first. (Specify the initial branching alternative as
3712	down with a call to way(-1). Each subsequent call to branch() performs the
3713	specified branch and advances the branch object state to the next branch
3714	alternative.)
3715	*/
3716	choice.possibleBranch->way(-1) ;
3717	choice.possibleBranch->branch() ;
3718	if (fabs(value - lowerValue) > integerTolerance) {
3719	solver->solveFromHotStart() ;
3720	/*
3721	We now have an estimate of objective degradation that we can use for strong
3722	branching. If we're over the cutoff, the variable is monotone up.
3723	If we actually made it to optimality, check for a solution, and if we have
3724	a good one, call setBestSolution to process it. Note that this may reduce the
3725	cutoff, so we check again to see if we can declare this variable monotone.
3726	*/
3727	if (solver->isProvenOptimal())
3728	iStatus = 0; // optimal
3729	else if (solver->isIterationLimitReached()
3730	&& !solver->isDualObjectiveLimitReached())
3731	iStatus = 2; // unknown
3732	else
3733	iStatus = 1; // infeasible
3734	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3735	choice.numItersDown = solver->getIterationCount();
3736	numberIterationsAllowed -= choice.numItersDown;
3737	objectiveChange = newObjectiveValue - objectiveValue_;
3738	if (!iStatus) {
3739	choice.finishedDown = true ;
3740	if (newObjectiveValue >= cutoff) {
3741	objectiveChange = 1.0e100; // say infeasible
3742	} else {
3743	// See if integer solution
3744	if (model->feasibleSolution(choice.numIntInfeasDown,
3745	choice.numObjInfeasDown)
3746	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
3747	model->setBestSolution(CBC_STRONGSOL,
3748	newObjectiveValue,
3749	solver->getColSolution()) ;
3750	model->setLastHeuristic(NULL);
3751	model->incrementUsed(solver->getColSolution());
3752	cutoff = model->getCutoff();
3753	if (newObjectiveValue >= cutoff) // new cutoff
3754	objectiveChange = 1.0e100 ;
3755	}
3756	}
3757	} else if (iStatus == 1) {
3758	objectiveChange = 1.0e100 ;
3759	} else {
3760	// Can't say much as we did not finish
3761	choice.finishedDown = false ;
3762	}
3763	choice.downMovement = objectiveChange ;
3764	}
3765	// restore bounds
3766	for ( j = 0; j < numberColumns; j++) {
3767	if (saveLower[j] != lower[j])
3768	solver->setColLower(j, saveLower[j]);
3769	if (saveUpper[j] != upper[j])
3770	solver->setColUpper(j, saveUpper[j]);
3771	}
3772	// repeat the whole exercise, forcing the variable up
3773	choice.possibleBranch->branch();
3774	if (fabs(value - upperValue) > integerTolerance) {
3775	solver->solveFromHotStart() ;
3776	/*
3777	We now have an estimate of objective degradation that we can use for strong
3778	branching. If we're over the cutoff, the variable is monotone up.
3779	If we actually made it to optimality, check for a solution, and if we have
3780	a good one, call setBestSolution to process it. Note that this may reduce the
3781	cutoff, so we check again to see if we can declare this variable monotone.
3782	*/
3783	if (solver->isProvenOptimal())
3784	iStatus = 0; // optimal
3785	else if (solver->isIterationLimitReached()
3786	&& !solver->isDualObjectiveLimitReached())
3787	iStatus = 2; // unknown
3788	else
3789	iStatus = 1; // infeasible
3790	newObjectiveValue = solver->getObjSense() * solver->getObjValue();
3791	choice.numItersUp = solver->getIterationCount();
3792	numberIterationsAllowed -= choice.numItersUp;
3793	objectiveChange = newObjectiveValue - objectiveValue_;
3794	if (!iStatus) {
3795	choice.finishedUp = true ;
3796	if (newObjectiveValue >= cutoff) {
3797	objectiveChange = 1.0e100; // say infeasible
3798	} else {
3799	// See if integer solution
3800	if (model->feasibleSolution(choice.numIntInfeasUp,
3801	choice.numObjInfeasUp)
3802	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
3803	model->setBestSolution(CBC_STRONGSOL,
3804	newObjectiveValue,
3805	solver->getColSolution()) ;
3806	model->setLastHeuristic(NULL);
3807	model->incrementUsed(solver->getColSolution());
3808	cutoff = model->getCutoff();
3809	if (newObjectiveValue >= cutoff) // new cutoff
3810	objectiveChange = 1.0e100 ;
3811	}
3812	}
3813	} else if (iStatus == 1) {
3814	objectiveChange = 1.0e100 ;
3815	} else {
3816	// Can't say much as we did not finish
3817	choice.finishedUp = false ;
3818	}
3819	choice.upMovement = objectiveChange ;
3820
3821	// restore bounds
3822	for ( j = 0; j < numberColumns; j++) {
3823	if (saveLower[j] != lower[j])
3824	solver->setColLower(j, saveLower[j]);
3825	if (saveUpper[j] != upper[j])
3826	solver->setColUpper(j, saveUpper[j]);
3827	}
3828	}
3829	// If objective goes above certain amount we can set bound
3830	int jInt = back[iColumn];
3831	newLower[jInt] = upperValue;
3832	if (choice.finishedDown)
3833	objLower[jInt] = choice.downMovement + objectiveValue_;
3834	else
3835	objLower[jInt] = objectiveValue_;
3836	newUpper[jInt] = lowerValue;
3837	if (choice.finishedUp)
3838	objUpper[jInt] = choice.upMovement + objectiveValue_;
3839	else
3840	objUpper[jInt] = objectiveValue_;
3841	objMin = CoinMin(CoinMin(objLower[jInt], objUpper[jInt]), objMin);
3842	/*
3843	End of evaluation for this candidate variable. Possibilities are:
3844	* Both sides below cutoff; this variable is a candidate for branching.
3845	* Both sides infeasible or above the objective cutoff: no further action
3846	here. Break from the evaluation loop and assume the node will be purged
3847	by the caller.
3848	* One side below cutoff: Install the branch (i.e., fix the variable). Break
3849	from the evaluation loop and assume the node will be reoptimised by the
3850	caller.
3851	*/
3852	if (choice.upMovement < 1.0e100) {
3853	if (choice.downMovement < 1.0e100) {
3854	objMax = CoinMax(CoinMax(objLower[jInt], objUpper[jInt]), objMax);
3855	// In case solution coming in was odd
3856	choice.upMovement = CoinMax(0.0, choice.upMovement);
3857	choice.downMovement = CoinMax(0.0, choice.downMovement);
3858	// feasible -
3859	model->messageHandler()->message(CBC_STRONG, *model->messagesPointer())
3860	<< iObject << iColumn
3861	<< choice.downMovement << choice.numIntInfeasDown
3862	<< choice.upMovement << choice.numIntInfeasUp
3863	<< value
3864	<< CoinMessageEol;
3865	} else {
3866	// up feasible, down infeasible
3867	anyAction = -1;
3868	if (!satisfied)
3869	needResolve = true;
3870	choice.fix = 1;
3871	numberToFix++;
3872	saveLower[iColumn] = upperValue;
3873	solver->setColLower(iColumn, upperValue);
3874	}
3875	} else {
3876	if (choice.downMovement < 1.0e100) {
3877	// down feasible, up infeasible
3878	anyAction = -1;
3879	if (!satisfied)
3880	needResolve = true;
3881	choice.fix = -1;
3882	numberToFix++;
3883	saveUpper[iColumn] = lowerValue;
3884	solver->setColUpper(iColumn, lowerValue);
3885	} else {
3886	// neither side feasible
3887	anyAction = -2;
3888	COIN_DETAIL_PRINT(printf("Both infeasible for choice %d sequence %d\n", i,
3889	model->object(choice.objectNumber)->columnNumber()));
3890	delete ws;
3891	ws = NULL;
3892	//solver->writeMps("bad");
3893	numberToFix = -1;
3894	delete choice.possibleBranch;
3895	choice.possibleBranch = NULL;
3896	break;
3897	}
3898	}
3899	delete choice.possibleBranch;
3900	if (numberIterationsAllowed <= 0)
3901	break;
3902	//printf("obj %d, col %d, down %g up %g value %g\n",iObject,iColumn,
3903	// choice.downMovement,choice.upMovement,value);
3904	}
3905	COIN_DETAIL_PRINT(printf("Best possible solution %g, can fix more if solution of %g found - looked at %d variables in %d iterations\n",
3906	objMin, objMax, iDo, model->numberAnalyzeIterations() - numberIterationsAllowed));
3907	model->setNumberAnalyzeIterations(numberIterationsAllowed);
3908	// Delete the snapshot
3909	solver->unmarkHotStart();
3910	// back to normal
3911	solver->setHintParam(OsiDoInBranchAndCut, true, OsiHintDo, NULL) ;
3912	solver->setIntParam(OsiMaxNumIterationHotStart, saveLimit);
3913	// restore basis
3914	solver->setWarmStart(ws);
3915	delete ws;
3916
3917	delete [] sort;
3918	delete [] whichObject;
3919	delete [] saveLower;
3920	delete [] saveUpper;
3921	delete [] back;
3922	// restore solution
3923	solver->setColSolution(saveSolution);
3924	# ifdef COIN_HAS_CLP
3925	if (osiclp)
3926	osiclp->setSpecialOptions(saveClpOptions);
3927	# endif
3928	model->reserveCurrentSolution(saveSolution);
3929	delete [] saveSolution;
3930	if (needResolve)
3931	solver->resolve();
3932	return numberToFix;
3933	}
3934
3935
3936	CbcNode::CbcNode(const CbcNode & rhs)
3937	: CoinTreeNode(rhs)
3938	{
3939	#ifdef CHECK_NODE
3940	printf("CbcNode %x Constructor from rhs %x\n", this, &rhs);
3941	#endif
3942	if (rhs.nodeInfo_)
3943	nodeInfo_ = rhs.nodeInfo_->clone();
3944	else
3945	nodeInfo_ = NULL;
3946	objectiveValue_ = rhs.objectiveValue_;
3947	guessedObjectiveValue_ = rhs.guessedObjectiveValue_;
3948	sumInfeasibilities_ = rhs.sumInfeasibilities_;
3949	if (rhs.branch_)
3950	branch_ = rhs.branch_->clone();
3951	else
3952	branch_ = NULL;
3953	depth_ = rhs.depth_;
3954	numberUnsatisfied_ = rhs.numberUnsatisfied_;
3955	nodeNumber_ = rhs.nodeNumber_;
3956	state_ = rhs.state_;
3957	if (nodeInfo_)
3958	assert ((state_&2) != 0);
3959	else
3960	assert ((state_&2) == 0);
3961	}
3962
3963	CbcNode &
3964	CbcNode::operator=(const CbcNode & rhs)
3965	{
3966	if (this != &rhs) {
3967	delete nodeInfo_;
3968	if (rhs.nodeInfo_)
3969	nodeInfo_ = rhs.nodeInfo_->clone();
3970	else
3971	nodeInfo_ = NULL;
3972	objectiveValue_ = rhs.objectiveValue_;
3973	guessedObjectiveValue_ = rhs.guessedObjectiveValue_;
3974	sumInfeasibilities_ = rhs.sumInfeasibilities_;
3975	if (rhs.branch_)
3976	branch_ = rhs.branch_->clone();
3977	else
3978	branch_ = NULL,
3979	depth_ = rhs.depth_;
3980	numberUnsatisfied_ = rhs.numberUnsatisfied_;
3981	nodeNumber_ = rhs.nodeNumber_;
3982	state_ = rhs.state_;
3983	if (nodeInfo_)
3984	assert ((state_&2) != 0);
3985	else
3986	assert ((state_&2) == 0);
3987	}
3988	return *this;
3989	}
3990	CbcNode::~CbcNode ()
3991	{
3992	#ifdef CHECK_NODE
3993	if (nodeInfo_) {
3994	printf("CbcNode %x Destructor nodeInfo %x (%d)\n",
3995	this, nodeInfo_, nodeInfo_->numberPointingToThis());
3996	//assert(nodeInfo_->numberPointingToThis()>=0);
3997	} else {
3998	printf("CbcNode %x Destructor nodeInfo %x (?)\n",
3999	this, nodeInfo_);
4000	}
4001	#endif
4002	if (nodeInfo_) {
4003	// was if (nodeInfo_&&(state_&2)!=0) {
4004	nodeInfo_->nullOwner();
4005	int numberToDelete = nodeInfo_->numberBranchesLeft();
4006	// CbcNodeInfo * parent = nodeInfo_->parent();
4007	//assert (nodeInfo_->numberPointingToThis()>0);
4008	if (nodeInfo_->decrement(numberToDelete) == 0 \|\| (state_&2) == 0) {
4009	if ((state_&2) == 0)
4010	nodeInfo_->nullParent();
4011	delete nodeInfo_;
4012	} else {
4013	//printf("node %x nodeinfo %x parent %x\n",this,nodeInfo_,nodeInfo_->parent());
4014	// anyway decrement parent
4015	//if (parent)
4016	///parent->decrement(1);
4017	}
4018	}
4019	delete branch_;
4020	}
4021	// Decrement active cut counts
4022	void
4023	CbcNode::decrementCuts(int change)
4024	{
4025	if (nodeInfo_)
4026	assert ((state_&2) != 0);
4027	else
4028	assert ((state_&2) == 0);
4029	if (nodeInfo_) {
4030	nodeInfo_->decrementCuts(change);
4031	}
4032	}
4033	void
4034	CbcNode::decrementParentCuts(CbcModel * model, int change)
4035	{
4036	if (nodeInfo_)
4037	assert ((state_&2) != 0);
4038	else
4039	assert ((state_&2) == 0);
4040	if (nodeInfo_) {
4041	nodeInfo_->decrementParentCuts(model, change);
4042	}
4043	}
4044
4045	/*
4046	Initialize reference counts (numberPointingToThis, numberBranchesLeft_)
4047	in the attached nodeInfo_.
4048	*/
4049	void
4050	CbcNode::initializeInfo()
4051	{
4052	assert(nodeInfo_ && branch_) ;
4053	nodeInfo_->initializeInfo(branch_->numberBranches());
4054	assert ((state_&2) != 0);
4055	assert (nodeInfo_->numberBranchesLeft() ==
4056	branch_->numberBranchesLeft());
4057	}
4058	// Nulls out node info
4059	void
4060	CbcNode::nullNodeInfo()
4061	{
4062	nodeInfo_ = NULL;
4063	// say not active
4064	state_ &= ~2;
4065	}
4066
4067	int
4068	CbcNode::branch(OsiSolverInterface * solver)
4069	{
4070	double changeInGuessed;
4071	assert (nodeInfo_->numberBranchesLeft() ==
4072	branch_->numberBranchesLeft());
4073	if (!solver)
4074	changeInGuessed = branch_->branch();
4075	else
4076	changeInGuessed = branch_->branch(solver);
4077	guessedObjectiveValue_ += changeInGuessed;
4078	//#define PRINTIT
4079	#ifdef PRINTIT
4080	int numberLeft = nodeInfo_->numberBranchesLeft();
4081	CbcNodeInfo * parent = nodeInfo_->parent();
4082	int parentNodeNumber = -1;
4083	CbcBranchingObject * object1 =
4084	dynamic_cast<CbcBranchingObject *>(branch_) ;
4085	//OsiObject * object = object1->
4086	//int sequence = object->columnNumber);
4087	int id = -1;
4088	double value = 0.0;
4089	if (object1) {
4090	id = object1->variable();
4091	value = object1->value();
4092	}
4093	printf("id %d value %g objvalue %g\n", id, value, objectiveValue_);
4094	if (parent)
4095	parentNodeNumber = parent->nodeNumber();
4096	printf("Node number %d, %s, way %d, depth %d, parent node number %d\n",
4097	nodeInfo_->nodeNumber(), (numberLeft == 2) ? "leftBranch" : "rightBranch",
4098	way(), depth_, parentNodeNumber);
4099	assert (parentNodeNumber != nodeInfo_->nodeNumber());
4100	#endif
4101	return nodeInfo_->branchedOn();
4102	}
4103	/* Active arm of the attached OsiBranchingObject.
4104
4105	In the simplest instance, coded -1 for the down arm of the branch, +1 for
4106	the up arm. But see OsiBranchingObject::way()
4107	Use nodeInfo--.numberBranchesLeft_ to see how active
4108
4109	Except that there is no OsiBranchingObject::way(), and this'll fail in any
4110	event because we have various OsiXXXBranchingObjects which aren't descended
4111	from CbcBranchingObjects. I think branchIndex() is the appropriate
4112	equivalent, but could be wrong. (lh, 061220)
4113
4114	071212: I'm finally getting back to cbc-generic and rescuing a lot of my
4115	annotation from branches/devel (which was killed in summer). I'm going to
4116	put back an assert(obj) just to see what happens. It's still present as of
4117	the most recent change to CbcNode (r833).
4118
4119	080104: Yep, we can arrive here with an OsiBranchingObject. Removed the
4120	assert, it's served its purpose.
4121
4122	080226: John finally noticed this problem and added a way() method to the
4123	OsiBranchingObject hierarchy. Removing my workaround.
4124
4125	*/
4126	int
4127	CbcNode::way() const
4128	{
4129	if (branch_) {
4130	CbcBranchingObject * obj =
4131	dynamic_cast <CbcBranchingObject *>(branch_) ;
4132	if (obj) {
4133	return obj->way();
4134	} else {
4135	OsiTwoWayBranchingObject * obj2 =
4136	dynamic_cast <OsiTwoWayBranchingObject *>(branch_) ;
4137	assert (obj2);
4138	return obj2->way();
4139	}
4140	} else {
4141	return 0;
4142	}
4143	}
4144	/* Create a branching object for the node
4145
4146	The routine scans the object list of the model and selects a set of
4147	unsatisfied objects as candidates for branching. The candidates are
4148	evaluated, and an appropriate branch object is installed.
4149
4150	The numberPassesLeft is decremented to stop fixing one variable each time
4151	and going on and on (e.g. for stock cutting, air crew scheduling)
4152
4153	If evaluation determines that an object is monotone or infeasible,
4154	the routine returns immediately. In the case of a monotone object,
4155	the branch object has already been called to modify the model.
4156
4157	Return value:
4158	<ul>
4159	<li> 0: A branching object has been installed
4160	<li> -1: A monotone object was discovered
4161	<li> -2: An infeasible object was discovered
4162	</ul>
4163	Branch state:
4164	<ul>
4165	<li> -1: start
4166	<li> -1: A monotone object was discovered
4167	<li> -2: An infeasible object was discovered
4168	</ul>
4169	*/
4170	int
4171	CbcNode::chooseOsiBranch (CbcModel * model,
4172	CbcNode * lastNode,
4173	OsiBranchingInformation * usefulInfo,
4174	int branchState)
4175	{
4176	int returnStatus = 0;
4177	if (lastNode)
4178	depth_ = lastNode->depth_ + 1;
4179	else
4180	depth_ = 0;
4181	OsiSolverInterface * solver = model->solver();
4182	objectiveValue_ = solver->getObjValue() * solver->getObjSense();
4183	usefulInfo->objectiveValue_ = objectiveValue_;
4184	usefulInfo->depth_ = depth_;
4185	const double * saveInfoSol = usefulInfo->solution_;
4186	double * saveSolution = new double[solver->getNumCols()];
4187	memcpy(saveSolution, solver->getColSolution(), solver->getNumCols()*sizeof(double));
4188	usefulInfo->solution_ = saveSolution;
4189	OsiChooseVariable * choose = model->branchingMethod()->chooseMethod();
4190	int numberUnsatisfied = -1;
4191	if (branchState < 0) {
4192	// initialize
4193	// initialize sum of "infeasibilities"
4194	sumInfeasibilities_ = 0.0;
4195	numberUnsatisfied = choose->setupList(usefulInfo, true);
4196	numberUnsatisfied_ = numberUnsatisfied;
4197	branchState = 0;
4198	if (numberUnsatisfied_ < 0) {
4199	// infeasible
4200	delete [] saveSolution;
4201	return -2;
4202	}
4203	}
4204	// unset best
4205	int best = -1;
4206	choose->setBestObjectIndex(-1);
4207	if (numberUnsatisfied) {
4208	if (branchState > 0 \|\| !choose->numberOnList()) {
4209	// we need to return at once - don't do strong branching or anything
4210	if (choose->numberOnList() \|\| !choose->numberStrong()) {
4211	best = choose->candidates()[0];
4212	choose->setBestObjectIndex(best);
4213	} else {
4214	// nothing on list - need to try again - keep any solution
4215	numberUnsatisfied = choose->setupList(usefulInfo, false);
4216	numberUnsatisfied_ = numberUnsatisfied;
4217	if (numberUnsatisfied) {
4218	best = choose->candidates()[0];
4219	choose->setBestObjectIndex(best);
4220	}
4221	}
4222	} else {
4223	// carry on with strong branching or whatever
4224	int returnCode = choose->chooseVariable(solver, usefulInfo, true);
4225	// update number of strong iterations etc
4226	model->incrementStrongInfo(choose->numberStrongDone(), choose->numberStrongIterations(),
4227	returnCode == -1 ? 0 : choose->numberStrongFixed(), returnCode == -1);
4228	if (returnCode > 1) {
4229	// has fixed some
4230	returnStatus = -1;
4231	} else if (returnCode == -1) {
4232	// infeasible
4233	returnStatus = -2;
4234	} else if (returnCode == 0) {
4235	// normal
4236	returnStatus = 0;
4237	numberUnsatisfied = 1;
4238	} else {
4239	// ones on list satisfied - double check
4240	numberUnsatisfied = choose->setupList(usefulInfo, false);
4241	numberUnsatisfied_ = numberUnsatisfied;
4242	if (numberUnsatisfied) {
4243	best = choose->candidates()[0];
4244	choose->setBestObjectIndex(best);
4245	}
4246	}
4247	}
4248	}
4249	delete branch_;
4250	branch_ = NULL;
4251	guessedObjectiveValue_ = COIN_DBL_MAX;//objectiveValue_; // for now
4252	if (!returnStatus) {
4253	if (numberUnsatisfied) {
4254	// create branching object
4255	const OsiObject * obj = model->solver()->object(choose->bestObjectIndex());
4256	//const OsiSolverInterface * solver = usefulInfo->solver_;
4257	branch_ = obj->createBranch(model->solver(), usefulInfo, obj->whichWay());
4258	}
4259	}
4260	usefulInfo->solution_ = saveInfoSol;
4261	delete [] saveSolution;
4262	// may have got solution
4263	if (choose->goodSolution()
4264	&& model->problemFeasibility()->feasible(model, -1) >= 0) {
4265	// yes
4266	double objValue = choose->goodObjectiveValue();
4267	model->setBestSolution(CBC_STRONGSOL,
4268	objValue,
4269	choose->goodSolution()) ;
4270	model->setLastHeuristic(NULL);
4271	model->incrementUsed(choose->goodSolution());
4272	choose->clearGoodSolution();
4273	}
4274	return returnStatus;
4275	}
4276	int
4277	CbcNode::chooseClpBranch (CbcModel * model,
4278	CbcNode * lastNode)
4279	{
4280	assert(lastNode);
4281	depth_ = lastNode->depth_ + 1;
4282	delete branch_;
4283	branch_ = NULL;
4284	OsiSolverInterface * solver = model->solver();
4285	//double saveObjectiveValue = solver->getObjValue();
4286	//double objectiveValue = CoinMax(solver->getObjSense()*saveObjectiveValue,objectiveValue_);
4287	const double * lower = solver->getColLower();
4288	const double * upper = solver->getColUpper();
4289	// point to useful information
4290	OsiBranchingInformation usefulInfo = model->usefulInformation();
4291	// and modify
4292	usefulInfo.depth_ = depth_;
4293	int i;
4294	//bool beforeSolution = model->getSolutionCount()==0;
4295	int numberObjects = model->numberObjects();
4296	int numberColumns = model->getNumCols();
4297	double * saveUpper = new double[numberColumns];
4298	double * saveLower = new double[numberColumns];
4299
4300	// Save solution in case heuristics need good solution later
4301
4302	double * saveSolution = new double[numberColumns];
4303	memcpy(saveSolution, solver->getColSolution(), numberColumns*sizeof(double));
4304	model->reserveCurrentSolution(saveSolution);
4305	for (i = 0; i < numberColumns; i++) {
4306	saveLower[i] = lower[i];
4307	saveUpper[i] = upper[i];
4308	}
4309	// Save basis
4310	CoinWarmStart * ws = solver->getWarmStart();
4311	numberUnsatisfied_ = 0;
4312	// initialize sum of "infeasibilities"
4313	sumInfeasibilities_ = 0.0;
4314	// Note looks as if off end (hidden one)
4315	OsiObject * object = model->modifiableObject(numberObjects);
4316	CbcGeneralDepth * thisOne = dynamic_cast <CbcGeneralDepth *> (object);
4317	assert (thisOne);
4318	OsiClpSolverInterface * clpSolver
4319	= dynamic_cast<OsiClpSolverInterface *> (solver);
4320	assert (clpSolver);
4321	ClpSimplex * simplex = clpSolver->getModelPtr();
4322	int preferredWay;
4323	double infeasibility = object->infeasibility(&usefulInfo, preferredWay);
4324	if (thisOne->whichSolution() >= 0) {
4325	ClpNode * nodeInfo = thisOne->nodeInfo(thisOne->whichSolution());
4326	nodeInfo->applyNode(simplex, 2);
4327	int saveLogLevel = simplex->logLevel();
4328	simplex->setLogLevel(0);
4329	simplex->dual();
4330	simplex->setLogLevel(saveLogLevel);
4331	double cutoff = model->getCutoff();
4332	bool goodSolution = true;
4333	if (simplex->status()) {
4334	//simplex->writeMps("bad7.mps",2);
4335	if (nodeInfo->objectiveValue() > cutoff - 1.0e-2)
4336	goodSolution = false;
4337	else
4338	assert (!simplex->status());
4339	}
4340	if (goodSolution) {
4341	double newObjectiveValue = solver->getObjSense() * solver->getObjValue();
4342	// See if integer solution
4343	int numInf;
4344	int numInf2;
4345	bool gotSol = model->feasibleSolution(numInf, numInf2);
4346	if (!gotSol) {
4347	COIN_DETAIL_PRINT(printf("numinf %d\n", numInf));
4348	double * sol = simplex->primalColumnSolution();
4349	for (int i = 0; i < numberColumns; i++) {
4350	if (simplex->isInteger(i)) {
4351	double value = floor(sol[i] + 0.5);
4352	if (fabs(value - sol[i]) > 1.0e-7) {
4353	COIN_DETAIL_PRINT(printf("%d value %g\n", i, sol[i]));
4354	if (fabs(value - sol[i]) < 1.0e-3) {
4355	sol[i] = value;
4356	}
4357	}
4358	}
4359	}
4360	simplex->writeMps("bad8.mps", 2);
4361	bool gotSol = model->feasibleSolution(numInf, numInf2);
4362	if (!gotSol)
4363	assert (gotSol);
4364	}
4365	model->setBestSolution(CBC_STRONGSOL,
4366	newObjectiveValue,
4367	solver->getColSolution()) ;
4368	model->setLastHeuristic(NULL);
4369	model->incrementUsed(solver->getColSolution());
4370	}
4371	}
4372	// restore bounds
4373	{
4374	for (int j = 0; j < numberColumns; j++) {
4375	if (saveLower[j] != lower[j])
4376	solver->setColLower(j, saveLower[j]);
4377	if (saveUpper[j] != upper[j])
4378	solver->setColUpper(j, saveUpper[j]);
4379	}
4380	}
4381	// restore basis
4382	solver->setWarmStart(ws);
4383	delete ws;
4384	int anyAction;
4385	//#define CHECK_PATH
4386	#ifdef CHECK_PATH
4387	extern int gotGoodNode_Z;
4388	if (gotGoodNode_Z >= 0)
4389	printf("good node %d %g\n", gotGoodNode_Z, infeasibility);
4390	#endif
4391	if (infeasibility > 0.0) {
4392	if (infeasibility == COIN_DBL_MAX) {
4393	anyAction = -2; // infeasible
4394	} else {
4395	branch_ = thisOne->createCbcBranch(solver, &usefulInfo, preferredWay);
4396	// Set to firat one (and change when re-pushing)
4397	CbcGeneralBranchingObject * branch =
4398	dynamic_cast <CbcGeneralBranchingObject *> (branch_);
4399	branch->state(objectiveValue_, sumInfeasibilities_,
4400	numberUnsatisfied_, 0);
4401	branch->setNode(this);
4402	anyAction = 0;
4403	}
4404	} else {
4405	anyAction = -1;
4406	}
4407	#ifdef CHECK_PATH
4408	gotGoodNode_Z = -1;
4409	#endif
4410	// Set guessed solution value
4411	guessedObjectiveValue_ = objectiveValue_ + 1.0e-5;
4412	delete [] saveLower;
4413	delete [] saveUpper;
4414
4415	// restore solution
4416	solver->setColSolution(saveSolution);
4417	delete [] saveSolution;
4418	return anyAction;
4419	}
4420	/* Double checks in case node can change its mind!
4421	Returns objective value
4422	Can change objective etc */
4423	double
4424	CbcNode::checkIsCutoff(double cutoff)
4425	{
4426	branch_->checkIsCutoff(cutoff);
4427	return objectiveValue_;
4428	}
4429

Note: See TracBrowser for help on using the repository browser.

Download in other formats: