Context Navigation

OsiClpSolverInterface.cpp @ 2395

Last change on this file since 2395 was 2386, checked in by unxusr, 2 years ago
formatting
File size: 368.5 KB

Line
1	// $Id$
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	#include <cassert>
7	#ifdef CBC_STATISTICS
8	extern int osi_crunch;
9	extern int osi_primal;
10	extern int osi_dual;
11	extern int osi_hot;
12	#endif
13	#include "CoinTime.hpp"
14	#include "CoinHelperFunctions.hpp"
15	#include "CoinIndexedVector.hpp"
16	#include "CoinModel.hpp"
17	#include "CoinMpsIO.hpp"
18	#include "CoinSort.hpp"
19	#include "ClpDualRowSteepest.hpp"
20	#include "ClpPrimalColumnSteepest.hpp"
21	#include "ClpPackedMatrix.hpp"
22	#include "ClpDualRowDantzig.hpp"
23	#include "ClpPrimalColumnDantzig.hpp"
24	#include "ClpFactorization.hpp"
25	#include "ClpObjective.hpp"
26	#include "ClpSimplex.hpp"
27	#include "ClpSimplexOther.hpp"
28	#include "ClpSimplexPrimal.hpp"
29	#include "ClpSimplexDual.hpp"
30	#include "ClpNonLinearCost.hpp"
31	#include "OsiClpSolverInterface.hpp"
32	#include "OsiBranchingObject.hpp"
33	#include "OsiCuts.hpp"
34	#include "OsiRowCut.hpp"
35	#include "OsiColCut.hpp"
36	#include "ClpPresolve.hpp"
37	#include "CoinLpIO.hpp"
38	//#define PRINT_TIME
39	#ifdef PRINT_TIME
40	static double totalTime = 0.0;
41	#endif
42	//#define SAVE_MODEL 1
43	#ifdef SAVE_MODEL
44	static int resolveTry = 0;
45	static int loResolveTry = 0;
46	static int hiResolveTry = 9999999;
47	#endif
48	//#############################################################################
49	// Solve methods
50	//#############################################################################
51	void OsiClpSolverInterface::initialSolve()
52	{
53	#define KEEP_SMALL
54	#ifdef KEEP_SMALL
55	if (smallModel_) {
56	delete[] spareArrays_;
57	spareArrays_ = NULL;
58	delete smallModel_;
59	smallModel_ = NULL;
60	}
61	#endif
62	if ((specialOptions_ & 2097152) != 0 \|\| (specialOptions_ & 4194304) != 0) {
63	bool takeHint;
64	OsiHintStrength strength;
65	int algorithm = 0;
66	getHintParam(OsiDoDualInInitial, takeHint, strength);
67	if (strength != OsiHintIgnore)
68	algorithm = takeHint ? -1 : 1;
69	if (algorithm > 0 \|\| (specialOptions_ & 4194304) != 0) {
70	// Gub
71	resolveGub((9 * modelPtr_->numberRows()) / 10);
72	return;
73	}
74	}
75	bool deleteSolver;
76	ClpSimplex *solver;
77	#ifdef PRINT_TIME
78	double time1 = CoinCpuTime();
79	#endif
80	int userFactorizationFrequency = modelPtr_->factorization()->maximumPivots();
81	int totalIterations = 0;
82	bool abortSearch = false;
83	ClpObjective *savedObjective = NULL;
84	double savedDualLimit = modelPtr_->dblParam_[ClpDualObjectiveLimit];
85	if (fakeObjective_) {
86	// Clear (no objective, 0-1 and in B&B)
87	modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions() & (~128));
88	// See if all with costs fixed
89	int numberColumns = modelPtr_->numberColumns_;
90	const double *obj = modelPtr_->objective();
91	const double *lower = modelPtr_->columnLower();
92	const double *upper = modelPtr_->columnUpper();
93	int i;
94	for (i = 0; i < numberColumns; i++) {
95	double objValue = obj[i];
96	if (objValue) {
97	if (lower[i] != upper[i])
98	break;
99	}
100	}
101	if (i == numberColumns) {
102	// Check (Clp fast dual)
103	if ((specialOptions_ & 524288) == 0) {
104	// Set fake
105	savedObjective = modelPtr_->objective_;
106	modelPtr_->objective_ = fakeObjective_;
107	modelPtr_->dblParam_[ClpDualObjectiveLimit] = COIN_DBL_MAX;
108	} else {
109	// Set (no objective, 0-1 and in B&B)
110	modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions() \| 128);
111	}
112	}
113	}
114	// Check (in branch and bound)
115	if ((specialOptions_ & 1024) == 0) {
116	solver = new ClpSimplex(true);
117	deleteSolver = true;
118	solver->borrowModel(*modelPtr_);
119	// See if user set factorization frequency
120	// borrowModel does not move
121	solver->factorization()->maximumPivots(userFactorizationFrequency);
122	} else {
123	solver = modelPtr_;
124	deleteSolver = false;
125	}
126	// Treat as if user simplex not enabled
127	int saveSolveType = solver->solveType();
128	bool doingPrimal = solver->algorithm() > 0;
129	if (saveSolveType == 2) {
130	disableSimplexInterface();
131	solver->setSolveType(1);
132	}
133	int saveOptions = solver->specialOptions();
134	solver->setSpecialOptions(saveOptions \| 64 \| 32768); // go as far as possible
135	// get original log levels
136	int saveMessageLevel = modelPtr_->logLevel();
137	int messageLevel = messageHandler()->logLevel();
138	int saveMessageLevel2 = messageLevel;
139	// Set message handler
140	if (!defaultHandler_)
141	solver->passInMessageHandler(handler_);
142	// But keep log level
143	solver->messageHandler()->setLogLevel(saveMessageLevel);
144	// set reasonable defaults
145	bool takeHint;
146	OsiHintStrength strength;
147	// Switch off printing if asked to
148	bool gotHint = (getHintParam(OsiDoReducePrint, takeHint, strength));
149	assert(gotHint);
150	if (strength != OsiHintIgnore && takeHint) {
151	if (messageLevel > 0)
152	messageLevel--;
153	}
154	if (messageLevel < saveMessageLevel)
155	solver->messageHandler()->setLogLevel(messageLevel);
156	// Allow for specialOptions_==1+8 forcing saving factorization
157	int startFinishOptions = 0;
158	if ((specialOptions_ & 9) == (1 + 8)) {
159	startFinishOptions = 1 + 2 + 4; // allow re-use of factorization
160	}
161	bool defaultHints = true;
162	{
163	int hint;
164	for (hint = OsiDoPresolveInInitial; hint < OsiLastHintParam; hint++) {
165	if (hint != OsiDoReducePrint && hint != OsiDoInBranchAndCut) {
166	bool yesNo;
167	OsiHintStrength strength;
168	getHintParam(static_cast< OsiHintParam >(hint), yesNo, strength);
169	if (yesNo) {
170	defaultHints = false;
171	break;
172	}
173	if (strength != OsiHintIgnore) {
174	defaultHints = false;
175	break;
176	}
177	}
178	}
179	}
180	ClpPresolve *pinfo = NULL;
181	/*
182	If basis then do primal (as user could do dual with resolve)
183	If not then see if dual feasible (and allow for gubs etc?)
184	*/
185	bool doPrimal = (basis_.numberBasicStructurals() > 0);
186	setBasis(basis_, solver);
187	bool inCbcOrOther = (modelPtr_->specialOptions() & 0x03000000) != 0;
188	if ((!defaultHints \|\| doPrimal) && !solveOptions_.getSpecialOption(6)) {
189	// scaling
190	// save initial state
191	const double *rowScale1 = solver->rowScale();
192	if (modelPtr_->solveType() == 1) {
193	gotHint = (getHintParam(OsiDoScale, takeHint, strength));
194	assert(gotHint);
195	if (strength == OsiHintIgnore \|\| takeHint) {
196	if (!solver->scalingFlag())
197	solver->scaling(3);
198	} else {
199	solver->scaling(0);
200	}
201	} else {
202	solver->scaling(0);
203	}
204	//solver->setDualBound(1.0e6);
205	//solver->setDualTolerance(1.0e-7);
206
207	//ClpDualRowSteepest steep;
208	//solver->setDualRowPivotAlgorithm(steep);
209	//solver->setPrimalTolerance(1.0e-8);
210	//ClpPrimalColumnSteepest steepP;
211	//solver->setPrimalColumnPivotAlgorithm(steepP);
212
213	// sort out hints;
214	// algorithm 0 whatever, -1 force dual, +1 force primal
215	int algorithm = 0;
216	gotHint = (getHintParam(OsiDoDualInInitial, takeHint, strength));
217	assert(gotHint);
218	if (strength != OsiHintIgnore)
219	algorithm = takeHint ? -1 : 1;
220	// crash 0 do lightweight if all slack, 1 do, -1 don't
221	int doCrash = 0;
222	gotHint = (getHintParam(OsiDoCrash, takeHint, strength));
223	assert(gotHint);
224	if (strength != OsiHintIgnore)
225	doCrash = takeHint ? 1 : -1;
226	// doPrimal set true if any structurals in basis so switch off crash
227	if (doPrimal)
228	doCrash = -1;
229
230	// presolve
231	gotHint = (getHintParam(OsiDoPresolveInInitial, takeHint, strength));
232	assert(gotHint);
233	if (strength != OsiHintIgnore && takeHint) {
234	pinfo = new ClpPresolve();
235	ClpSimplex model2 = pinfo->presolvedModel(solver, 1.0e-8);
236	if (!model2) {
237	// problem found to be infeasible - whats best?
238	model2 = solver;
239	delete pinfo;
240	pinfo = NULL;
241	} else {
242	model2->setSpecialOptions(solver->specialOptions());
243	}
244
245	// change from 200 (unless changed)
246	if (modelPtr_->factorization()->maximumPivots() == 200)
247	model2->factorization()->maximumPivots(100 + model2->numberRows() / 50);
248	else
249	model2->factorization()->maximumPivots(userFactorizationFrequency);
250	int savePerturbation = model2->perturbation();
251	if (savePerturbation == 100)
252	model2->setPerturbation(50);
253	if (!doPrimal) {
254	// faster if bounds tightened
255	//int numberInfeasibilities = model2->tightenPrimalBounds();
256	model2->tightenPrimalBounds();
257	// look further
258	bool crashResult = false;
259	if (doCrash > 0)
260	crashResult = (solver->crash(1000.0, 1) > 0);
261	else if (doCrash == 0 && algorithm > 0)
262	crashResult = (solver->crash(1000.0, 1) > 0);
263	doPrimal = crashResult;
264	}
265	if (algorithm < 0)
266	doPrimal = false;
267	else if (algorithm > 0)
268	doPrimal = true;
269	if (!doPrimal) {
270	//if (numberInfeasibilities)
271	//std::cout<<"** Analysis indicates model infeasible"
272	// <<std::endl;
273	// up dual bound for safety
274	//model2->setDualBound(1.0e11);
275	disasterHandler_->setOsiModel(this);
276	if (inCbcOrOther) {
277	disasterHandler_->setSimplex(model2);
278	disasterHandler_->setWhereFrom(4);
279	model2->setDisasterHandler(disasterHandler_);
280	}
281	model2->dual(0);
282	totalIterations += model2->numberIterations();
283	if (inCbcOrOther) {
284	if (disasterHandler_->inTrouble()) {
285	#ifdef COIN_DEVELOP
286	printf("dual trouble a\n");
287	#endif
288	if (disasterHandler_->typeOfDisaster()) {
289	// We want to abort
290	abortSearch = true;
291	goto disaster;
292	}
293	// try just going back in
294	disasterHandler_->setPhase(1);
295	model2->dual();
296	totalIterations += model2->numberIterations();
297	if (disasterHandler_->inTrouble()) {
298	#ifdef COIN_DEVELOP
299	printf("dual trouble b\n");
300	#endif
301	if (disasterHandler_->typeOfDisaster()) {
302	// We want to abort
303	abortSearch = true;
304	goto disaster;
305	}
306	// try primal with original basis
307	disasterHandler_->setPhase(2);
308	setBasis(basis_, model2);
309	model2->primal();
310	totalIterations += model2->numberIterations();
311	}
312	if (disasterHandler_->inTrouble()) {
313	#ifdef COIN_DEVELOP
314	printf("disaster - treat as infeasible\n");
315	#endif
316	if (disasterHandler_->typeOfDisaster()) {
317	// We want to abort
318	abortSearch = true;
319	goto disaster;
320	}
321	model2->setProblemStatus(1);
322	}
323	}
324	// reset
325	model2->setDisasterHandler(NULL);
326	}
327	// check if clp thought it was in a loop
328	if (model2->status() == 3 && !model2->hitMaximumIterations()) {
329	// switch algorithm
330	disasterHandler_->setOsiModel(this);
331	if (inCbcOrOther) {
332	disasterHandler_->setSimplex(model2);
333	disasterHandler_->setWhereFrom(6);
334	model2->setDisasterHandler(disasterHandler_);
335	}
336	model2->primal();
337	totalIterations += model2->numberIterations();
338	if (inCbcOrOther) {
339	if (disasterHandler_->inTrouble()) {
340	#ifdef COIN_DEVELOP
341	printf("primal trouble a\n");
342	#endif
343	if (disasterHandler_->typeOfDisaster()) {
344	// We want to abort
345	abortSearch = true;
346	goto disaster;
347	}
348	// try just going back in (but with dual)
349	disasterHandler_->setPhase(1);
350	model2->dual();
351	totalIterations += model2->numberIterations();
352	if (disasterHandler_->inTrouble()) {
353	#ifdef COIN_DEVELOP
354	printf("primal trouble b\n");
355	#endif
356	if (disasterHandler_->typeOfDisaster()) {
357	// We want to abort
358	abortSearch = true;
359	goto disaster;
360	}
361	// try primal with original basis
362	disasterHandler_->setPhase(2);
363	setBasis(basis_, model2);
364	model2->dual();
365	totalIterations += model2->numberIterations();
366	}
367	if (disasterHandler_->inTrouble()) {
368	#ifdef COIN_DEVELOP
369	printf("disaster - treat as infeasible\n");
370	#endif
371	if (disasterHandler_->typeOfDisaster()) {
372	// We want to abort
373	abortSearch = true;
374	goto disaster;
375	}
376	model2->setProblemStatus(1);
377	}
378	}
379	// reset
380	model2->setDisasterHandler(NULL);
381	}
382	}
383	} else {
384	// up infeasibility cost for safety
385	//model2->setInfeasibilityCost(1.0e10);
386	disasterHandler_->setOsiModel(this);
387	if (inCbcOrOther) {
388	disasterHandler_->setSimplex(model2);
389	disasterHandler_->setWhereFrom(6);
390	model2->setDisasterHandler(disasterHandler_);
391	}
392	model2->primal(1);
393	totalIterations += model2->numberIterations();
394	if (inCbcOrOther) {
395	if (disasterHandler_->inTrouble()) {
396	#ifdef COIN_DEVELOP
397	printf("primal trouble a\n");
398	#endif
399	if (disasterHandler_->typeOfDisaster()) {
400	// We want to abort
401	abortSearch = true;
402	goto disaster;
403	}
404	// try just going back in (but with dual)
405	disasterHandler_->setPhase(1);
406	model2->dual();
407	totalIterations += model2->numberIterations();
408	if (disasterHandler_->inTrouble()) {
409	#ifdef COIN_DEVELOP
410	printf("primal trouble b\n");
411	#endif
412	if (disasterHandler_->typeOfDisaster()) {
413	// We want to abort
414	abortSearch = true;
415	goto disaster;
416	}
417	// try primal with original basis
418	disasterHandler_->setPhase(2);
419	setBasis(basis_, model2);
420	model2->dual();
421	totalIterations += model2->numberIterations();
422	}
423	if (disasterHandler_->inTrouble()) {
424	#ifdef COIN_DEVELOP
425	printf("disaster - treat as infeasible\n");
426	#endif
427	if (disasterHandler_->typeOfDisaster()) {
428	// We want to abort
429	abortSearch = true;
430	goto disaster;
431	}
432	model2->setProblemStatus(1);
433	}
434	}
435	// reset
436	model2->setDisasterHandler(NULL);
437	}
438	// check if clp thought it was in a loop
439	if (model2->status() == 3 && !model2->hitMaximumIterations()) {
440	// switch algorithm
441	disasterHandler_->setOsiModel(this);
442	if (inCbcOrOther) {
443	disasterHandler_->setSimplex(model2);
444	disasterHandler_->setWhereFrom(4);
445	model2->setDisasterHandler(disasterHandler_);
446	}
447	model2->dual(0);
448	totalIterations += model2->numberIterations();
449	if (inCbcOrOther) {
450	if (disasterHandler_->inTrouble()) {
451	#ifdef COIN_DEVELOP
452	printf("dual trouble a\n");
453	#endif
454	if (disasterHandler_->typeOfDisaster()) {
455	// We want to abort
456	abortSearch = true;
457	goto disaster;
458	}
459	// try just going back in
460	disasterHandler_->setPhase(1);
461	model2->dual();
462	totalIterations += model2->numberIterations();
463	if (disasterHandler_->inTrouble()) {
464	#ifdef COIN_DEVELOP
465	printf("dual trouble b\n");
466	#endif
467	if (disasterHandler_->typeOfDisaster()) {
468	// We want to abort
469	abortSearch = true;
470	goto disaster;
471	}
472	// try primal with original basis
473	disasterHandler_->setPhase(2);
474	setBasis(basis_, model2);
475	model2->primal();
476	totalIterations += model2->numberIterations();
477	}
478	if (disasterHandler_->inTrouble()) {
479	#ifdef COIN_DEVELOP
480	printf("disaster - treat as infeasible\n");
481	#endif
482	if (disasterHandler_->typeOfDisaster()) {
483	// We want to abort
484	abortSearch = true;
485	goto disaster;
486	}
487	model2->setProblemStatus(1);
488	}
489	}
490	// reset
491	model2->setDisasterHandler(NULL);
492	}
493	}
494	}
495	model2->setPerturbation(savePerturbation);
496	if (model2 != solver) {
497	int presolvedStatus = model2->status();
498	pinfo->postsolve(true);
499	delete pinfo;
500	pinfo = NULL;
501
502	delete model2;
503	int oldStatus = solver->status();
504	solver->setProblemStatus(presolvedStatus);
505	if (solver->logLevel() == 63) // for gcc 4.6 bug
506	printf("pstat %d stat %d\n", presolvedStatus, oldStatus);
507	//printf("Resolving from postsolved model\n");
508	// later try without (1) and check duals before solve
509	if (presolvedStatus != 3
510	&& (presolvedStatus \|\| oldStatus == -1)) {
511	if (!inCbcOrOther \|\| presolvedStatus != 1) {
512	disasterHandler_->setOsiModel(this);
513	if (inCbcOrOther) {
514	disasterHandler_->setSimplex(solver); // as "borrowed"
515	disasterHandler_->setWhereFrom(6);
516	solver->setDisasterHandler(disasterHandler_);
517	}
518	solver->primal(1);
519	totalIterations += solver->numberIterations();
520	if (inCbcOrOther) {
521	if (disasterHandler_->inTrouble()) {
522	#ifdef COIN_DEVELOP
523	printf("primal trouble a\n");
524	#endif
525	if (disasterHandler_->typeOfDisaster()) {
526	// We want to abort
527	abortSearch = true;
528	goto disaster;
529	}
530	// try just going back in (but with dual)
531	disasterHandler_->setPhase(1);
532	solver->dual();
533	totalIterations += solver->numberIterations();
534	if (disasterHandler_->inTrouble()) {
535	#ifdef COIN_DEVELOP
536	printf("primal trouble b\n");
537	#endif
538	if (disasterHandler_->typeOfDisaster()) {
539	// We want to abort
540	abortSearch = true;
541	goto disaster;
542	}
543	// try primal with original basis
544	disasterHandler_->setPhase(2);
545	setBasis(basis_, solver);
546	solver->dual();
547	totalIterations += solver->numberIterations();
548	}
549	if (disasterHandler_->inTrouble()) {
550	#ifdef COIN_DEVELOP
551	printf("disaster - treat as infeasible\n");
552	#endif
553	if (disasterHandler_->typeOfDisaster()) {
554	// We want to abort
555	abortSearch = true;
556	goto disaster;
557	}
558	solver->setProblemStatus(1);
559	}
560	}
561	// reset
562	solver->setDisasterHandler(NULL);
563	}
564	}
565	}
566	}
567	lastAlgorithm_ = 1; // primal
568	//if (solver->numberIterations())
569	//printf("****** iterated %d\n",solver->numberIterations());
570	} else {
571	// do we want crash
572	if (doCrash > 0)
573	solver->crash(1000.0, 2);
574	else if (doCrash == 0)
575	solver->crash(1000.0, 0);
576	if (algorithm < 0)
577	doPrimal = false;
578	else if (algorithm > 0)
579	doPrimal = true;
580	disasterHandler_->setOsiModel(this);
581	disasterHandler_->setSimplex(solver); // as "borrowed"
582	bool inCbcOrOther = (modelPtr_->specialOptions() & 0x03000000) != 0;
583	if (!doPrimal)
584	disasterHandler_->setWhereFrom(4);
585	else
586	disasterHandler_->setWhereFrom(6);
587	if (inCbcOrOther)
588	solver->setDisasterHandler(disasterHandler_);
589	if (!doPrimal) {
590	//printf("doing dual\n");
591	solver->dual(0);
592	totalIterations += solver->numberIterations();
593	if (inCbcOrOther) {
594	if (disasterHandler_->inTrouble()) {
595	#ifdef COIN_DEVELOP
596	printf("dual trouble a\n");
597	#endif
598	if (disasterHandler_->typeOfDisaster()) {
599	// We want to abort
600	abortSearch = true;
601	goto disaster;
602	}
603	// try just going back in
604	disasterHandler_->setPhase(1);
605	solver->dual();
606	totalIterations += solver->numberIterations();
607	if (disasterHandler_->inTrouble()) {
608	#ifdef COIN_DEVELOP
609	printf("dual trouble b\n");
610	#endif
611	if (disasterHandler_->typeOfDisaster()) {
612	// We want to abort
613	abortSearch = true;
614	goto disaster;
615	}
616	// try primal with original basis
617	disasterHandler_->setPhase(2);
618	setBasis(basis_, solver);
619	solver->primal();
620	totalIterations += solver->numberIterations();
621	}
622	if (disasterHandler_->inTrouble()) {
623	#ifdef COIN_DEVELOP
624	printf("disaster - treat as infeasible\n");
625	#endif
626	if (disasterHandler_->typeOfDisaster()) {
627	// We want to abort
628	abortSearch = true;
629	goto disaster;
630	}
631	solver->setProblemStatus(1);
632	}
633	}
634	// reset
635	solver->setDisasterHandler(NULL);
636	}
637	lastAlgorithm_ = 2; // dual
638	// check if clp thought it was in a loop
639	if (solver->status() == 3 && !solver->hitMaximumIterations()) {
640	// switch algorithm
641	solver->primal(0);
642	totalIterations += solver->numberIterations();
643	lastAlgorithm_ = 1; // primal
644	}
645	} else {
646	//printf("doing primal\n");
647	solver->primal(1);
648	totalIterations += solver->numberIterations();
649	if (inCbcOrOther) {
650	if (disasterHandler_->inTrouble()) {
651	#ifdef COIN_DEVELOP
652	printf("primal trouble a\n");
653	#endif
654	if (disasterHandler_->typeOfDisaster()) {
655	// We want to abort
656	abortSearch = true;
657	goto disaster;
658	}
659	// try just going back in (but with dual)
660	disasterHandler_->setPhase(1);
661	solver->dual();
662	totalIterations += solver->numberIterations();
663	if (disasterHandler_->inTrouble()) {
664	#ifdef COIN_DEVELOP
665	printf("primal trouble b\n");
666	#endif
667	if (disasterHandler_->typeOfDisaster()) {
668	// We want to abort
669	abortSearch = true;
670	goto disaster;
671	}
672	// try primal with original basis
673	disasterHandler_->setPhase(2);
674	setBasis(basis_, solver);
675	solver->dual();
676	totalIterations += solver->numberIterations();
677	}
678	if (disasterHandler_->inTrouble()) {
679	#ifdef COIN_DEVELOP
680	printf("disaster - treat as infeasible\n");
681	#endif
682	if (disasterHandler_->typeOfDisaster()) {
683	// We want to abort
684	abortSearch = true;
685	goto disaster;
686	}
687	solver->setProblemStatus(1);
688	}
689	}
690	// reset
691	solver->setDisasterHandler(NULL);
692	}
693	lastAlgorithm_ = 1; // primal
694	// check if clp thought it was in a loop
695	if (solver->status() == 3 && !solver->hitMaximumIterations()) {
696	// switch algorithm
697	solver->dual(0);
698	totalIterations += solver->numberIterations();
699	lastAlgorithm_ = 2; // dual
700	}
701	}
702	}
703	// If scaled feasible but unscaled infeasible take action
704	if (!solver->status() && cleanupScaling_) {
705	solver->cleanup(cleanupScaling_);
706	}
707	basis_ = getBasis(solver);
708	//basis_.print();
709	const double *rowScale2 = solver->rowScale();
710	solver->setSpecialOptions(saveOptions);
711	if (!rowScale1 && rowScale2) {
712	// need to release memory
713	if (!solver->savedRowScale_) {
714	solver->setRowScale(NULL);
715	solver->setColumnScale(NULL);
716	} else {
717	solver->rowScale_ = NULL;
718	solver->columnScale_ = NULL;
719	}
720	}
721	} else {
722	// User doing nothing and all slack basis
723	ClpSolve options = solveOptions_;
724	bool yesNo;
725	OsiHintStrength strength;
726	getHintParam(OsiDoInBranchAndCut, yesNo, strength);
727	if (yesNo) {
728	solver->setSpecialOptions(solver->specialOptions() \| 1024);
729	}
730	solver->initialSolve(options);
731	totalIterations += solver->numberIterations();
732	lastAlgorithm_ = 2; // say dual
733	// If scaled feasible but unscaled infeasible take action
734	if (!solver->status() && cleanupScaling_) {
735	solver->cleanup(cleanupScaling_);
736	}
737	basis_ = getBasis(solver);
738	//basis_.print();
739	}
740	solver->messageHandler()->setLogLevel(saveMessageLevel);
741	disaster:
742	if (deleteSolver) {
743	solver->returnModel(*modelPtr_);
744	delete solver;
745	}
746	if (startFinishOptions) {
747	int save = modelPtr_->logLevel();
748	if (save < 2)
749	modelPtr_->setLogLevel(0);
750	modelPtr_->dual(0, startFinishOptions);
751	totalIterations += modelPtr_->numberIterations();
752	modelPtr_->setLogLevel(save);
753	}
754	if (saveSolveType == 2) {
755	enableSimplexInterface(doingPrimal);
756	}
757	if (savedObjective) {
758	// fix up
759	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
760	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
761	modelPtr_->objective_ = savedObjective;
762	if (!modelPtr_->problemStatus_) {
763	CoinZeroN(modelPtr_->dual_, modelPtr_->numberRows_);
764	CoinZeroN(modelPtr_->reducedCost_, modelPtr_->numberColumns_);
765	if (modelPtr_->dj_ && (modelPtr_->whatsChanged_ & 1) != 0)
766	CoinZeroN(modelPtr_->dj_, modelPtr_->numberColumns_ + modelPtr_->numberRows_);
767	modelPtr_->computeObjectiveValue();
768	}
769	}
770	modelPtr_->setNumberIterations(totalIterations);
771	handler_->setLogLevel(saveMessageLevel2);
772	if (modelPtr_->problemStatus_ == 3 && lastAlgorithm_ == 2)
773	modelPtr_->computeObjectiveValue();
774	// mark so we can pick up objective value quickly
775	modelPtr_->upperIn_ = 0.0;
776	#ifdef PRINT_TIME
777	time1 = CoinCpuTime() - time1;
778	totalTime += time1;
779	#endif
780	assert(!modelPtr_->disasterHandler());
781	if (lastAlgorithm_ < 1 \|\| lastAlgorithm_ > 2)
782	lastAlgorithm_ = 1;
783	if (abortSearch) {
784	lastAlgorithm_ = -911;
785	modelPtr_->setProblemStatus(4);
786	}
787	modelPtr_->whatsChanged_ \|= 0x30000;
788	#if 0
789	// delete scaled matrix and rowcopy for safety
790	delete modelPtr_->scaledMatrix_;
791	modelPtr_->scaledMatrix_=NULL;
792	delete modelPtr_->rowCopy_;
793	modelPtr_->rowCopy_=NULL;
794	#endif
795	#ifdef PRINT_TIME
796	std::cout << time1 << " seconds - total " << totalTime << std::endl;
797	#endif
798	delete pinfo;
799	}
800	//-----------------------------------------------------------------------------
801	void OsiClpSolverInterface::resolve()
802	{
803	#ifdef COIN_DEVELOP
804	{
805	int i;
806	int n = getNumCols();
807	const double *lower = getColLower();
808	const double *upper = getColUpper();
809	for (i = 0; i < n; i++) {
810	assert(lower[i] < 1.0e12);
811	assert(upper[i] > -1.0e12);
812	}
813	n = getNumRows();
814	lower = getRowLower();
815	upper = getRowUpper();
816	for (i = 0; i < n; i++) {
817	assert(lower[i] < 1.0e12);
818	assert(upper[i] > -1.0e12);
819	}
820	}
821	#endif
822	if ((stuff_.solverOptions_ & 65536) != 0) {
823	modelPtr_->fastDual2(&stuff_);
824	return;
825	}
826	if ((specialOptions_ & 2097152) != 0 \|\| (specialOptions_ & 4194304) != 0) {
827	bool takeHint;
828	OsiHintStrength strength;
829	int algorithm = 0;
830	getHintParam(OsiDoDualInResolve, takeHint, strength);
831	if (strength != OsiHintIgnore)
832	algorithm = takeHint ? -1 : 1;
833	if (algorithm > 0 \|\| (specialOptions_ & 4194304) != 0) {
834	// Gub
835	resolveGub((9 * modelPtr_->numberRows()) / 10);
836	return;
837	}
838	}
839	//void pclp(char *);
840	//pclp("res");
841	bool takeHint;
842	OsiHintStrength strength;
843	bool gotHint = (getHintParam(OsiDoInBranchAndCut, takeHint, strength));
844	assert(gotHint);
845	// mark so we can pick up objective value quickly
846	modelPtr_->upperIn_ = 0.0;
847	if ((specialOptions_ & 4096) != 0) {
848	// Quick check to see if optimal
849	modelPtr_->checkSolutionInternal();
850	if (modelPtr_->problemStatus() == 0) {
851	modelPtr_->setNumberIterations(0);
852	return;
853	}
854	}
855	int totalIterations = 0;
856	bool abortSearch = false;
857	ClpObjective *savedObjective = NULL;
858	double savedDualLimit = modelPtr_->dblParam_[ClpDualObjectiveLimit];
859	if (fakeObjective_) {
860	modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions() & (~128));
861	// See if all with costs fixed
862	int numberColumns = modelPtr_->numberColumns_;
863	const double *obj = modelPtr_->objective();
864	const double *lower = modelPtr_->columnLower();
865	const double *upper = modelPtr_->columnUpper();
866	int i;
867	for (i = 0; i < numberColumns; i++) {
868	double objValue = obj[i];
869	if (objValue) {
870	if (lower[i] != upper[i])
871	break;
872	}
873	}
874	if (i == numberColumns) {
875	if ((specialOptions_ & 524288) == 0) {
876	// Set fake
877	savedObjective = modelPtr_->objective_;
878	modelPtr_->objective_ = fakeObjective_;
879	modelPtr_->dblParam_[ClpDualObjectiveLimit] = COIN_DBL_MAX;
880	} else {
881	modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions() \| 128);
882	}
883	}
884	}
885	// If using Clp initialSolve and primal - just do here
886	gotHint = (getHintParam(OsiDoDualInResolve, takeHint, strength));
887	assert(gotHint);
888	if (strength != OsiHintIgnore && !takeHint && solveOptions_.getSpecialOption(6)) {
889	ClpSolve options = solveOptions_;
890	// presolve
891	getHintParam(OsiDoPresolveInResolve, takeHint, strength);
892	if (strength != OsiHintIgnore && !takeHint)
893	options.setPresolveType(ClpSolve::presolveOff);
894	int saveOptions = modelPtr_->specialOptions();
895	getHintParam(OsiDoInBranchAndCut, takeHint, strength);
896	if (takeHint) {
897	modelPtr_->setSpecialOptions(modelPtr_->specialOptions() \| 1024);
898	}
899	setBasis(basis_, modelPtr_);
900	modelPtr_->initialSolve(options);
901	lastAlgorithm_ = 1; // say primal
902	// If scaled feasible but unscaled infeasible take action
903	if (!modelPtr_->status() && cleanupScaling_) {
904	modelPtr_->cleanup(cleanupScaling_);
905	}
906	modelPtr_->setSpecialOptions(saveOptions); // restore
907	basis_ = getBasis(modelPtr_);
908	}
909	int saveSolveType = modelPtr_->solveType();
910	bool doingPrimal = modelPtr_->algorithm() > 0;
911	if (saveSolveType == 2) {
912	disableSimplexInterface();
913	}
914	int saveOptions = modelPtr_->specialOptions();
915	int startFinishOptions = 0;
916	if (specialOptions_ != 0x80000000) {
917	if ((specialOptions_ & 1) == 0) {
918	startFinishOptions = 0;
919	modelPtr_->setSpecialOptions(saveOptions \| (64 \| 1024 \| 32768));
920	} else {
921	startFinishOptions = 1 + 4;
922	if ((specialOptions_ & 8) != 0)
923	startFinishOptions += 2; // allow re-use of factorization
924	if ((specialOptions_ & 4) == 0 \|\| !takeHint)
925	modelPtr_->setSpecialOptions(saveOptions \| (64 \| 128 \| 512 \| 1024 \| 4096 \| 32768));
926	else
927	modelPtr_->setSpecialOptions(saveOptions \| (64 \| 128 \| 512 \| 1024 \| 2048 \| 4096 \| 32768));
928	}
929	} else {
930	modelPtr_->setSpecialOptions(saveOptions \| 64 \| 32768);
931	}
932	//printf("options %d size %d\n",modelPtr_->specialOptions(),modelPtr_->numberColumns());
933	//modelPtr_->setSolveType(1);
934	// Set message handler to have same levels etc
935	int saveMessageLevel = modelPtr_->logLevel();
936	int messageLevel = messageHandler()->logLevel();
937	bool oldDefault;
938	CoinMessageHandler *saveHandler = NULL;
939	if (!defaultHandler_)
940	saveHandler = modelPtr_->pushMessageHandler(handler_, oldDefault);
941	//printf("basis before dual\n");
942	//basis_.print();
943	setBasis(basis_, modelPtr_);
944	#ifdef SAVE_MODEL
945	resolveTry++;
946	#if SAVE_MODEL > 1
947	if (resolveTry >= loResolveTry && resolveTry <= hiResolveTry) {
948	char fileName[20];
949	sprintf(fileName, "save%d.mod", resolveTry);
950	modelPtr_->saveModel(fileName);
951	}
952	#endif
953	#endif
954	// set reasonable defaults
955	// Switch off printing if asked to
956	gotHint = (getHintParam(OsiDoReducePrint, takeHint, strength));
957	assert(gotHint);
958	if (strength != OsiHintIgnore && takeHint) {
959	if (messageLevel > 0)
960	messageLevel--;
961	}
962	if (messageLevel < modelPtr_->messageHandler()->logLevel())
963	modelPtr_->messageHandler()->setLogLevel(messageLevel);
964	// See if user set factorization frequency
965	int userFactorizationFrequency = modelPtr_->factorization()->maximumPivots();
966	// scaling
967	if (modelPtr_->solveType() == 1) {
968	gotHint = (getHintParam(OsiDoScale, takeHint, strength));
969	assert(gotHint);
970	if (strength == OsiHintIgnore \|\| takeHint) {
971	if (!modelPtr_->scalingFlag())
972	modelPtr_->scaling(3);
973	} else {
974	modelPtr_->scaling(0);
975	}
976	} else {
977	modelPtr_->scaling(0);
978	}
979	// sort out hints;
980	// algorithm -1 force dual, +1 force primal
981	int algorithm = -1;
982	gotHint = (getHintParam(OsiDoDualInResolve, takeHint, strength));
983	assert(gotHint);
984	if (strength != OsiHintIgnore)
985	algorithm = takeHint ? -1 : 1;
986	//modelPtr_->saveModel("save.bad");
987	// presolve
988	gotHint = (getHintParam(OsiDoPresolveInResolve, takeHint, strength));
989	assert(gotHint);
990	if (strength != OsiHintIgnore && takeHint) {
991	#ifdef KEEP_SMALL
992	if (smallModel_) {
993	delete[] spareArrays_;
994	spareArrays_ = NULL;
995	delete smallModel_;
996	smallModel_ = NULL;
997	}
998	#endif
999	ClpPresolve pinfo;
1000	if ((specialOptions_ & 128) != 0) {
1001	specialOptions_ &= ~128;
1002	}
1003	if ((modelPtr_->specialOptions() & 1024) != 0) {
1004	pinfo.setDoDual(false);
1005	pinfo.setDoTripleton(false);
1006	pinfo.setDoDupcol(false);
1007	pinfo.setDoDuprow(false);
1008	pinfo.setDoSingletonColumn(false);
1009	}
1010	ClpSimplex model2 = pinfo.presolvedModel(modelPtr_, 1.0e-8);
1011	if (!model2) {
1012	// problem found to be infeasible - whats best?
1013	model2 = modelPtr_;
1014	}
1015	// return number of rows
1016	int stats = reinterpret_cast< int >(getApplicationData());
1017	if (stats) {
1018	stats[0] = model2->numberRows();
1019	stats[1] = model2->numberColumns();
1020	}
1021	//printf("rows %d -> %d, columns %d -> %d\n",
1022	// modelPtr_->numberRows(),model2->numberRows(),
1023	// modelPtr_->numberColumns(),model2->numberColumns());
1024	// change from 200
1025	if (modelPtr_->factorization()->maximumPivots() == 200)
1026	model2->factorization()->maximumPivots(100 + model2->numberRows() / 50);
1027	else
1028	model2->factorization()->maximumPivots(userFactorizationFrequency);
1029	if (algorithm < 0) {
1030	model2->dual();
1031	totalIterations += model2->numberIterations();
1032	// check if clp thought it was in a loop
1033	if (model2->status() == 3 && !model2->hitMaximumIterations()) {
1034	// switch algorithm
1035	model2->primal();
1036	totalIterations += model2->numberIterations();
1037	}
1038	} else {
1039	model2->primal(1);
1040	totalIterations += model2->numberIterations();
1041	// check if clp thought it was in a loop
1042	if (model2->status() == 3 && !model2->hitMaximumIterations()) {
1043	// switch algorithm
1044	model2->dual();
1045	totalIterations += model2->numberIterations();
1046	}
1047	}
1048	if (model2 != modelPtr_) {
1049	int finalStatus = model2->status();
1050	pinfo.postsolve(true);
1051
1052	delete model2;
1053	// later try without (1) and check duals before solve
1054	if (finalStatus != 3 && (finalStatus \|\| modelPtr_->status() == -1)) {
1055	modelPtr_->primal(1);
1056	totalIterations += modelPtr_->numberIterations();
1057	lastAlgorithm_ = 1; // primal
1058	//if (modelPtr_->numberIterations())
1059	//printf("****** iterated %d\n",modelPtr_->numberIterations());
1060	}
1061	}
1062	} else {
1063	//modelPtr_->setLogLevel(63);
1064	//modelPtr_->setDualTolerance(1.0e-7);
1065	if (false && modelPtr_->scalingFlag_ > 0 && !modelPtr_->rowScale_ && !modelPtr_->rowCopy_ && matrixByRow_) {
1066	assert(matrixByRow_->getNumElements() == modelPtr_->clpMatrix()->getNumElements());
1067	modelPtr_->setNewRowCopy(new ClpPackedMatrix(*matrixByRow_));
1068	}
1069	if (algorithm < 0) {
1070	//writeMps("try1");
1071	int savePerturbation = modelPtr_->perturbation();
1072	if ((specialOptions_ & 2) != 0)
1073	modelPtr_->setPerturbation(100);
1074	//modelPtr_->messageHandler()->setLogLevel(1);
1075	//writeMpsNative("bad",NULL,NULL,2,1,1.0);
1076	disasterHandler_->setOsiModel(this);
1077	bool inCbcOrOther = (modelPtr_->specialOptions() & 0x03000000) != 0;
1078	#if 0
1079	// See how many integers fixed
1080	bool skipCrunch=true;
1081	const char * integerInformation = modelPtr_->integerType_;
1082	if (integerInformation) {
1083	int numberColumns = modelPtr_->numberColumns_;
1084	const double * lower = modelPtr_->columnLower();
1085	const double * upper = modelPtr_->columnUpper();
1086	int target=CoinMax(1,numberColumns/10000);
1087	for (int i=0;i<numberColumns;i++) {
1088	if (integerInformation[i]) {
1089	if (lower[i]==upper[i]) {
1090	target--;
1091	if (!target) {
1092	skipCrunch=false;
1093	break;
1094	}
1095	}
1096	}
1097	}
1098	}
1099	#endif
1100	if ((specialOptions_ & 1) == 0 \|\| (specialOptions_ & 2048) != 0 \|\| (modelPtr_->specialOptions_ & 2097152) != 0 /\|\|skipCrunch/) {
1101	disasterHandler_->setWhereFrom(0); // dual
1102	if (inCbcOrOther)
1103	modelPtr_->setDisasterHandler(disasterHandler_);
1104	bool specialScale;
1105	if ((specialOptions_ & 131072) != 0 && !modelPtr_->rowScale_) {
1106	modelPtr_->rowScale_ = rowScale_.array();
1107	modelPtr_->columnScale_ = columnScale_.array();
1108	specialScale = true;
1109	} else {
1110	specialScale = false;
1111	}
1112	#ifdef KEEP_SMALL
1113	if (smallModel_) {
1114	delete[] spareArrays_;
1115	spareArrays_ = NULL;
1116	delete smallModel_;
1117	smallModel_ = NULL;
1118	}
1119	#endif
1120	#ifdef CBC_STATISTICS
1121	osi_dual++;
1122	#endif
1123	modelPtr_->dual(0, startFinishOptions);
1124	totalIterations += modelPtr_->numberIterations();
1125	if (specialScale) {
1126	modelPtr_->rowScale_ = NULL;
1127	modelPtr_->columnScale_ = NULL;
1128	}
1129	} else {
1130	#ifdef CBC_STATISTICS
1131	osi_crunch++;
1132	#endif
1133	crunch();
1134	totalIterations += modelPtr_->numberIterations();
1135	if (modelPtr_->problemStatus() == 4)
1136	goto disaster;
1137	// should have already been fixed if problems
1138	inCbcOrOther = false;
1139	}
1140	if (inCbcOrOther) {
1141	if (disasterHandler_->inTrouble()) {
1142	if (disasterHandler_->typeOfDisaster()) {
1143	// We want to abort
1144	abortSearch = true;
1145	goto disaster;
1146	}
1147	// try just going back in
1148	disasterHandler_->setPhase(1);
1149	modelPtr_->dual();
1150	totalIterations += modelPtr_->numberIterations();
1151	if (disasterHandler_->inTrouble()) {
1152	if (disasterHandler_->typeOfDisaster()) {
1153	// We want to abort
1154	abortSearch = true;
1155	goto disaster;
1156	}
1157	// try primal with original basis
1158	disasterHandler_->setPhase(2);
1159	setBasis(basis_, modelPtr_);
1160	modelPtr_->primal();
1161	totalIterations += modelPtr_->numberIterations();
1162	}
1163	if (disasterHandler_->inTrouble()) {
1164	#ifdef COIN_DEVELOP
1165	printf("disaster - treat as infeasible\n");
1166	#endif
1167	if (disasterHandler_->typeOfDisaster()) {
1168	// We want to abort
1169	abortSearch = true;
1170	goto disaster;
1171	}
1172	modelPtr_->setProblemStatus(1);
1173	}
1174	}
1175	// reset
1176	modelPtr_->setDisasterHandler(NULL);
1177	}
1178	if (modelPtr_->problemStatus() == 4) {
1179	// bad bounds?
1180	modelPtr_->setProblemStatus(1);
1181	}
1182	if (!modelPtr_->problemStatus() && 0) {
1183	int numberColumns = modelPtr_->numberColumns();
1184	const double *columnLower = modelPtr_->columnLower();
1185	const double *columnUpper = modelPtr_->columnUpper();
1186	int nBad = 0;
1187	for (int i = 0; i < numberColumns; i++) {
1188	if (columnLower[i] == columnUpper[i] && modelPtr_->getColumnStatus(i) == ClpSimplex::basic) {
1189	nBad++;
1190	modelPtr_->setColumnStatus(i, ClpSimplex::isFixed);
1191	}
1192	}
1193	if (nBad) {
1194	modelPtr_->primal(1);
1195	totalIterations += modelPtr_->numberIterations();
1196	printf("%d fixed basic - %d iterations\n", nBad, modelPtr_->numberIterations());
1197	}
1198	}
1199	assert(modelPtr_->objectiveValue() < 1.0e100);
1200	modelPtr_->setPerturbation(savePerturbation);
1201	lastAlgorithm_ = 2; // dual
1202	// check if clp thought it was in a loop
1203	if (modelPtr_->status() == 3 && !modelPtr_->hitMaximumIterations()) {
1204	modelPtr_->setSpecialOptions(saveOptions);
1205	// switch algorithm
1206	//modelPtr_->messageHandler()->setLogLevel(63);
1207	// Allow for catastrophe
1208	int saveMax = modelPtr_->maximumIterations();
1209	int numberIterations = modelPtr_->numberIterations();
1210	int numberRows = modelPtr_->numberRows();
1211	int numberColumns = modelPtr_->numberColumns();
1212	if (modelPtr_->maximumIterations() > 100000 + numberIterations)
1213	modelPtr_->setMaximumIterations(numberIterations + 1000 + 2 * numberRows + numberColumns);
1214	modelPtr_->primal(0, startFinishOptions);
1215	totalIterations += modelPtr_->numberIterations();
1216	modelPtr_->setMaximumIterations(saveMax);
1217	lastAlgorithm_ = 1; // primal
1218	if (modelPtr_->status() == 3 && !modelPtr_->hitMaximumIterations()) {
1219	#ifdef COIN_DEVELOP
1220	printf("in trouble - try all slack\n");
1221	#endif
1222	CoinWarmStartBasis allSlack;
1223	setBasis(allSlack, modelPtr_);
1224	modelPtr_->dual();
1225	totalIterations += modelPtr_->numberIterations();
1226	if (modelPtr_->status() == 3 && !modelPtr_->hitMaximumIterations()) {
1227	if (modelPtr_->numberPrimalInfeasibilities()) {
1228	#ifdef COIN_DEVELOP
1229	printf("Real real trouble - treat as infeasible\n");
1230	#endif
1231	modelPtr_->setProblemStatus(1);
1232	} else {
1233	#ifdef COIN_DEVELOP
1234	printf("Real real trouble - treat as optimal\n");
1235	#endif
1236	modelPtr_->setProblemStatus(0);
1237	}
1238	}
1239	}
1240	}
1241	assert(modelPtr_->objectiveValue() < 1.0e100);
1242	} else {
1243	#ifdef KEEP_SMALL
1244	if (smallModel_) {
1245	delete[] spareArrays_;
1246	spareArrays_ = NULL;
1247	delete smallModel_;
1248	smallModel_ = NULL;
1249	}
1250	#endif
1251	//printf("doing primal\n");
1252	#ifdef CBC_STATISTICS
1253	osi_primal++;
1254	#endif
1255	modelPtr_->primal(1, startFinishOptions);
1256	totalIterations += modelPtr_->numberIterations();
1257	lastAlgorithm_ = 1; // primal
1258	// check if clp thought it was in a loop
1259	if (modelPtr_->status() == 3 && !modelPtr_->hitMaximumIterations()) {
1260	// switch algorithm
1261	modelPtr_->dual();
1262	totalIterations += modelPtr_->numberIterations();
1263	lastAlgorithm_ = 2; // dual
1264	}
1265	}
1266	}
1267	// If scaled feasible but unscaled infeasible take action
1268	//if (!modelPtr_->status()&&cleanupScaling_) {
1269	if (cleanupScaling_) {
1270	modelPtr_->cleanup(cleanupScaling_);
1271	}
1272	basis_ = getBasis(modelPtr_);
1273	disaster:
1274	//printf("basis after dual\n");
1275	//basis_.print();
1276	if (!defaultHandler_)
1277	modelPtr_->popMessageHandler(saveHandler, oldDefault);
1278	modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
1279	if (saveSolveType == 2) {
1280	int saveStatus = modelPtr_->problemStatus_;
1281	enableSimplexInterface(doingPrimal);
1282	modelPtr_->problemStatus_ = saveStatus;
1283	}
1284	#ifdef COIN_DEVELOP_x
1285	extern bool doingDoneBranch;
1286	if (doingDoneBranch) {
1287	if (modelPtr_->numberIterations())
1288	printf("***** done %d iterations after general\n", modelPtr_->numberIterations());
1289	doingDoneBranch = false;
1290	}
1291	#endif
1292	modelPtr_->setNumberIterations(totalIterations);
1293	//modelPtr_->setSolveType(saveSolveType);
1294	if (abortSearch) {
1295	lastAlgorithm_ = -911;
1296	modelPtr_->setProblemStatus(4);
1297	}
1298	if (savedObjective) {
1299	// fix up
1300	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
1301	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
1302	modelPtr_->objective_ = savedObjective;
1303	if (!modelPtr_->problemStatus_) {
1304	CoinZeroN(modelPtr_->dual_, modelPtr_->numberRows_);
1305	CoinZeroN(modelPtr_->reducedCost_, modelPtr_->numberColumns_);
1306	if (modelPtr_->dj_ && (modelPtr_->whatsChanged_ & 1) != 0)
1307	CoinZeroN(modelPtr_->dj_, modelPtr_->numberColumns_ + modelPtr_->numberRows_);
1308	modelPtr_->computeObjectiveValue();
1309	}
1310	}
1311	modelPtr_->setSpecialOptions(saveOptions); // restore
1312	if (modelPtr_->problemStatus_ == 3 && lastAlgorithm_ == 2)
1313	modelPtr_->computeObjectiveValue();
1314	if (lastAlgorithm_ < 1 \|\| lastAlgorithm_ > 2)
1315	lastAlgorithm_ = 1;
1316	#ifdef SAVE_MODEL
1317	if (resolveTry >= loResolveTry && resolveTry <= hiResolveTry) {
1318	printf("resolve %d took %d iterations - algorithm %d\n", resolveTry, modelPtr_->numberIterations(), lastAlgorithm_);
1319	}
1320	#endif
1321	// Make sure whatsChanged not out of sync
1322	if (!modelPtr_->columnUpperWork_)
1323	modelPtr_->whatsChanged_ &= ~0xffff;
1324	modelPtr_->whatsChanged_ \|= 0x30000;
1325	}
1326	#include "ClpSimplexOther.hpp"
1327	// Resolve an LP relaxation after problem modification (try GUB)
1328	void OsiClpSolverInterface::resolveGub(int needed)
1329	{
1330	bool takeHint;
1331	OsiHintStrength strength;
1332	// Switch off printing if asked to
1333	getHintParam(OsiDoReducePrint, takeHint, strength);
1334	int saveMessageLevel = modelPtr_->logLevel();
1335	if (strength != OsiHintIgnore && takeHint) {
1336	int messageLevel = messageHandler()->logLevel();
1337	if (messageLevel > 0)
1338	modelPtr_->messageHandler()->setLogLevel(messageLevel - 1);
1339	else
1340	modelPtr_->messageHandler()->setLogLevel(0);
1341	}
1342	//modelPtr_->messageHandler()->setLogLevel(1);
1343	setBasis(basis_, modelPtr_);
1344	// find gub
1345	int numberRows = modelPtr_->numberRows();
1346	int *which = new int[numberRows];
1347	int numberColumns = modelPtr_->numberColumns();
1348	int *whichC = new int[numberColumns + numberRows];
1349	ClpSimplex model2 = static_cast< ClpSimplexOther >(modelPtr_)->gubVersion(which, whichC,
1350	needed, 100);
1351	if (model2) {
1352	// move in solution
1353	static_cast< ClpSimplexOther * >(model2)->setGubBasis(*modelPtr_,
1354	which, whichC);
1355	model2->setLogLevel(CoinMin(1, model2->logLevel()));
1356	ClpPrimalColumnSteepest steepest(5);
1357	model2->setPrimalColumnPivotAlgorithm(steepest);
1358	//double time1 = CoinCpuTime();
1359	model2->primal();
1360	//printf("Primal took %g seconds\n",CoinCpuTime()-time1);
1361	static_cast< ClpSimplexOther * >(model2)->getGubBasis(*modelPtr_,
1362	which, whichC);
1363	int totalIterations = model2->numberIterations();
1364	delete model2;
1365	//modelPtr_->setLogLevel(63);
1366	modelPtr_->primal(1);
1367	modelPtr_->setNumberIterations(totalIterations + modelPtr_->numberIterations());
1368	} else {
1369	modelPtr_->dual();
1370	}
1371	delete[] which;
1372	delete[] whichC;
1373	basis_ = getBasis(modelPtr_);
1374	modelPtr_->messageHandler()->setLogLevel(saveMessageLevel);
1375	}
1376	// Sort of lexicographic resolve
1377	void OsiClpSolverInterface::lexSolve()
1378	{
1379	#if 1
1380	abort();
1381	#else
1382	((ClpSimplexPrimal *)modelPtr_)->lexSolve();
1383	printf("itA %d\n", modelPtr_->numberIterations());
1384	modelPtr_->primal();
1385	printf("itB %d\n", modelPtr_->numberIterations());
1386	basis_ = getBasis(modelPtr_);
1387	#endif
1388	}
1389
1390	/* Sets up solver for repeated use by Osi interface.
1391	The normal usage does things like keeping factorization around so can be used.
1392	Will also do things like keep scaling and row copy of matrix if
1393	matrix does not change.
1394	adventure:
1395	0 - safe stuff as above
1396	1 - will take more risks - if it does not work then bug which will be fixed
1397	2 - don't bother doing most extreme termination checks e.g. don't bother
1398	re-factorizing if less than 20 iterations.
1399	3 - Actually safer than 1 (mainly just keeps factorization)
1400
1401	printOut - -1 always skip round common messages instead of doing some work
1402	0 skip if normal defaults
1403	1 leaves
1404	*/
1405	void OsiClpSolverInterface::setupForRepeatedUse(int senseOfAdventure, int printOut)
1406	{
1407	// First try
1408	switch (senseOfAdventure) {
1409	case 0:
1410	specialOptions_ = 8;
1411	break;
1412	case 1:
1413	specialOptions_ = 1 + 2 + 8;
1414	break;
1415	case 2:
1416	specialOptions_ = 1 + 2 + 4 + 8;
1417	break;
1418	case 3:
1419	specialOptions_ = 1 + 8;
1420	break;
1421	}
1422	//#define NO_CRUNCH2
1423	#ifdef NO_CRUNCH2
1424	specialOptions_ &= ~1;
1425	#endif
1426	bool stopPrinting = false;
1427	if (printOut < 0) {
1428	stopPrinting = true;
1429	} else if (!printOut) {
1430	bool takeHint;
1431	OsiHintStrength strength;
1432	getHintParam(OsiDoReducePrint, takeHint, strength);
1433	int messageLevel = messageHandler()->logLevel();
1434	if (strength != OsiHintIgnore && takeHint)
1435	messageLevel--;
1436	stopPrinting = (messageLevel <= 0);
1437	}
1438	#ifndef COIN_DEVELOP
1439	if (stopPrinting) {
1440	CoinMessages *messagesPointer = modelPtr_->messagesPointer();
1441	// won't even build messages
1442	messagesPointer->setDetailMessages(100, 10000, reinterpret_cast< int * >(NULL));
1443	}
1444	#endif
1445	}
1446	#ifndef NDEBUG
1447	// For errors to make sure print to screen
1448	// only called in debug mode
1449	static void indexError(int index,
1450	std::string methodName)
1451	{
1452	std::cerr << "Illegal index " << index << " in OsiClpSolverInterface::" << methodName << std::endl;
1453	throw CoinError("Illegal index", methodName, "OsiClpSolverInterface");
1454	}
1455	#endif
1456	//#############################################################################
1457	// Parameter related methods
1458	//#############################################################################
1459
1460	bool OsiClpSolverInterface::setIntParam(OsiIntParam key, int value)
1461	{
1462	return modelPtr_->setIntParam(static_cast< ClpIntParam >(key), value);
1463	}
1464
1465	//-----------------------------------------------------------------------------
1466
1467	bool OsiClpSolverInterface::setDblParam(OsiDblParam key, double value)
1468	{
1469	if (key != OsiLastDblParam) {
1470	if (key == OsiDualObjectiveLimit \|\| key == OsiPrimalObjectiveLimit)
1471	value *= modelPtr_->optimizationDirection();
1472	return modelPtr_->setDblParam(static_cast< ClpDblParam >(key), value);
1473	} else {
1474	return false;
1475	}
1476	}
1477
1478	//-----------------------------------------------------------------------------
1479
1480	bool OsiClpSolverInterface::setStrParam(OsiStrParam key, const std::string &value)
1481	{
1482	assert(key != OsiSolverName);
1483	if (key != OsiLastStrParam) {
1484	return modelPtr_->setStrParam(static_cast< ClpStrParam >(key), value);
1485	} else {
1486	return false;
1487	}
1488	}
1489
1490	//-----------------------------------------------------------------------------
1491
1492	bool OsiClpSolverInterface::getIntParam(OsiIntParam key, int &value) const
1493	{
1494	return modelPtr_->getIntParam(static_cast< ClpIntParam >(key), value);
1495	}
1496
1497	//-----------------------------------------------------------------------------
1498
1499	bool OsiClpSolverInterface::getDblParam(OsiDblParam key, double &value) const
1500	{
1501	if (key != OsiLastDblParam) {
1502	bool condition = modelPtr_->getDblParam(static_cast< ClpDblParam >(key), value);
1503	if (key == OsiDualObjectiveLimit \|\| key == OsiPrimalObjectiveLimit)
1504	value *= modelPtr_->optimizationDirection();
1505	return condition;
1506	} else {
1507	return false;
1508	}
1509	}
1510
1511	//-----------------------------------------------------------------------------
1512
1513	bool OsiClpSolverInterface::getStrParam(OsiStrParam key, std::string &value) const
1514	{
1515	if (key == OsiSolverName) {
1516	value = "clp";
1517	return true;
1518	}
1519	if (key != OsiLastStrParam) {
1520	return modelPtr_->getStrParam(static_cast< ClpStrParam >(key), value);
1521	} else {
1522	return false;
1523	}
1524	}
1525
1526	//#############################################################################
1527	// Methods returning info on how the solution process terminated
1528	//#############################################################################
1529
1530	bool OsiClpSolverInterface::isAbandoned() const
1531	{
1532	// not sure about -1 (should not happen)
1533	return (modelPtr_->status() == 4 \|\| modelPtr_->status() == -1 \|\| (modelPtr_->status() == 1 && modelPtr_->secondaryStatus() == 8));
1534	}
1535
1536	bool OsiClpSolverInterface::isProvenOptimal() const
1537	{
1538
1539	const int stat = modelPtr_->status();
1540	return (stat == 0);
1541	}
1542
1543	bool OsiClpSolverInterface::isProvenPrimalInfeasible() const
1544	{
1545
1546	const int stat = modelPtr_->status();
1547	if (stat != 1)
1548	return false;
1549	return true;
1550	}
1551
1552	bool OsiClpSolverInterface::isProvenDualInfeasible() const
1553	{
1554	const int stat = modelPtr_->status();
1555	return stat == 2;
1556	}
1557	/*
1558	NOTE - Coding if limit > 1.0e30 says that 1.0e29 is loose bound
1559	so all maximization tests are changed
1560	*/
1561	bool OsiClpSolverInterface::isPrimalObjectiveLimitReached() const
1562	{
1563	double limit = 0.0;
1564	modelPtr_->getDblParam(ClpPrimalObjectiveLimit, limit);
1565	if (fabs(limit) > 1e30) {
1566	// was not ever set
1567	return false;
1568	}
1569
1570	const double obj = modelPtr_->objectiveValue();
1571	int maxmin = static_cast< int >(modelPtr_->optimizationDirection());
1572
1573	switch (lastAlgorithm_) {
1574	case 0: // no simplex was needed
1575	return maxmin > 0 ? (obj < limit) /minim/ : (-obj < limit) /maxim/;
1576	case 2: // dual simplex
1577	if (modelPtr_->status() == 0) // optimal
1578	return maxmin > 0 ? (obj < limit) /minim/ : (-obj < limit) /maxim/;
1579	return false;
1580	case 1: // primal simplex
1581	return maxmin > 0 ? (obj < limit) /minim/ : (-obj < limit) /maxim/;
1582	}
1583	return false; // fake return
1584	}
1585
1586	bool OsiClpSolverInterface::isDualObjectiveLimitReached() const
1587	{
1588
1589	const int stat = modelPtr_->status();
1590	if (stat == 1)
1591	return true;
1592	else if (stat < 0)
1593	return false;
1594	double limit = 0.0;
1595	modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
1596	if (fabs(limit) > 1e30) {
1597	// was not ever set
1598	return false;
1599	}
1600
1601	const double obj = modelPtr_->objectiveValue();
1602	int maxmin = static_cast< int >(modelPtr_->optimizationDirection());
1603
1604	switch (lastAlgorithm_) {
1605	case 0: // no simplex was needed
1606	return maxmin > 0 ? (obj > limit) /minim/ : (-obj > limit) /maxim/;
1607	case 1: // primal simplex
1608	if (stat == 0) // optimal
1609	return maxmin > 0 ? (obj > limit) /minim/ : (-obj > limit) /maxim/;
1610	return false;
1611	case 2: // dual simplex
1612	if (stat != 0 && stat != 3)
1613	// over dual limit
1614	return true;
1615	return maxmin > 0 ? (obj > limit) /minim/ : (-obj > limit) /maxim/;
1616	}
1617	return false; // fake return
1618	}
1619
1620	bool OsiClpSolverInterface::isIterationLimitReached() const
1621	{
1622	const int status = modelPtr_->status();
1623	const int secondaryStatus = modelPtr_->secondaryStatus();
1624	return (status == 3 && secondaryStatus != 9);
1625	}
1626
1627	//#############################################################################
1628	// WarmStart related methods
1629	//#############################################################################
1630	CoinWarmStart *OsiClpSolverInterface::getEmptyWarmStart() const
1631	{
1632	return (static_cast< CoinWarmStart * >(new CoinWarmStartBasis()));
1633	}
1634
1635	CoinWarmStart *OsiClpSolverInterface::getWarmStart() const
1636	{
1637
1638	return new CoinWarmStartBasis(basis_);
1639	}
1640	/* Get warm start information.
1641	Return warm start information for the current state of the solver
1642	interface. If there is no valid warm start information, an empty warm
1643	start object wil be returned. This does not necessarily create an
1644	object - may just point to one. must Delete set true if user
1645	should delete returned object.
1646	OsiClp version always returns pointer and false.
1647	*/
1648	CoinWarmStart *
1649	OsiClpSolverInterface::getPointerToWarmStart(bool &mustDelete)
1650	{
1651	mustDelete = false;
1652	return &basis_;
1653	}
1654
1655	//-----------------------------------------------------------------------------
1656
1657	bool OsiClpSolverInterface::setWarmStart(const CoinWarmStart *warmstart)
1658	{
1659	modelPtr_->whatsChanged_ &= 0xffff;
1660	const CoinWarmStartBasis ws = dynamic_cast< const CoinWarmStartBasis >(warmstart);
1661	if (ws) {
1662	basis_ = CoinWarmStartBasis(*ws);
1663	return true;
1664	} else if (!warmstart) {
1665	// create from current basis
1666	basis_ = getBasis(modelPtr_);
1667	return true;
1668	} else {
1669	return false;
1670	}
1671	}
1672
1673	//#############################################################################
1674	// Hotstart related methods (primarily used in strong branching)
1675	//#############################################################################
1676	void OsiClpSolverInterface::markHotStart()
1677	{
1678	#ifdef CBC_STATISTICS
1679	osi_hot++;
1680	#endif
1681	//printf("HotStart options %x changed %x, small %x spare %x\n",
1682	// specialOptions_,modelPtr_->whatsChanged_,
1683	// smallModel_,spareArrays_);
1684	modelPtr_->setProblemStatus(0);
1685	saveData_.perturbation_ = 0;
1686	saveData_.specialOptions_ = modelPtr_->specialOptions_;
1687	modelPtr_->specialOptions_ \|= 0x1000000;
1688	modelPtr_->specialOptions_ = saveData_.specialOptions_;
1689	ClpObjective *savedObjective = NULL;
1690	double savedDualLimit = modelPtr_->dblParam_[ClpDualObjectiveLimit];
1691	if (fakeObjective_) {
1692	modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions() & (~128));
1693	// See if all with costs fixed
1694	int numberColumns = modelPtr_->numberColumns_;
1695	const double *obj = modelPtr_->objective();
1696	const double *lower = modelPtr_->columnLower();
1697	const double *upper = modelPtr_->columnUpper();
1698	int i;
1699	for (i = 0; i < numberColumns; i++) {
1700	double objValue = obj[i];
1701	if (objValue) {
1702	if (lower[i] != upper[i])
1703	break;
1704	}
1705	}
1706	if (i == numberColumns) {
1707	if ((specialOptions_ & 524288) == 0) {
1708	// Set fake
1709	savedObjective = modelPtr_->objective_;
1710	modelPtr_->objective_ = fakeObjective_;
1711	modelPtr_->dblParam_[ClpDualObjectiveLimit] = COIN_DBL_MAX;
1712	saveData_.perturbation_ = 1;
1713	} else {
1714	modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions() \| 128);
1715	}
1716	}
1717	}
1718	#define CLEAN_HOT_START
1719	#ifdef CLEAN_HOT_START
1720	if ((specialOptions_ & 65536) != 0) {
1721	//specialOptions_ \|= 65536;
1722	saveData_.scalingFlag_ = modelPtr_->logLevel();
1723	if (modelPtr_->logLevel() < 2)
1724	modelPtr_->setLogLevel(0);
1725	assert((specialOptions_ & 128) == 0);
1726	// space for save arrays
1727	int numberColumns = modelPtr_->numberColumns();
1728	int numberRows = modelPtr_->numberRows();
1729	// Get space for strong branching
1730	int size = static_cast< int >((1 + 4 * (numberRows + numberColumns)) * sizeof(double));
1731	// and for save of original column bounds
1732	size += static_cast< int >(2 * numberColumns * sizeof(double));
1733	size += static_cast< int >((1 + 4 * numberRows + 2 * numberColumns) * sizeof(int));
1734	size += numberRows + numberColumns;
1735	assert(spareArrays_ == NULL);
1736	spareArrays_ = new char[size];
1737	//memset(spareArrays_,0x20,size);
1738	// Setup for strong branching
1739	assert(factorization_ == NULL);
1740	if ((specialOptions_ & 131072) != 0) {
1741	assert(lastNumberRows_ >= 0);
1742	if (modelPtr_->rowScale_ != rowScale_.array()) {
1743	assert(modelPtr_->columnScale_ != columnScale_.array());
1744	delete[] modelPtr_->rowScale_;
1745	modelPtr_->rowScale_ = NULL;
1746	delete[] modelPtr_->columnScale_;
1747	modelPtr_->columnScale_ = NULL;
1748	if (lastNumberRows_ == modelPtr_->numberRows()) {
1749	// use scaling
1750	modelPtr_->rowScale_ = rowScale_.array();
1751	modelPtr_->columnScale_ = columnScale_.array();
1752	} else {
1753	specialOptions_ &= ~131072;
1754	}
1755	}
1756	lastNumberRows_ = -1 - lastNumberRows_;
1757	}
1758	factorization_ = static_cast< ClpSimplexDual * >(modelPtr_)->setupForStrongBranching(spareArrays_, numberRows,
1759	numberColumns, true);
1760	double arrayD = reinterpret_cast< double >(spareArrays_);
1761	arrayD[0] = modelPtr_->objectiveValue() * modelPtr_->optimizationDirection();
1762	double *saveSolution = arrayD + 1;
1763	double *saveLower = saveSolution + (numberRows + numberColumns);
1764	double *saveUpper = saveLower + (numberRows + numberColumns);
1765	double *saveObjective = saveUpper + (numberRows + numberColumns);
1766	double *saveLowerOriginal = saveObjective + (numberRows + numberColumns);
1767	double *saveUpperOriginal = saveLowerOriginal + numberColumns;
1768	CoinMemcpyN(modelPtr_->columnLower(), numberColumns, saveLowerOriginal);
1769	CoinMemcpyN(modelPtr_->columnUpper(), numberColumns, saveUpperOriginal);
1770	#if 0
1771	if (whichRange_&&whichRange_[0]) {
1772	// get ranging information
1773	int numberToDo = whichRange_[0];
1774	int * which = new int [numberToDo];
1775	// Convert column numbers
1776	int * backColumn = whichColumn+numberColumns;
1777	for (int i=0;i<numberToDo;i++) {
1778	int iColumn = whichRange_[i+1];
1779	which[i]=backColumn[iColumn];
1780	}
1781	double * downRange=new double [numberToDo];
1782	double * upRange=new double [numberToDo];
1783	int * whichDown = new int [numberToDo];
1784	int * whichUp = new int [numberToDo];
1785	modelPtr_->gutsOfSolution(NULL,NULL,false);
1786	// Tell code we can increase costs in some cases
1787	modelPtr_->setCurrentDualTolerance(0.0);
1788	((ClpSimplexOther *) modelPtr_)->dualRanging(numberToDo,which,
1789	upRange, whichUp, downRange, whichDown);
1790	delete [] whichDown;
1791	delete [] whichUp;
1792	delete [] which;
1793	rowActivity_=upRange;
1794	columnActivity_=downRange;
1795	}
1796	#endif
1797	if (savedObjective) {
1798	// fix up
1799	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
1800	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
1801	modelPtr_->objective_ = savedObjective;
1802	if (!modelPtr_->problemStatus_) {
1803	CoinZeroN(modelPtr_->dual_, modelPtr_->numberRows_);
1804	CoinZeroN(modelPtr_->reducedCost_, modelPtr_->numberColumns_);
1805	if (modelPtr_->dj_ && (modelPtr_->whatsChanged_ & 1) != 0)
1806	CoinZeroN(modelPtr_->dj_, modelPtr_->numberColumns_ + modelPtr_->numberRows_);
1807	modelPtr_->computeObjectiveValue();
1808	}
1809	}
1810	return;
1811	}
1812	#endif
1813	if ((specialOptions_ & 8192) == 0 && false) { // \|\|(specialOptions_&65536)!=0) {
1814	delete ws_;
1815	ws_ = dynamic_cast< CoinWarmStartBasis * >(getWarmStart());
1816	int numberRows = modelPtr_->numberRows();
1817	rowActivity_ = new double[numberRows];
1818	CoinMemcpyN(modelPtr_->primalRowSolution(), numberRows, rowActivity_);
1819	int numberColumns = modelPtr_->numberColumns();
1820	columnActivity_ = new double[numberColumns];
1821	CoinMemcpyN(modelPtr_->primalColumnSolution(), numberColumns, columnActivity_);
1822	} else {
1823	#if 0
1824	int saveLevel = modelPtr_->logLevel();
1825	modelPtr_->setLogLevel(0);
1826	//modelPtr_->dual();
1827	OsiClpSolverInterface::resolve();
1828	if (modelPtr_->numberIterations()>0)
1829	printf("**** iterated large %d\n",modelPtr_->numberIterations());
1830	//else
1831	//printf("no iterations\n");
1832	modelPtr_->setLogLevel(saveLevel);
1833	#endif
1834	// called from CbcNode
1835	int numberColumns = modelPtr_->numberColumns();
1836	int numberRows = modelPtr_->numberRows();
1837	// Get space for crunch and strong branching (too much)
1838	int size = static_cast< int >((1 + 4 * (numberRows + numberColumns)) * sizeof(double));
1839	// and for save of original column bounds
1840	size += static_cast< int >(2 * numberColumns * sizeof(double));
1841	size += static_cast< int >((1 + 4 * numberRows + 2 * numberColumns) * sizeof(int));
1842	size += numberRows + numberColumns;
1843	#ifdef KEEP_SMALL
1844	if (smallModel_ && (modelPtr_->whatsChanged_ & 0x30000) != 0x30000) {
1845	//printf("Bounds changed ? %x\n",modelPtr_->whatsChanged_);
1846	delete smallModel_;
1847	smallModel_ = NULL;
1848	}
1849	if (!smallModel_) {
1850	delete[] spareArrays_;
1851	spareArrays_ = NULL;
1852	}
1853	#endif
1854	if (spareArrays_ == NULL) {
1855	//if (smallModel_)
1856	//printf("small model %x\n",smallModel_);
1857	delete smallModel_;
1858	smallModel_ = NULL;
1859	spareArrays_ = new char[size];
1860	//memset(spareArrays_,0x20,size);
1861	} else {
1862	double arrayD = reinterpret_cast< double >(spareArrays_);
1863	double *saveSolution = arrayD + 1;
1864	double *saveLower = saveSolution + (numberRows + numberColumns);
1865	double *saveUpper = saveLower + (numberRows + numberColumns);
1866	double *saveObjective = saveUpper + (numberRows + numberColumns);
1867	double *saveLowerOriginal = saveObjective + (numberRows + numberColumns);
1868	double *saveUpperOriginal = saveLowerOriginal + numberColumns;
1869	double *lowerOriginal = modelPtr_->columnLower();
1870	double *upperOriginal = modelPtr_->columnUpper();
1871	arrayD = saveUpperOriginal + numberColumns;
1872	int savePivot = reinterpret_cast< int >(arrayD);
1873	int *whichRow = savePivot + numberRows;
1874	int whichColumn = whichRow + 3 numberRows;
1875	int nSame = 0;
1876	int nSub = 0;
1877	for (int i = 0; i < numberColumns; i++) {
1878	double lo = lowerOriginal[i];
1879	//char * xx = (char *) (saveLowerOriginal+i);
1880	//assert (xx[0]!=0x20\|\|xx[1]!=0x20);
1881	//xx = (char *) (saveUpperOriginal+i);
1882	//assert (xx[0]!=0x20\|\|xx[1]!=0x20);
1883	double loOld = saveLowerOriginal[i];
1884	//assert (!loOld\|\|fabs(loOld)>1.0e-30);
1885	double up = upperOriginal[i];
1886	double upOld = saveUpperOriginal[i];
1887	if (lo >= loOld && up <= upOld) {
1888	if (lo == loOld && up == upOld) {
1889	nSame++;
1890	} else {
1891	nSub++;
1892	//if (!isInteger(i))
1893	//nSub+=10;
1894	}
1895	}
1896	}
1897	//printf("Mark Hot %d bounds same, %d interior, %d bad\n",
1898	// nSame,nSub,numberColumns-nSame-nSub);
1899	if (nSame < numberColumns) {
1900	if (nSame + nSub < numberColumns) {
1901	delete smallModel_;
1902	smallModel_ = NULL;
1903	} else {
1904	// we can fix up (but should we if large number fixed?)
1905	assert(smallModel_);
1906	double *lowerSmall = smallModel_->columnLower();
1907	double *upperSmall = smallModel_->columnUpper();
1908	int numberColumns2 = smallModel_->numberColumns();
1909	for (int i = 0; i < numberColumns2; i++) {
1910	int iColumn = whichColumn[i];
1911	lowerSmall[i] = lowerOriginal[iColumn];
1912	upperSmall[i] = upperOriginal[iColumn];
1913	}
1914	}
1915	}
1916	}
1917	double arrayD = reinterpret_cast< double >(spareArrays_);
1918	arrayD[0] = modelPtr_->objectiveValue() * modelPtr_->optimizationDirection();
1919	double *saveSolution = arrayD + 1;
1920	double *saveLower = saveSolution + (numberRows + numberColumns);
1921	double *saveUpper = saveLower + (numberRows + numberColumns);
1922	double *saveObjective = saveUpper + (numberRows + numberColumns);
1923	double *saveLowerOriginal = saveObjective + (numberRows + numberColumns);
1924	double *saveUpperOriginal = saveLowerOriginal + numberColumns;
1925	arrayD = saveUpperOriginal + numberColumns;
1926	int savePivot = reinterpret_cast< int >(arrayD);
1927	int *whichRow = savePivot + numberRows;
1928	int whichColumn = whichRow + 3 numberRows;
1929	int arrayI = whichColumn + 2 numberColumns;
1930	//unsigned char * saveStatus = (unsigned char *) (arrayI+1);
1931	// Use dual region
1932	double *rhs = modelPtr_->dualRowSolution();
1933	int nBound = 0;
1934	ClpSimplex *small;
1935	#ifndef KEEP_SMALL
1936	assert(!smallModel_);
1937	small = static_cast< ClpSimplexOther * >(modelPtr_)->crunch(rhs, whichRow, whichColumn, nBound, true);
1938	bool needSolveInSetupHotStart = true;
1939	#else
1940	bool needSolveInSetupHotStart = true;
1941	if (!smallModel_) {
1942	#ifndef NDEBUG
1943	CoinFillN(whichRow, 3 * numberRows + 2 * numberColumns, -1);
1944	#endif
1945	small = static_cast< ClpSimplexOther * >(modelPtr_)->crunch(rhs, whichRow, whichColumn, nBound, true);
1946	#ifndef NDEBUG
1947	int nCopy = 3 * numberRows + 2 * numberColumns;
1948	for (int i = 0; i < nCopy; i++)
1949	assert(whichRow[i] >= -CoinMax(numberRows, numberColumns) && whichRow[i] < CoinMax(numberRows, numberColumns));
1950	#endif
1951	smallModel_ = small;
1952	//int hotIts = small->intParam_[ClpMaxNumIterationHotStart];
1953	//if (5*small->factorization_->maximumPivots()>
1954	// 4*hotIts)
1955	//small->factorization_->maximumPivots(hotIts+1);
1956	} else {
1957	assert((modelPtr_->whatsChanged_ & 0x30000) == 0x30000);
1958	//delete [] spareArrays_;
1959	//spareArrays_ = NULL;
1960	assert(spareArrays_);
1961	int nCopy = 3 * numberRows + 2 * numberColumns;
1962	nBound = whichRow[nCopy];
1963	#ifndef NDEBUG
1964	for (int i = 0; i < nCopy; i++)
1965	assert(whichRow[i] >= -CoinMax(numberRows, numberColumns) && whichRow[i] < CoinMax(numberRows, numberColumns));
1966	#endif
1967	needSolveInSetupHotStart = false;
1968	small = smallModel_;
1969	}
1970	#endif
1971	if (small) {
1972	small->specialOptions_ \|= 262144;
1973	small->specialOptions_ &= ~65536;
1974	}
1975	if (small && (specialOptions_ & 131072) != 0) {
1976	assert(lastNumberRows_ >= 0);
1977	int numberRows2 = small->numberRows();
1978	int numberColumns2 = small->numberColumns();
1979	double rowScale2 = new double[2 numberRows2];
1980	const double *rowScale = rowScale_.array();
1981	double *inverseScale2 = rowScale2 + numberRows2;
1982	const double *inverseScale = rowScale + modelPtr_->numberRows_;
1983	int i;
1984	for (i = 0; i < numberRows2; i++) {
1985	int iRow = whichRow[i];
1986	rowScale2[i] = rowScale[iRow];
1987	inverseScale2[i] = inverseScale[iRow];
1988	}
1989	small->setRowScale(rowScale2);
1990	double columnScale2 = new double[2 numberColumns2];
1991	const double *columnScale = columnScale_.array();
1992	inverseScale2 = columnScale2 + numberColumns2;
1993	inverseScale = columnScale + modelPtr_->numberColumns_;
1994	for (i = 0; i < numberColumns2; i++) {
1995	int iColumn = whichColumn[i];
1996	columnScale2[i] = columnScale[iColumn];
1997	inverseScale2[i] = inverseScale[iColumn];
1998	}
1999	small->setColumnScale(columnScale2);
2000	}
2001	if (!small) {
2002	// should never be infeasible .... but
2003	delete[] spareArrays_;
2004	spareArrays_ = NULL;
2005	delete ws_;
2006	ws_ = dynamic_cast< CoinWarmStartBasis * >(getWarmStart());
2007	int numberRows = modelPtr_->numberRows();
2008	rowActivity_ = new double[numberRows];
2009	CoinMemcpyN(modelPtr_->primalRowSolution(), numberRows, rowActivity_);
2010	int numberColumns = modelPtr_->numberColumns();
2011	columnActivity_ = new double[numberColumns];
2012	CoinMemcpyN(modelPtr_->primalColumnSolution(), numberColumns, columnActivity_);
2013	modelPtr_->setProblemStatus(1);
2014	if (savedObjective) {
2015	// fix up
2016	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
2017	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2018	modelPtr_->objective_ = savedObjective;
2019	}
2020	return;
2021	}
2022	int clpOptions = modelPtr_->specialOptions();
2023	clpOptions &= ~65536;
2024	if ((specialOptions_ & 1) == 0) {
2025	small->setSpecialOptions(clpOptions \| (64 \| 1024));
2026	} else {
2027	if ((specialOptions_ & 4) == 0)
2028	small->setSpecialOptions(clpOptions \| (64 \| 128 \| 512 \| 1024 \| 4096));
2029	else
2030	small->setSpecialOptions(clpOptions \| (64 \| 128 \| 512 \| 1024 \| 2048 \| 4096));
2031	}
2032	arrayI[0] = nBound;
2033	assert(smallModel_ == NULL \|\| small == smallModel_);
2034	if ((specialOptions_ & 256) != 0 \|\| 1) {
2035	// only need to do this on second pass in CbcNode
2036	if (modelPtr_->logLevel() < 2)
2037	small->setLogLevel(0);
2038	small->specialOptions_ \|= 262144;
2039	small->moreSpecialOptions_ = modelPtr_->moreSpecialOptions_;
2040	#define SETUP_HOT
2041	#ifndef SETUP_HOT
2042	small->dual();
2043	#else
2044	assert(factorization_ == NULL);
2045	//needSolveInSetupHotStart=true;
2046	ClpFactorization factorization = static_cast< ClpSimplexDual >(small)->setupForStrongBranching(spareArrays_,
2047	numberRows, numberColumns,
2048	needSolveInSetupHotStart);
2049	#endif
2050	if (small->numberIterations() > 0 && small->logLevel() > 2)
2051	printf("**** iterated small %d\n", small->numberIterations());
2052	//small->setLogLevel(0);
2053	// Could be infeasible if forced one way (and other way stopped on iterations)
2054	// could also be stopped on iterations
2055	if (small->status()) {
2056	#ifndef KEEP_SMALL
2057	if (small != modelPtr_)
2058	delete small;
2059	//delete smallModel_;
2060	//smallModel_=NULL;
2061	assert(!smallModel_);
2062	#else
2063	assert(small == smallModel_);
2064	if (smallModel_ != modelPtr_) {
2065	delete smallModel_;
2066	}
2067	smallModel_ = NULL;
2068	#endif
2069	delete[] spareArrays_;
2070	spareArrays_ = NULL;
2071	delete ws_;
2072	ws_ = dynamic_cast< CoinWarmStartBasis * >(getWarmStart());
2073	int numberRows = modelPtr_->numberRows();
2074	rowActivity_ = new double[numberRows];
2075	CoinMemcpyN(modelPtr_->primalRowSolution(), numberRows, rowActivity_);
2076	int numberColumns = modelPtr_->numberColumns();
2077	columnActivity_ = new double[numberColumns];
2078	CoinMemcpyN(modelPtr_->primalColumnSolution(), numberColumns, columnActivity_);
2079	modelPtr_->setProblemStatus(1);
2080	if (savedObjective) {
2081	// fix up
2082	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
2083	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2084	modelPtr_->objective_ = savedObjective;
2085	}
2086	return;
2087	} else {
2088	// update model
2089	static_cast< ClpSimplexOther * >(modelPtr_)->afterCrunch(*small, whichRow, whichColumn, nBound);
2090	}
2091	#ifndef SETUP_HOT
2092	assert(factorization_ == NULL);
2093	factorization_ = static_cast< ClpSimplexDual * >(small)->setupForStrongBranching(spareArrays_, numberRows,
2094	numberColumns, false);
2095	#else
2096	assert(factorization != NULL \|\| small->problemStatus_);
2097	factorization_ = factorization;
2098	if (factorization_ == NULL)
2099	factorization_ = static_cast< ClpSimplexDual * >(small)->setupForStrongBranching(spareArrays_, numberRows, numberColumns, false);
2100	#endif
2101	} else {
2102	assert(factorization_ == NULL);
2103	factorization_ = static_cast< ClpSimplexDual * >(small)->setupForStrongBranching(spareArrays_, numberRows,
2104	numberColumns, false);
2105	}
2106	smallModel_ = small;
2107	if (modelPtr_->logLevel() < 2)
2108	smallModel_->setLogLevel(0);
2109	// Setup for strong branching
2110	int numberColumns2 = smallModel_->numberColumns();
2111	CoinMemcpyN(modelPtr_->columnLower(), numberColumns, saveLowerOriginal);
2112	CoinMemcpyN(modelPtr_->columnUpper(), numberColumns, saveUpperOriginal);
2113	const double *smallLower = smallModel_->columnLower();
2114	const double *smallUpper = smallModel_->columnUpper();
2115	// But modify if bounds changed in small
2116	for (int i = 0; i < numberColumns2; i++) {
2117	int iColumn = whichColumn[i];
2118	saveLowerOriginal[iColumn] = CoinMax(saveLowerOriginal[iColumn],
2119	smallLower[i]);
2120	saveUpperOriginal[iColumn] = CoinMin(saveUpperOriginal[iColumn],
2121	smallUpper[i]);
2122	}
2123	if (whichRange_ && whichRange_[0]) {
2124	// get ranging information
2125	int numberToDo = whichRange_[0];
2126	int *which = new int[numberToDo];
2127	// Convert column numbers
2128	int *backColumn = whichColumn + numberColumns;
2129	for (int i = 0; i < numberToDo; i++) {
2130	int iColumn = whichRange_[i + 1];
2131	which[i] = backColumn[iColumn];
2132	}
2133	double *downRange = new double[numberToDo];
2134	double *upRange = new double[numberToDo];
2135	int *whichDown = new int[numberToDo];
2136	int *whichUp = new int[numberToDo];
2137	smallModel_->setFactorization(*factorization_);
2138	smallModel_->gutsOfSolution(NULL, NULL, false);
2139	// Tell code we can increase costs in some cases
2140	smallModel_->setCurrentDualTolerance(0.0);
2141	static_cast< ClpSimplexOther * >(smallModel_)->dualRanging(numberToDo, which, upRange, whichUp, downRange, whichDown);
2142	delete[] whichDown;
2143	delete[] whichUp;
2144	delete[] which;
2145	rowActivity_ = upRange;
2146	columnActivity_ = downRange;
2147	}
2148	}
2149	if (savedObjective) {
2150	// fix up
2151	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
2152	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2153	modelPtr_->objective_ = savedObjective;
2154	if (!modelPtr_->problemStatus_) {
2155	CoinZeroN(modelPtr_->dual_, modelPtr_->numberRows_);
2156	CoinZeroN(modelPtr_->reducedCost_, modelPtr_->numberColumns_);
2157	if (modelPtr_->dj_ && (modelPtr_->whatsChanged_ & 1) != 0)
2158	CoinZeroN(modelPtr_->dj_, modelPtr_->numberColumns_ + modelPtr_->numberRows_);
2159	modelPtr_->computeObjectiveValue();
2160	}
2161	}
2162	}
2163
2164	void OsiClpSolverInterface::solveFromHotStart()
2165	{
2166	#ifdef KEEP_SMALL
2167	if (!spareArrays_) {
2168	assert(!smallModel_);
2169	assert(modelPtr_->problemStatus_ == 1);
2170	return;
2171	}
2172	#endif
2173	ClpObjective *savedObjective = NULL;
2174	double savedDualLimit = modelPtr_->dblParam_[ClpDualObjectiveLimit];
2175	if (saveData_.perturbation_) {
2176	// Set fake
2177	savedObjective = modelPtr_->objective_;
2178	modelPtr_->objective_ = fakeObjective_;
2179	modelPtr_->dblParam_[ClpDualObjectiveLimit] = COIN_DBL_MAX;
2180	}
2181	int numberRows = modelPtr_->numberRows();
2182	int numberColumns = modelPtr_->numberColumns();
2183	modelPtr_->getIntParam(ClpMaxNumIteration, itlimOrig_);
2184	int itlim;
2185	modelPtr_->getIntParam(ClpMaxNumIterationHotStart, itlim);
2186	// Is there an extra copy of scaled bounds
2187	int extraCopy = (modelPtr_->maximumRows_ > 0) ? modelPtr_->maximumRows_ + modelPtr_->maximumColumns_ : 0;
2188	#ifdef CLEAN_HOT_START
2189	if ((specialOptions_ & 65536) != 0) {
2190	double arrayD = reinterpret_cast< double >(spareArrays_);
2191	double saveObjectiveValue = arrayD[0];
2192	double *saveSolution = arrayD + 1;
2193	int number = numberRows + numberColumns;
2194	CoinMemcpyN(saveSolution, number, modelPtr_->solutionRegion());
2195	double *saveLower = saveSolution + (numberRows + numberColumns);
2196	CoinMemcpyN(saveLower, number, modelPtr_->lowerRegion());
2197	double *saveUpper = saveLower + (numberRows + numberColumns);
2198	CoinMemcpyN(saveUpper, number, modelPtr_->upperRegion());
2199	double *saveObjective = saveUpper + (numberRows + numberColumns);
2200	CoinMemcpyN(saveObjective, number, modelPtr_->costRegion());
2201	double *saveLowerOriginal = saveObjective + (numberRows + numberColumns);
2202	double *saveUpperOriginal = saveLowerOriginal + numberColumns;
2203	arrayD = saveUpperOriginal + numberColumns;
2204	int savePivot = reinterpret_cast< int >(arrayD);
2205	CoinMemcpyN(savePivot, numberRows, modelPtr_->pivotVariable());
2206	int *whichRow = savePivot + numberRows;
2207	int whichColumn = whichRow + 3 numberRows;
2208	int arrayI = whichColumn + 2 numberColumns;
2209	unsigned char saveStatus = reinterpret_cast< unsigned char >(arrayI + 1);
2210	CoinMemcpyN(saveStatus, number, modelPtr_->statusArray());
2211	modelPtr_->setFactorization(*factorization_);
2212	double *columnLower = modelPtr_->columnLower();
2213	double *columnUpper = modelPtr_->columnUpper();
2214	// make sure whatsChanged_ has 1 set
2215	modelPtr_->setWhatsChanged(511);
2216	double *lowerInternal = modelPtr_->lowerRegion();
2217	double *upperInternal = modelPtr_->upperRegion();
2218	double rhsScale = modelPtr_->rhsScale();
2219	const double *columnScale = NULL;
2220	if (modelPtr_->scalingFlag() > 0)
2221	columnScale = modelPtr_->columnScale();
2222	// and do bounds in case dual needs them
2223	int iColumn;
2224	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
2225	if (columnLower[iColumn] > saveLowerOriginal[iColumn]) {
2226	double value = columnLower[iColumn];
2227	value *= rhsScale;
2228	if (columnScale)
2229	value /= columnScale[iColumn];
2230	lowerInternal[iColumn] = value;
2231	if (extraCopy)
2232	lowerInternal[iColumn + extraCopy] = value;
2233	}
2234	if (columnUpper[iColumn] < saveUpperOriginal[iColumn]) {
2235	double value = columnUpper[iColumn];
2236	value *= rhsScale;
2237	if (columnScale)
2238	value /= columnScale[iColumn];
2239	upperInternal[iColumn] = value;
2240	if (extraCopy)
2241	upperInternal[iColumn + extraCopy] = value;
2242	}
2243	}
2244	// Start of fast iterations
2245	bool alwaysFinish = ((specialOptions_ & 32) == 0) ? true : false;
2246	//modelPtr_->setLogLevel(1);
2247	modelPtr_->setIntParam(ClpMaxNumIteration, itlim);
2248	int saveNumberFake = (static_cast< ClpSimplexDual * >(modelPtr_))->numberFake_;
2249	int status = (static_cast< ClpSimplexDual * >(modelPtr_))->fastDual(alwaysFinish);
2250	(static_cast< ClpSimplexDual * >(modelPtr_))->numberFake_ = saveNumberFake;
2251
2252	int problemStatus = modelPtr_->problemStatus();
2253	double objectiveValue = modelPtr_->objectiveValue() * modelPtr_->optimizationDirection();
2254	CoinAssert(modelPtr_->problemStatus() \|\| modelPtr_->objectiveValue() < 1.0e50);
2255	// make sure plausible
2256	double obj = CoinMax(objectiveValue, saveObjectiveValue);
2257	if (problemStatus == 10 \|\| problemStatus < 0) {
2258	// was trying to clean up or something odd
2259	if (problemStatus == 10) {
2260	obj = saveObjectiveValue;
2261	lastAlgorithm_ = 1; // so won't fail on cutoff (in CbcNode)
2262	}
2263	status = 1;
2264	}
2265	if (status) {
2266	// not finished - might be optimal
2267	modelPtr_->checkPrimalSolution(modelPtr_->solutionRegion(0),
2268	modelPtr_->solutionRegion(1));
2269	//modelPtr_->gutsOfSolution(NULL,NULL,0);
2270	//if (problemStatus==3)
2271	//modelPtr_->computeObjectiveValue();
2272	objectiveValue = modelPtr_->objectiveValue() * modelPtr_->optimizationDirection();
2273	obj = CoinMax(objectiveValue, saveObjectiveValue);
2274	if (!modelPtr_->numberDualInfeasibilities()) {
2275	double limit = 0.0;
2276	modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2277	if (modelPtr_->secondaryStatus() == 1 && !problemStatus && obj < limit) {
2278	obj = limit;
2279	problemStatus = 3;
2280	}
2281	if (!modelPtr_->numberPrimalInfeasibilities() && obj < limit) {
2282	problemStatus = 0;
2283	} else if (problemStatus == 10) {
2284	problemStatus = 3;
2285	obj = saveObjectiveValue;
2286	} else if (!modelPtr_->numberPrimalInfeasibilities()) {
2287	problemStatus = 1; // infeasible
2288	}
2289	} else {
2290	// can't say much
2291	//if (problemStatus==3)
2292	//modelPtr_->computeObjectiveValue();
2293	lastAlgorithm_ = 1; // so won't fail on cutoff (in CbcNode)
2294	obj = saveObjectiveValue;
2295	problemStatus = 3;
2296	}
2297	} else if (!problemStatus) {
2298	if (modelPtr_->isDualObjectiveLimitReached())
2299	problemStatus = 1; // infeasible
2300	}
2301	if (status && !problemStatus) {
2302	problemStatus = 3; // can't be sure
2303	lastAlgorithm_ = 1;
2304	}
2305	if (problemStatus < 0)
2306	problemStatus = 3;
2307	modelPtr_->setProblemStatus(problemStatus);
2308	modelPtr_->setObjectiveValue(obj * modelPtr_->optimizationDirection());
2309	double *solution = modelPtr_->primalColumnSolution();
2310	const double *solution2 = modelPtr_->solutionRegion();
2311	// could just do changed bounds - also try double size scale so can use * not /
2312	if (!columnScale) {
2313	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
2314	solution[iColumn] = solution2[iColumn];
2315	}
2316	} else {
2317	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
2318	solution[iColumn] = solution2[iColumn] * columnScale[iColumn];
2319	}
2320	}
2321	CoinMemcpyN(saveLowerOriginal, numberColumns, columnLower);
2322	CoinMemcpyN(saveUpperOriginal, numberColumns, columnUpper);
2323	#if 0
2324	// could combine with loop above
2325	if (modelPtr_==modelPtr_)
2326	modelPtr_->computeObjectiveValue();
2327	if (status&&!problemStatus) {
2328	memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2329	modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2330	modelPtr_->checkSolutionInternal();
2331	//modelPtr_->setLogLevel(1);
2332	//modelPtr_->allSlackBasis();
2333	//modelPtr_->primal(1);
2334	//memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2335	//modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2336	//modelPtr_->checkSolutionInternal();
2337	assert (!modelPtr_->problemStatus());
2338	}
2339	#endif
2340	// and back bounds
2341	CoinMemcpyN(saveLower, number, modelPtr_->lowerRegion());
2342	CoinMemcpyN(saveUpper, number, modelPtr_->upperRegion());
2343	if (extraCopy) {
2344	CoinMemcpyN(saveLower, number, modelPtr_->lowerRegion() + extraCopy);
2345	CoinMemcpyN(saveUpper, number, modelPtr_->upperRegion() + extraCopy);
2346	}
2347	modelPtr_->setIntParam(ClpMaxNumIteration, itlimOrig_);
2348	if (savedObjective) {
2349	// fix up
2350	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
2351	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2352	modelPtr_->objective_ = savedObjective;
2353	if (!modelPtr_->problemStatus_) {
2354	CoinZeroN(modelPtr_->dual_, modelPtr_->numberRows_);
2355	CoinZeroN(modelPtr_->reducedCost_, modelPtr_->numberColumns_);
2356	if (modelPtr_->dj_ && (modelPtr_->whatsChanged_ & 1) != 0)
2357	CoinZeroN(modelPtr_->dj_, modelPtr_->numberColumns_ + modelPtr_->numberRows_);
2358	modelPtr_->computeObjectiveValue();
2359	}
2360	}
2361	return;
2362	}
2363	#endif
2364	if (smallModel_ == NULL) {
2365	setWarmStart(ws_);
2366	CoinMemcpyN(rowActivity_, numberRows, modelPtr_->primalRowSolution());
2367	CoinMemcpyN(columnActivity_, numberColumns, modelPtr_->primalColumnSolution());
2368	modelPtr_->setIntParam(ClpMaxNumIteration, CoinMin(itlim, 9999));
2369	resolve();
2370	} else {
2371	assert(spareArrays_);
2372	double arrayD = reinterpret_cast< double >(spareArrays_);
2373	double saveObjectiveValue = arrayD[0];
2374	double *saveSolution = arrayD + 1;
2375	// double check arrays exist (? for nonlinear)
2376	//if (!smallModel_->solutionRegion())
2377	//smallModel_->createRim(63);
2378	int numberRows2 = smallModel_->numberRows();
2379	int numberColumns2 = smallModel_->numberColumns();
2380	int number = numberRows2 + numberColumns2;
2381	CoinMemcpyN(saveSolution, number, smallModel_->solutionRegion());
2382	double *saveLower = saveSolution + (numberRows + numberColumns);
2383	CoinMemcpyN(saveLower, number, smallModel_->lowerRegion());
2384	double *saveUpper = saveLower + (numberRows + numberColumns);
2385	CoinMemcpyN(saveUpper, number, smallModel_->upperRegion());
2386	double *saveObjective = saveUpper + (numberRows + numberColumns);
2387	CoinMemcpyN(saveObjective, number, smallModel_->costRegion());
2388	double *saveLowerOriginal = saveObjective + (numberRows + numberColumns);
2389	double *saveUpperOriginal = saveLowerOriginal + numberColumns;
2390	arrayD = saveUpperOriginal + numberColumns;
2391	int savePivot = reinterpret_cast< int >(arrayD);
2392	CoinMemcpyN(savePivot, numberRows2, smallModel_->pivotVariable());
2393	int *whichRow = savePivot + numberRows;
2394	int whichColumn = whichRow + 3 numberRows;
2395	int arrayI = whichColumn + 2 numberColumns;
2396	unsigned char saveStatus = reinterpret_cast< unsigned char >(arrayI + 1);
2397	CoinMemcpyN(saveStatus, number, smallModel_->statusArray());
2398	/* If factorization_ NULL then infeasible
2399	not really sure how could have slipped through.
2400	But this can't make situation worse */
2401	if (factorization_)
2402	smallModel_->setFactorization(*factorization_);
2403	//int * backColumn = whichColumn+numberColumns;
2404	const double *lowerBig = modelPtr_->columnLower();
2405	const double *upperBig = modelPtr_->columnUpper();
2406	// make sure whatsChanged_ has 1 set
2407	smallModel_->setWhatsChanged(511);
2408	double *lowerSmall = smallModel_->lowerRegion();
2409	double *upperSmall = smallModel_->upperRegion();
2410	double *solutionSmall = smallModel_->solutionRegion();
2411	double *lowerSmallReal = smallModel_->columnLower();
2412	double *upperSmallReal = smallModel_->columnUpper();
2413	int i;
2414	double rhsScale = smallModel_->rhsScale();
2415	const double *columnScale = NULL;
2416	const double *rowScale = NULL;
2417	if (smallModel_->scalingFlag() > 0) {
2418	columnScale = smallModel_->columnScale();
2419	rowScale = smallModel_->rowScale();
2420	}
2421	// and do bounds in case dual needs them
2422	// may be infeasible
2423	for (i = 0; i < numberColumns2; i++) {
2424	int iColumn = whichColumn[i];
2425	if (lowerBig[iColumn] > saveLowerOriginal[iColumn]) {
2426	double value = lowerBig[iColumn];
2427	lowerSmallReal[i] = value;
2428	value *= rhsScale;
2429	if (columnScale)
2430	value /= columnScale[i];
2431	lowerSmall[i] = value;
2432	if (smallModel_->getColumnStatus(i) == ClpSimplex::atLowerBound)
2433	solutionSmall[i] = value;
2434	}
2435	if (upperBig[iColumn] < saveUpperOriginal[iColumn]) {
2436	double value = upperBig[iColumn];
2437	upperSmallReal[i] = value;
2438	value *= rhsScale;
2439	if (columnScale)
2440	value /= columnScale[i];
2441	upperSmall[i] = value;
2442	if (smallModel_->getColumnStatus(i) == ClpSimplex::atUpperBound)
2443	solutionSmall[i] = value;
2444	}
2445	if (upperSmall[i] < lowerSmall[i] - 1.0e-8)
2446	break;
2447	}
2448	/* If factorization_ NULL then infeasible
2449	not really sure how could have slipped through.
2450	But this can't make situation worse */
2451	bool infeasible = (i < numberColumns2 \|\| !factorization_);
2452	// Start of fast iterations
2453	bool alwaysFinish = ((specialOptions_ & 32) == 0) ? true : false;
2454	//smallModel_->setLogLevel(1);
2455	smallModel_->setIntParam(ClpMaxNumIteration, itlim);
2456	int saveNumberFake = (static_cast< ClpSimplexDual * >(smallModel_))->numberFake_;
2457	int status;
2458	if (!infeasible) {
2459	if ((modelPtr_->specialOptions() & 0x011200000) == 0x11200000) {
2460	smallModel_->specialOptions_ \|= 2097152;
2461	//smallModel_->specialOptions_&=~2097152;
2462	delete[] modelPtr_->ray_;
2463	delete[] smallModel_->ray_;
2464	modelPtr_->ray_ = NULL;
2465	smallModel_->ray_ = NULL;
2466	}
2467	status = static_cast< ClpSimplexDual * >(smallModel_)->fastDual(alwaysFinish);
2468	if ((modelPtr_->specialOptions() & 0x011200000) == 0x11200000 && smallModel_->ray_) {
2469	if (smallModel_->sequenceOut_ < smallModel_->numberColumns_)
2470	modelPtr_->sequenceOut_ = whichColumn[smallModel_->sequenceOut_];
2471	else
2472	modelPtr_->sequenceOut_ = whichRow[smallModel_->sequenceOut_ - smallModel_->numberColumns_] + modelPtr_->numberColumns_;
2473	if (smallModel_->rowScale_) {
2474	// scale ray
2475	double scaleFactor = 1.0;
2476	if (smallModel_->sequenceOut_ < smallModel_->numberColumns_)
2477	scaleFactor = smallModel_->columnScale_[smallModel_->sequenceOut_];
2478	for (int i = 0; i < smallModel_->numberRows_; i++) {
2479	smallModel_->ray_[i] = smallModel_->rowScale_[i] scaleFactor;
2480	}
2481	}
2482	}
2483	} else {
2484	status = 0;
2485	smallModel_->setProblemStatus(1);
2486	}
2487	(static_cast< ClpSimplexDual * >(smallModel_))->numberFake_ = saveNumberFake;
2488	if (smallModel_->numberIterations() == -98) {
2489	printf("rrrrrrrrrrrr\n");
2490	smallModel_->checkPrimalSolution(smallModel_->solutionRegion(0),
2491	smallModel_->solutionRegion(1));
2492	//smallModel_->gutsOfSolution(NULL,NULL,0);
2493	//if (problemStatus==3)
2494	//smallModel_->computeObjectiveValue();
2495	printf("robj %g\n", smallModel_->objectiveValue() * modelPtr_->optimizationDirection());
2496	writeMps("rr.mps");
2497	smallModel_->writeMps("rr_small.mps");
2498	ClpSimplex temp = *smallModel_;
2499	printf("small\n");
2500	temp.setLogLevel(63);
2501	temp.dual();
2502	double limit = 0.0;
2503	modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2504	if (temp.problemStatus() == 0 && temp.objectiveValue() < limit) {
2505	printf("inf obj %g, true %g - offsets %g %g\n", smallModel_->objectiveValue(),
2506	temp.objectiveValue(),
2507	smallModel_->objectiveOffset(), temp.objectiveOffset());
2508	}
2509	printf("big\n");
2510	temp = *modelPtr_;
2511	temp.dual();
2512	if (temp.problemStatus() == 0 && temp.objectiveValue() < limit) {
2513	printf("inf obj %g, true %g - offsets %g %g\n", smallModel_->objectiveValue(),
2514	temp.objectiveValue(),
2515	smallModel_->objectiveOffset(), temp.objectiveOffset());
2516	}
2517	}
2518	int problemStatus = smallModel_->problemStatus();
2519	double objectiveValue = smallModel_->objectiveValue() * modelPtr_->optimizationDirection();
2520	CoinAssert(smallModel_->problemStatus() \|\| smallModel_->objectiveValue() < 1.0e50);
2521	// make sure plausible
2522	double obj = CoinMax(objectiveValue, saveObjectiveValue);
2523	if (problemStatus == 10 \|\| problemStatus < 0) {
2524	// was trying to clean up or something odd
2525	if (problemStatus == 10)
2526	lastAlgorithm_ = 1; // so won't fail on cutoff (in CbcNode)
2527	status = 1;
2528	}
2529	if (status) {
2530	// not finished - might be optimal
2531	smallModel_->checkPrimalSolution(smallModel_->solutionRegion(0),
2532	smallModel_->solutionRegion(1));
2533	//smallModel_->gutsOfSolution(NULL,NULL,0);
2534	//if (problemStatus==3)
2535	//smallModel_->computeObjectiveValue();
2536	objectiveValue = smallModel_->objectiveValue() * modelPtr_->optimizationDirection();
2537	if (problemStatus != 10)
2538	obj = CoinMax(objectiveValue, saveObjectiveValue);
2539	if (!smallModel_->numberDualInfeasibilities()) {
2540	double limit = 0.0;
2541	modelPtr_->getDblParam(ClpDualObjectiveLimit, limit);
2542	if (smallModel_->secondaryStatus() == 1 && !problemStatus && obj < limit) {
2543	#if 0
2544	// switch off
2545	ClpSimplex temp = *smallModel_;
2546	temp.dual();
2547	if (temp.problemStatus()==0&&temp.objectiveValue()<limit) {
2548	printf("inf obj %g, true %g - offsets %g %g\n",smallModel_->objectiveValue(),
2549	temp.objectiveValue(),
2550	smallModel_->objectiveOffset(),temp.objectiveOffset());
2551	}
2552	lastAlgorithm_=1;
2553	obj=limit;
2554	problemStatus=10;
2555	#else
2556	obj = limit;
2557	problemStatus = 3;
2558	#endif
2559	}
2560	if (!smallModel_->numberPrimalInfeasibilities() && obj < limit) {
2561	problemStatus = 0;
2562	#if 0
2563	ClpSimplex temp = *smallModel_;
2564	temp.dual();
2565	if (temp.numberIterations())
2566	printf("temp iterated\n");
2567	assert (temp.problemStatus()==0&&temp.objectiveValue()<limit);
2568	#endif
2569	} else if (problemStatus == 10) {
2570	problemStatus = 3;
2571	} else if (!smallModel_->numberPrimalInfeasibilities()) {
2572	problemStatus = 1; // infeasible
2573	}
2574	} else {
2575	// can't say much
2576	//if (problemStatus==3)
2577	//smallModel_->computeObjectiveValue();
2578	lastAlgorithm_ = 1; // so won't fail on cutoff (in CbcNode)
2579	problemStatus = 3;
2580	}
2581	} else if (!problemStatus) {
2582	if (smallModel_->isDualObjectiveLimitReached())
2583	problemStatus = 1; // infeasible
2584	}
2585	if (status && !problemStatus) {
2586	problemStatus = 3; // can't be sure
2587	lastAlgorithm_ = 1;
2588	}
2589	if (problemStatus < 0)
2590	problemStatus = 3;
2591	modelPtr_->setProblemStatus(problemStatus);
2592	modelPtr_->setObjectiveValue(obj * modelPtr_->optimizationDirection());
2593	modelPtr_->setSumDualInfeasibilities(smallModel_->sumDualInfeasibilities());
2594	modelPtr_->setNumberDualInfeasibilities(smallModel_->numberDualInfeasibilities());
2595	modelPtr_->setSumPrimalInfeasibilities(smallModel_->sumPrimalInfeasibilities());
2596	modelPtr_->setNumberPrimalInfeasibilities(smallModel_->numberPrimalInfeasibilities());
2597	double *solution = modelPtr_->primalColumnSolution();
2598	const double *solution2 = smallModel_->solutionRegion();
2599	double *djs = modelPtr_->dualColumnSolution();
2600	if (!columnScale) {
2601	for (i = 0; i < numberColumns2; i++) {
2602	int iColumn = whichColumn[i];
2603	solution[iColumn] = solution2[i];
2604	lowerSmallReal[i] = saveLowerOriginal[iColumn];
2605	upperSmallReal[i] = saveUpperOriginal[iColumn];
2606	}
2607	} else {
2608	for (i = 0; i < numberColumns2; i++) {
2609	int iColumn = whichColumn[i];
2610	solution[iColumn] = solution2[i] * columnScale[i];
2611	lowerSmallReal[i] = saveLowerOriginal[iColumn];
2612	upperSmallReal[i] = saveUpperOriginal[iColumn];
2613	}
2614	}
2615	// compute duals and djs
2616	double *dual = modelPtr_->dualRowSolution();
2617	const double *dual2 = smallModel_->dualRowSolution();
2618	if (!rowScale) {
2619	for (i = 0; i < numberRows2; i++) {
2620	int iRow = whichRow[i];
2621	dual[iRow] = dual2[i];
2622	}
2623	} else {
2624	for (i = 0; i < numberRows2; i++) {
2625	int iRow = whichRow[i];
2626	dual[iRow] = dual2[i] * rowScale[i];
2627	}
2628	}
2629	memcpy(djs, modelPtr_->objective(), numberColumns * sizeof(double));
2630	modelPtr_->clpMatrix()->transposeTimes(-1.0, dual, djs);
2631	// could combine with loop above
2632	if (modelPtr_ == smallModel_)
2633	modelPtr_->computeObjectiveValue();
2634	if (problemStatus == 1 && smallModel_->ray_) {
2635	delete[] modelPtr_->ray_;
2636	// get ray to full problem
2637	int numberRows = modelPtr_->numberRows();
2638	int numberRows2 = smallModel_->numberRows();
2639	int nCopy = 3 * numberRows + 2 * numberColumns;
2640	int nBound = whichRow[nCopy];
2641	double *ray = new double[numberRows];
2642	memset(ray, 0, numberRows * sizeof(double));
2643	for (int i = 0; i < numberRows2; i++) {
2644	int iRow = whichRow[i];
2645	ray[iRow] = smallModel_->ray_[i];
2646	}
2647	// Column copy of matrix
2648	const double *element = getMatrixByCol()->getElements();
2649	const int *row = getMatrixByCol()->getIndices();
2650	const CoinBigIndex *columnStart = getMatrixByCol()->getVectorStarts();
2651	const int *columnLength = getMatrixByCol()->getVectorLengths();
2652	// translate
2653	for (int jRow = nBound; jRow < 2 * numberRows; jRow++) {
2654	int iRow = whichRow[jRow];
2655	int iColumn = whichRow[jRow + numberRows];
2656	if (modelPtr_->getColumnStatus(iColumn) == ClpSimplex::basic) {
2657	double value = 0.0;
2658	double sum = 0.0;
2659	for (CoinBigIndex j = columnStart[iColumn];
2660	j < columnStart[iColumn] + columnLength[iColumn]; j++) {
2661	if (iRow == row[j]) {
2662	value = element[j];
2663	} else {
2664	sum += ray[row[j]] * element[j];
2665	}
2666	}
2667	ray[iRow] = -sum / value;
2668	}
2669	}
2670	for (int i = 0; i < modelPtr_->numberColumns_; i++) {
2671	if (modelPtr_->getStatus(i) != ClpSimplex::basic && modelPtr_->columnLower_[i] == modelPtr_->columnUpper_[i])
2672	modelPtr_->setStatus(i, ClpSimplex::isFixed);
2673	}
2674	modelPtr_->ray_ = ray;
2675	//delete [] smallModel_->ray_;
2676	//smallModel_->ray_=NULL;
2677	modelPtr_->directionOut_ = smallModel_->directionOut_;
2678	lastAlgorithm_ = 2; // dual
2679	}
2680	#if 1
2681	if (status && !problemStatus) {
2682	memset(modelPtr_->primalRowSolution(), 0, numberRows * sizeof(double));
2683	modelPtr_->clpMatrix()->times(1.0, solution, modelPtr_->primalRowSolution());
2684	modelPtr_->checkSolutionInternal();
2685	//modelPtr_->setLogLevel(1);
2686	//modelPtr_->allSlackBasis();
2687	//modelPtr_->primal(1);
2688	//memset(modelPtr_->primalRowSolution(),0,numberRows*sizeof(double));
2689	//modelPtr_->clpMatrix()->times(1.0,solution,modelPtr_->primalRowSolution());
2690	//modelPtr_->checkSolutionInternal();
2691	assert(!modelPtr_->problemStatus());
2692	}
2693	#endif
2694	modelPtr_->setNumberIterations(smallModel_->numberIterations());
2695	// and back bounds
2696	CoinMemcpyN(saveLower, number, smallModel_->lowerRegion());
2697	CoinMemcpyN(saveUpper, number, smallModel_->upperRegion());
2698	}
2699	if (savedObjective) {
2700	// fix up
2701	modelPtr_->dblParam_[ClpDualObjectiveLimit] = savedDualLimit;
2702	//modelPtr_->setMoreSpecialOptions(modelPtr_->moreSpecialOptions()&(~32));
2703	modelPtr_->objective_ = savedObjective;
2704	if (!modelPtr_->problemStatus_) {
2705	CoinZeroN(modelPtr_->dual_, modelPtr_->numberRows_);
2706	CoinZeroN(modelPtr_->reducedCost_, modelPtr_->numberColumns_);
2707	if (modelPtr_->dj_ && (modelPtr_->whatsChanged_ & 1) != 0)
2708	CoinZeroN(modelPtr_->dj_, modelPtr_->numberColumns_ + modelPtr_->numberRows_);
2709	modelPtr_->computeObjectiveValue();
2710	}
2711	}
2712	modelPtr_->setIntParam(ClpMaxNumIteration, itlimOrig_);
2713	}
2714
2715	void OsiClpSolverInterface::unmarkHotStart()
2716	{
2717	#ifdef CLEAN_HOT_START
2718	if ((specialOptions_ & 65536) != 0) {
2719	modelPtr_->setLogLevel(saveData_.scalingFlag_);
2720	modelPtr_->deleteRim(0);
2721	if (lastNumberRows_ < 0) {
2722	specialOptions_ \|= 131072;
2723	lastNumberRows_ = -1 - lastNumberRows_;
2724	if (modelPtr_->rowScale_) {
2725	if (modelPtr_->rowScale_ != rowScale_.array()) {
2726	delete[] modelPtr_->rowScale_;
2727	delete[] modelPtr_->columnScale_;
2728	}
2729	modelPtr_->rowScale_ = NULL;
2730	modelPtr_->columnScale_ = NULL;
2731	}
2732	}
2733	delete factorization_;
2734	delete[] spareArrays_;
2735	smallModel_ = NULL;
2736	spareArrays_ = NULL;
2737	factorization_ = NULL;
2738	delete[] rowActivity_;
2739	delete[] columnActivity_;
2740	rowActivity_ = NULL;
2741	columnActivity_ = NULL;
2742	return;
2743	}
2744	#endif
2745	if (smallModel_ == NULL) {
2746	setWarmStart(ws_);
2747	int numberRows = modelPtr_->numberRows();
2748	int numberColumns = modelPtr_->numberColumns();
2749	CoinMemcpyN(rowActivity_, numberRows, modelPtr_->primalRowSolution());
2750	CoinMemcpyN(columnActivity_, numberColumns, modelPtr_->primalColumnSolution());
2751	delete ws_;
2752	ws_ = NULL;
2753	} else {
2754	#ifndef KEEP_SMALL
2755	if (smallModel_ != modelPtr_)
2756	delete smallModel_;
2757	smallModel_ = NULL;
2758	#else
2759	if (smallModel_ == modelPtr_) {
2760	smallModel_ = NULL;
2761	} else if (smallModel_) {
2762	if (!spareArrays_) {
2763	delete smallModel_;
2764	smallModel_ = NULL;
2765	delete factorization_;
2766	factorization_ = NULL;
2767	} else {
2768	static_cast< ClpSimplexDual * >(smallModel_)->cleanupAfterStrongBranching(factorization_);
2769	if ((smallModel_->specialOptions_ & 4096) == 0) {
2770	delete factorization_;
2771	}
2772	}
2773	}
2774	#endif
2775	//delete [] spareArrays_;
2776	//spareArrays_=NULL;
2777	factorization_ = NULL;
2778	}
2779	delete[] rowActivity_;
2780	delete[] columnActivity_;
2781	rowActivity_ = NULL;
2782	columnActivity_ = NULL;
2783	// Make sure whatsChanged not out of sync
2784	if (!modelPtr_->columnUpperWork_)
2785	modelPtr_->whatsChanged_ &= ~0xffff;
2786	modelPtr_->specialOptions_ = saveData_.specialOptions_;
2787	}
2788
2789	#ifdef CONFLICT_CUTS
2790	// Return a conflict analysis cut from small model
2791	OsiRowCut *
2792	OsiClpSolverInterface::smallModelCut(const double originalLower, const double originalUpper,
2793	int numberRowsAtContinuous, const int *whichGenerator,
2794	int typeCut)
2795	{
2796	if (smallModel_ && smallModel_->ray_) {
2797	//printf("Could do small cut\n");
2798	int numberRows = modelPtr_->numberRows();
2799	int numberRows2 = smallModel_->numberRows();
2800	int numberColumns = modelPtr_->numberColumns();
2801	int numberColumns2 = smallModel_->numberColumns();
2802	double arrayD = reinterpret_cast< double >(spareArrays_);
2803	double *saveSolution = arrayD + 1;
2804	double *saveLower = saveSolution + (numberRows + numberColumns);
2805	double *saveUpper = saveLower + (numberRows + numberColumns);
2806	double *saveObjective = saveUpper + (numberRows + numberColumns);
2807	double *saveLowerOriginal = saveObjective + (numberRows + numberColumns);
2808	double *saveUpperOriginal = saveLowerOriginal + numberColumns;
2809	//double * lowerOriginal = modelPtr_->columnLower();
2810	//double * upperOriginal = modelPtr_->columnUpper();
2811	int savePivot = reinterpret_cast< int >(saveUpperOriginal + numberColumns);
2812	int *whichRow = savePivot + numberRows;
2813	int whichColumn = whichRow + 3 numberRows;
2814	int nCopy = 3 * numberRows + 2 * numberColumns;
2815	int nBound = whichRow[nCopy];
2816	if (smallModel_->sequenceOut_ >= 0 && smallModel_->sequenceOut_ < numberColumns2)
2817	modelPtr_->sequenceOut_ = whichColumn[smallModel_->sequenceOut_];
2818	else
2819	modelPtr_->sequenceOut_ = modelPtr_->numberColumns_ + whichRow[smallModel_->sequenceOut_];
2820	unsigned char *saveStatus = CoinCopyOfArray(modelPtr_->status_,
2821	numberRows + numberColumns);
2822	// get ray to full problem
2823	for (int i = 0; i < numberColumns2; i++) {
2824	int iColumn = whichColumn[i];
2825	modelPtr_->setStatus(iColumn, smallModel_->getStatus(i));
2826	}
2827	double *ray = new double[numberRows + numberColumns2 + numberColumns];
2828	char *mark = new char[numberRows];
2829	memset(ray, 0, (numberRows + numberColumns2 + numberColumns) * sizeof(double));
2830	double *smallFarkas = ray + numberRows;
2831	double *farkas = smallFarkas + numberColumns2;
2832	double *saveScale = smallModel_->rowScale_;
2833	smallModel_->rowScale_ = NULL;
2834	smallModel_->transposeTimes(1.0, smallModel_->ray_, smallFarkas);
2835	smallModel_->rowScale_ = saveScale;
2836	for (int i = 0; i < numberColumns2; i++)
2837	farkas[whichColumn[i]] = smallFarkas[i];
2838	memset(mark, 0, numberRows);
2839	for (int i = 0; i < numberRows2; i++) {
2840	int iRow = whichRow[i];
2841	modelPtr_->setRowStatus(iRow, smallModel_->getRowStatus(i));
2842	ray[iRow] = smallModel_->ray_[i];
2843	mark[iRow] = 1;
2844	}
2845	// Column copy of matrix
2846	const double *element = getMatrixByCol()->getElements();
2847	const int *row = getMatrixByCol()->getIndices();
2848	const CoinBigIndex *columnStart = getMatrixByCol()->getVectorStarts();
2849	const int *columnLength = getMatrixByCol()->getVectorLengths();
2850	// pick up small pivot row
2851	int pivotRow = smallModel_->spareIntArray_[3];
2852	// translate
2853	if (pivotRow >= 0)
2854	pivotRow = whichRow[pivotRow];
2855	modelPtr_->spareIntArray_[3] = pivotRow;
2856	for (int jRow = nBound; jRow < 2 * numberRows; jRow++) {
2857	int iRow = whichRow[jRow];
2858	int iColumn = whichRow[jRow + numberRows];
2859	if (modelPtr_->getColumnStatus(iColumn) == ClpSimplex::basic) {
2860	double value = 0.0;
2861	double sum = 0.0;
2862	for (CoinBigIndex j = columnStart[iColumn];
2863	j < columnStart[iColumn] + columnLength[iColumn]; j++) {
2864	if (iRow == row[j]) {
2865	value = element[j];
2866	} else if (mark[row[j]]) {
2867	sum += ray[row[j]] * element[j];
2868	}
2869	}
2870	double target = farkas[iColumn];
2871	if (iRow != pivotRow) {
2872	ray[iRow] = (target - sum) / value;
2873	} else {
2874	printf("what now - direction %d wanted %g sum %g value %g\n",
2875	smallModel_->directionOut_, ray[iRow],
2876	sum, value);
2877	}
2878	mark[iRow] = 1;
2879	}
2880	}
2881	delete[] mark;
2882	for (int i = 0; i < modelPtr_->numberColumns_; i++) {
2883	if (modelPtr_->getStatus(i) != ClpSimplex::basic && modelPtr_->columnLower_[i] == modelPtr_->columnUpper_[i])
2884	modelPtr_->setStatus(i, ClpSimplex::isFixed);
2885	}
2886	#if 0
2887	{
2888	int nBad=0;
2889	double * fBig = new double [2numberColumns+2numberRows];
2890	double * fSmall = fBig+numberColumns;
2891	double * effectiveRhs = fSmall+numberColumns;
2892	double * effectiveRhs2 = effectiveRhs+numberRows;
2893	memset(fBig,0,2numberColumnssizeof(double));
2894	{
2895	memcpy(fBig,modelPtr_->columnLower_,numberColumns*sizeof(double));
2896	const double * columnLower = modelPtr_->columnLower_;
2897	const double * columnUpper = modelPtr_->columnUpper_;
2898	for (int i=0;i<numberColumns2;i++) {
2899	int iBig=whichColumn[i];
2900	if (smallModel_->getStatus(i)==ClpSimplex::atLowerBound\|\|
2901	smallModel_->getStatus(i)==ClpSimplex::isFixed\|\|
2902	columnLower[iBig]==columnUpper[iBig]) {
2903	fSmall[i]=columnLower[iBig];
2904	} else if (smallModel_->getStatus(i)==ClpSimplex::atUpperBound) {
2905	fSmall[i]=columnUpper[iBig];
2906	} else if (i==smallModel_->sequenceOut_) {
2907	if (smallModel_->directionOut_<0) {
2908	// above upper bound
2909	fSmall[i]=columnUpper[iBig];
2910	} else {
2911	// below lower bound
2912	fSmall[i]=columnLower[iBig];
2913	}
2914	} else if (smallModel_->columnLower_[i]==smallModel_->columnUpper_[i]) {
2915	fSmall[i]=smallModel_->columnLower_[i];
2916	}
2917	fBig[whichColumn[i]]=fSmall[i];
2918	}
2919	{
2920	// only works for one test case
2921	memcpy(effectiveRhs2,modelPtr_->rowUpper_,numberRows*sizeof(double));
2922	double * saveScale = modelPtr_->rowScale_;
2923	ClpPackedMatrix * saveMatrix = modelPtr_->scaledMatrix_;
2924	modelPtr_->rowScale_=NULL;
2925	modelPtr_->scaledMatrix_=NULL;
2926	modelPtr_->times(-1.0,fBig,effectiveRhs2);
2927	modelPtr_->scaledMatrix_=saveMatrix;
2928	modelPtr_->rowScale_=saveScale;
2929	memset(fBig,0,numberColumns*sizeof(double));
2930	}
2931	const double * rowLower = smallModel_->rowLower_;
2932	const double * rowUpper = smallModel_->rowUpper_;
2933	int pivotRow = smallModel_->spareIntArray_[3];
2934	for (int i=0;i<numberRows2;i++) {
2935	if (smallModel_->getRowStatus(i)==ClpSimplex::atLowerBound\|\|
2936	rowUpper[i]==rowLower[i]\|\|
2937	smallModel_->getRowStatus(i)==ClpSimplex::isFixed) {
2938	effectiveRhs[i]=rowLower[i];
2939	//if (smallModel_->getRowStatus(i)==ClpSimplex::atLowerBound)
2940	//assert (ray[i]>-1.0e-3);
2941	} else if (smallModel_->getRowStatus(i)==ClpSimplex::atUpperBound) {
2942	effectiveRhs[i]=rowUpper[i];
2943	//assert (ray[i]<1.0e-3);
2944	} else {
2945	if (smallModel_->getRowStatus(i)!=ClpSimplex::basic) {
2946	assert (rowUpper[i]>1.0e30\|\|rowLower[i]<-1.0e30); // eventually skip
2947	if (rowUpper[i]<1.0e30)
2948	effectiveRhs[i]=rowUpper[i];
2949	else
2950	effectiveRhs[i]=rowLower[i];
2951	}
2952	}
2953	if (smallModel_->getRowStatus(i)==ClpSimplex::basic) {
2954	effectiveRhs[i]=0.0;
2955	// we should leave pivot row in and use direction for bound
2956	if (fabs(smallModel_->ray_[i])>1.0e-8) {
2957	printf("Basic small slack value %g on %d - pivotRow %d\n",smallModel_->ray_[i],i,pivotRow);
2958	if (i==pivotRow) {
2959	if (smallModel_->directionOut_<0)
2960	effectiveRhs[i]=rowUpper[i];
2961	else
2962	effectiveRhs[i]=rowLower[i];
2963	} else {
2964	//smallModel_->ray_[i]=0.0;
2965	}
2966	}
2967	}
2968	}
2969	//ClpPackedMatrix * saveMatrix = smallModel_->scaledMatrix_;
2970	double * saveScale = smallModel_->rowScale_;
2971	smallModel_->rowScale_=NULL;
2972	smallModel_->scaledMatrix_=NULL;
2973	smallModel_->times(-1.0,fSmall,effectiveRhs);
2974	double bSum=0.0;
2975	for (int i=0;i<numberRows2;i++) {
2976	bSum += effectiveRhs[i]*smallModel_->ray_[i];
2977	}
2978	printf("Small bsum %g\n",bSum);
2979	smallModel_->rowScale_=saveScale;
2980	}
2981	double * saveScale = smallModel_->rowScale_;
2982	smallModel_->rowScale_=NULL;
2983	memset(fSmall,0,numberColumns*sizeof(double));
2984	smallModel_->transposeTimes(1.0,smallModel_->ray_,fSmall);
2985	smallModel_->rowScale_=saveScale;
2986	ClpPackedMatrix * saveMatrix = modelPtr_->scaledMatrix_;
2987	saveScale = modelPtr_->rowScale_;
2988	modelPtr_->rowScale_=NULL;
2989	modelPtr_->scaledMatrix_=NULL;
2990	modelPtr_->transposeTimes(1.0,ray,fBig);
2991	modelPtr_->scaledMatrix_=saveMatrix;
2992	modelPtr_->rowScale_=saveScale;
2993	for (int i=0;i<numberColumns2;i++) {
2994	int iBig=whichColumn[i];
2995	if (fabs(fBig[iBig]-fSmall[i])>1.0e-4) {
2996	printf("small %d %g big %d %g\n",
2997	i,fSmall[i],iBig,fBig[iBig]);
2998	nBad++;
2999	}
3000	}
3001	delete [] fBig;
3002	if (nBad)
3003	printf("Bad %d odd %d\n",nBad,2*numberRows-nBound);
3004	}
3005	#endif
3006	modelPtr_->ray_ = ray;
3007	lastAlgorithm_ = 2;
3008	modelPtr_->directionOut_ = smallModel_->directionOut_;
3009	OsiRowCut *cut = modelCut(originalLower, originalUpper,
3010	numberRowsAtContinuous, whichGenerator, typeCut);
3011	smallModel_->deleteRay();
3012	memcpy(modelPtr_->status_, saveStatus, numberRows + numberColumns);
3013	delete[] saveStatus;
3014	if (cut) {
3015	//printf("got small cut\n");
3016	//delete cut;
3017	//cut=NULL;
3018	}
3019	return cut;
3020	} else {
3021	return NULL;
3022	}
3023	}
3024	static int debugMode = 0;
3025	// Return a conflict analysis cut from model
3026	// If type is 0 then genuine cut, if 1 then only partially processed
3027	OsiRowCut *
3028	OsiClpSolverInterface::modelCut(const double originalLower, const double originalUpper,
3029	int numberRowsAtContinuous, const int *whichGenerator,
3030	int typeCut)
3031	{
3032	OsiRowCut *cut = NULL;
3033	//return NULL;
3034	if (modelPtr_->ray_) {
3035	if (lastAlgorithm_ == 2) {
3036	//printf("Could do normal cut\n");
3037	// could use existing arrays
3038	int numberRows = modelPtr_->numberRows_;
3039	int numberColumns = modelPtr_->numberColumns_;
3040	double farkas = new double[2 numberColumns + numberRows];
3041	double *bound = farkas + numberColumns;
3042	double *effectiveRhs = bound + numberColumns;
3043	/const/ double *ray = modelPtr_->ray_;
3044	// have to get rid of local cut rows
3045	whichGenerator -= numberRowsAtContinuous;
3046	int badRows = 0;
3047	for (int iRow = numberRowsAtContinuous; iRow < numberRows; iRow++) {
3048	int iType = whichGenerator[iRow];
3049	if ((iType >= 0 && iType < 20000)) {
3050	if (fabs(ray[iRow]) > 1.0e-10) {
3051	badRows++;
3052	}
3053	ray[iRow] = 0.0;
3054	}
3055	}
3056	ClpSimplex debugModel;
3057	if ((debugMode & 4) != 0)
3058	debugModel = *modelPtr_;
3059	if (badRows && (debugMode & 1) != 0)
3060	printf("%d rows from local cuts\n", badRows);
3061	// get farkas row
3062	ClpPackedMatrix *saveMatrix = modelPtr_->scaledMatrix_;
3063	double *saveScale = modelPtr_->rowScale_;
3064	modelPtr_->rowScale_ = NULL;
3065	modelPtr_->scaledMatrix_ = NULL;
3066	memset(farkas, 0, (2 * numberColumns + numberRows) * sizeof(double));
3067	modelPtr_->transposeTimes(-1.0, ray, farkas);
3068	//const char * integerInformation = modelPtr_->integerType_;
3069	//assert (integerInformation);
3070
3071	int sequenceOut = modelPtr_->sequenceOut_;
3072	// Put nonzero bounds in bound
3073	const double *columnLower = modelPtr_->columnLower_;
3074	const double *columnUpper = modelPtr_->columnUpper_;
3075	const double *columnValue = modelPtr_->columnActivity_;
3076	int numberBad = 0;
3077	int nNonzeroBasic = 0;
3078	for (int i = 0; i < numberColumns; i++) {
3079	double value = farkas[i];
3080	double boundValue = 0.0;
3081	if (modelPtr_->getStatus(i) == ClpSimplex::basic) {
3082	// treat as zero if small
3083	if (fabs(value) < 1.0e-8) {
3084	value = 0.0;
3085	farkas[i] = 0.0;
3086	}
3087	if (value) {
3088	//printf("basic %d direction %d farkas %g\n",
3089	//i,modelPtr_->directionOut_,value);
3090	nNonzeroBasic++;
3091	if (value < 0.0)
3092	boundValue = columnLower[i];
3093	else
3094	boundValue = columnUpper[i];
3095	}
3096	} else if (fabs(value) > 1.0e-10) {
3097	if (value < 0.0) {
3098	if (columnValue[i] > columnLower[i] + 1.0e-5 && value < -1.0e-7) {
3099	//printf("bad %d alpha %g not at lower\n",i,value);
3100	numberBad++;
3101	}
3102	boundValue = columnLower[i];
3103	} else {
3104	if (columnValue[i] < columnUpper[i] - 1.0e-5 && value > 1.0e-7) {
3105	//printf("bad %d alpha %g not at upper\n",i,value);
3106	numberBad++;
3107	}
3108	boundValue = columnUpper[i];
3109	}
3110	}
3111	bound[i] = boundValue;
3112	if (fabs(boundValue) > 1.0e10)
3113	numberBad++;
3114	}
3115	const double *rowLower = modelPtr_->rowLower_;
3116	const double *rowUpper = modelPtr_->rowUpper_;
3117	//int pivotRow = modelPtr_->spareIntArray_[3];
3118	//bool badPivot=pivotRow<0;
3119	for (int i = 0; i < numberRows; i++) {
3120	double value = ray[i];
3121	double rhsValue = 0.0;
3122	if (modelPtr_->getRowStatus(i) == ClpSimplex::basic) {
3123	// treat as zero if small
3124	if (fabs(value) < 1.0e-8) {
3125	value = 0.0;
3126	ray[i] = 0.0;
3127	}
3128	if (value) {
3129	//printf("row basic %d direction %d ray %g\n",
3130	// i,modelPtr_->directionOut_,value);
3131	nNonzeroBasic++;
3132	if (value < 0.0)
3133	rhsValue = rowLower[i];
3134	else
3135	rhsValue = rowUpper[i];
3136	}
3137	} else if (fabs(value) > 1.0e-10) {
3138	if (value < 0.0)
3139	rhsValue = rowLower[i];
3140	else
3141	rhsValue = rowUpper[i];
3142	}
3143	effectiveRhs[i] = rhsValue;
3144	}
3145	{
3146	double bSum = 0.0;
3147	for (int i = 0; i < numberRows; i++) {
3148	bSum += effectiveRhs[i] * ray[i];
3149	}
3150	//printf("before bounds - bSum %g\n",bSum);
3151	}
3152	modelPtr_->times(-1.0, bound, effectiveRhs);
3153	double bSum = 0.0;
3154	for (int i = 0; i < numberRows; i++) {
3155	bSum += effectiveRhs[i] * ray[i];
3156	}
3157	#if 0
3158	double bSum2=0.0;
3159	#if 1
3160	{
3161	int iSequence = modelPtr_->sequenceOut_;
3162	double lower,value,upper;
3163	if (iSequence<numberColumns) {
3164	lower=columnLower[iSequence];
3165	value=columnValue[iSequence];
3166	upper=columnUpper[iSequence];
3167	} else {
3168	iSequence -= numberColumns;
3169	lower=rowLower[iSequence];
3170	value=modelPtr_->rowActivity_[iSequence];
3171	upper=rowUpper[iSequence];
3172	}
3173	if (value<lower)
3174	bSum2=value-lower;
3175	else if (value>upper)
3176	bSum2=-(value-upper);
3177	else
3178	printf("OUCHXX %g %g %g\n",
3179	lower,value,upper);
3180	}
3181	#else
3182	if (modelPtr_->valueOut_<modelPtr_->lowerOut_)
3183	bSum2=modelPtr_->valueOut_-modelPtr_->lowerOut_;
3184	else if (modelPtr_->valueOut_>modelPtr_->upperOut_)
3185	bSum2=-(modelPtr_->valueOut_-modelPtr_->upperOut_);
3186	else
3187	printf("OUCHXX %g %g %g\n",
3188	modelPtr_->lowerOut_,modelPtr_->valueOut_,modelPtr_->upperOut_);
3189	#endif
3190	#if 0
3191	if (fabs(bSum-bSum2)>1.0e-5)
3192	printf("XX ");
3193	printf("after bounds - bSum %g - bSum2 %g\n",bSum,bSum2);
3194	#endif
3195	#endif
3196	modelPtr_->scaledMatrix_ = saveMatrix;
3197	modelPtr_->rowScale_ = saveScale;
3198	if (numberBad \|\| bSum > -1.0e-4 /\|\|nNonzeroBasic>1/ /\|\|bSum2>-1.0e-4/) {
3199	#if PRINT_CONFLICT > 1 //ndef NDEBUG
3200	printf("bad BOUND bSum %g (bSum2 %g) - %d bad - %d basic\n",
3201	bSum, bSum2, numberBad, nNonzeroBasic);
3202	#endif
3203	} else {
3204	if (numberColumns < 0)
3205	debugMode = -numberColumns;
3206	if ((debugMode & 4) != 0) {
3207	int *tempC = new int[numberColumns];
3208	double *temp = new double[numberColumns];
3209	int n = 0;
3210	for (int i = 0; i < numberColumns; i++) {
3211	if (fabs(farkas[i]) > 1.0e-12) {
3212	temp[n] = farkas[i];
3213	tempC[n++] = i;
3214	}
3215	}
3216	debugModel.addRow(n, tempC, temp, bSum, bSum);
3217	delete[] tempC;
3218	delete[] temp;
3219	}
3220	if ((debugMode & 2) != 0) {
3221	ClpSimplex dual = *modelPtr_;
3222	dual.setLogLevel(63);
3223	dual.scaling(0);
3224	dual.dual();
3225	assert(dual.problemStatus_ == 1);
3226	if (dual.ray_) {
3227	double *farkas2 = dual.reducedCost_;
3228	// Put nonzero bounds in farkas
3229	const double *columnLower = dual.columnLower_;
3230	const double *columnUpper = dual.columnUpper_;
3231	for (int i = 0; i < numberColumns; i++) {
3232	farkas2[i] = 0.0;
3233	if (dual.getStatus(i) == ClpSimplex::atLowerBound \|\| columnLower[i] == columnUpper[i]) {
3234	farkas2[i] = columnLower[i];
3235	} else if (dual.getStatus(i) == ClpSimplex::atUpperBound) {
3236	farkas2[i] = columnUpper[i];
3237	} else if (i == dual.sequenceOut_) {
3238	if (dual.directionOut_ < 0) {
3239	// above upper bound
3240	farkas2[i] = columnUpper[i];
3241	} else {
3242	// below lower bound
3243	farkas2[i] = columnLower[i];
3244	}
3245	} else if (columnLower[i] == columnUpper[i]) {
3246	farkas2[i] = columnLower[i];
3247	}
3248	}
3249	double *effectiveRhs2 = dual.rowActivity_;
3250	const double *rowLower = dual.rowLower_;
3251	const double *rowUpper = dual.rowUpper_;
3252	int pivotRow = dual.spareIntArray_[3];
3253	for (int i = 0; i < numberRows; i++) {
3254	if (dual.getRowStatus(i) == ClpSimplex::atLowerBound \|\| rowUpper[i] == rowLower[i] \|\| dual.getRowStatus(i) == ClpSimplex::isFixed) {
3255	effectiveRhs2[i] = rowLower[i];
3256	} else if (dual.getRowStatus(i) == ClpSimplex::atUpperBound) {
3257	effectiveRhs2[i] = rowUpper[i];
3258	} else {
3259	if (dual.getRowStatus(i) != ClpSimplex::basic) {
3260	assert(rowUpper[i] > 1.0e30 \|\| rowLower[i] < -1.0e30); // eventually skip
3261	if (rowUpper[i] < 1.0e30)
3262	effectiveRhs2[i] = rowUpper[i];
3263	else
3264	effectiveRhs2[i] = rowLower[i];
3265	}
3266	}
3267	if (dual.getRowStatus(i) == ClpSimplex::basic) {
3268	effectiveRhs2[i] = 0.0;
3269	// we should leave pivot row in and use direction for bound
3270	if (fabs(dual.ray_[i]) > 1.0e-8) {
3271	printf("Basic slack value %g on %d - pivotRow %d\n", ray[i], i, pivotRow);
3272	if (i == pivotRow) {
3273	if (dual.directionOut_ < 0)
3274	effectiveRhs[i] = rowUpper[i];
3275	else
3276	effectiveRhs[i] = rowLower[i];
3277	} else {
3278	dual.ray_[i] = 0.0;
3279	}
3280	}
3281	}
3282	}
3283	dual.times(-1.0, farkas2, effectiveRhs2);
3284	double bSum2 = 0.0;
3285	for (int i = 0; i < numberRows; i++) {
3286	bSum2 += effectiveRhs2[i] * dual.ray_[i];
3287	}
3288	printf("Alternate bSum %g\n", bSum2);
3289	memset(farkas2, 0, numberColumns * sizeof(double));
3290	dual.transposeTimes(-1.0, dual.ray_, farkas2);
3291	int nBad = 0;
3292	double largest = -1.0;
3293	double smallest = 1.0e30;
3294	for (int i = 0; i < numberColumns; i++) {
3295	//if (fabs(farkas[i])>1.0e-12)
3296	//printf("%d ptr farkas %g dual farkas %g\n",i,farkas[i],farkas2[i]);
3297	if (fabs(farkas[i] - farkas2[i]) > 1.0e-7) {
3298	nBad++;
3299	largest = CoinMax(largest, fabs(farkas[i] - farkas2[i]));
3300	smallest = CoinMin(smallest, fabs(farkas[i] - farkas2[i]));
3301	//printf("%d ptr farkas %g dual farkas %g\n",i,farkas[i],farkas2[i]);
3302	}
3303	}
3304	if (nBad)
3305	printf("%d farkas difference %g to %g\n", nBad, smallest, largest);
3306	dual.primal();
3307	assert(dual.problemStatus_ == 1);
3308	assert(!nBad);
3309	}
3310	}
3311	const char *integerInformation = modelPtr_->integerType_;
3312	assert(integerInformation);
3313	int *conflict = new int[numberColumns];
3314	double *sort = new double[numberColumns];
3315	double relax = 0.0;
3316	int nConflict = 0;
3317	int nOriginal = 0;
3318	int nFixed = 0;
3319	for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
3320	double thisRelax = 0.0;
3321	if (integerInformation[iColumn]) {
3322	if ((debugMode & 1) != 0)
3323	printf("%d status %d %g <= %g <=%g (orig %g, %g) farkas %g\n",
3324	iColumn, modelPtr_->getStatus(iColumn), columnLower[iColumn],
3325	modelPtr_->columnActivity_[iColumn], columnUpper[iColumn],
3326	originalLower[iColumn], originalUpper[iColumn],
3327	farkas[iColumn]);
3328	double gap = originalUpper[iColumn] - originalLower[iColumn];
3329	if (!gap)
3330	continue;
3331	if (gap == columnUpper[iColumn] - columnLower[iColumn])
3332	nOriginal++;
3333	if (columnUpper[iColumn] == columnLower[iColumn])
3334	nFixed++;
3335	if (fabs(farkas[iColumn]) < 1.0e-15) {
3336	farkas[iColumn] = 0.0;
3337	continue;
3338	}
3339	// temp
3340	if (gap >= 20000.0 && false) {
3341	// can't use
3342	if (farkas[iColumn] < 0.0) {
3343	assert(originalLower[iColumn] - columnLower[iColumn] <= 0.0);
3344	// farkas is negative - relax lower bound all way
3345	relax += farkas[iColumn] * (originalLower[iColumn] - columnLower[iColumn]);
3346	} else {
3347	assert(originalUpper[iColumn] - columnUpper[iColumn] >= 0.0);
3348	// farkas is positive - relax upper bound all way
3349	relax += farkas[iColumn] * (originalUpper[iColumn] - columnUpper[iColumn]);
3350	}
3351	continue;
3352	}
3353	if (originalLower[iColumn] == columnLower[iColumn]) {
3354	// may need to be careful if odd basic (due to local cut)
3355	if (farkas[iColumn] > 0.0 /*&&(modelPtr_->getStatus(iColumn)==ClpSimplex::atUpperBound
3356	\|\|modelPtr_->getStatus(iColumn)==ClpSimplex::isFixed
3357	\|\|iColumn==sequenceOut)*/
3358	) {
3359	// farkas is positive - add to list
3360	gap = originalUpper[iColumn] - columnUpper[iColumn];
3361	if (gap) {
3362	sort[nConflict] = -farkas[iColumn] * gap;
3363	conflict[nConflict++] = iColumn;
3364	}
3365	//assert (gap>columnUpper[iColumn]-columnLower[iColumn]);
3366	}
3367	} else if (originalUpper[iColumn] == columnUpper[iColumn]) {
3368	// may need to be careful if odd basic (due to local cut)
3369	if (farkas[iColumn] < 0.0 /*&&(modelPtr_->getStatus(iColumn)==ClpSimplex::atLowerBound
3370	\|\|modelPtr_->getStatus(iColumn)==ClpSimplex::isFixed
3371	\|\|iColumn==sequenceOut)*/
3372	) {
3373	// farkas is negative - add to list
3374	gap = columnLower[iColumn] - originalLower[iColumn];
3375	if (gap) {
3376	sort[nConflict] = farkas[iColumn] * gap;
3377	conflict[nConflict++] = iColumn;
3378	}
3379	//assert (gap>columnUpper[iColumn]-columnLower[iColumn]);
3380	}
3381	} else {
3382	// can't use
3383	if (farkas[iColumn] < 0.0) {
3384	assert(originalLower[iColumn] - columnLower[iColumn] <= 0.0);
3385	// farkas is negative - relax lower bound all way
3386	thisRelax = farkas[iColumn] * (originalLower[iColumn] - columnLower[iColumn]);
3387	} else {
3388	assert(originalUpper[iColumn] - columnUpper[iColumn] >= 0.0);
3389	// farkas is positive - relax upper bound all way
3390	thisRelax = farkas[iColumn] * (originalUpper[iColumn] - columnUpper[iColumn]);
3391	}
3392	}
3393	} else {
3394	// not integer - but may have been got at
3395	double gap = originalUpper[iColumn] - originalLower[iColumn];
3396	if (gap > columnUpper[iColumn] - columnLower[iColumn]) {
3397	// can't use
3398	if (farkas[iColumn] < 0.0) {
3399	assert(originalLower[iColumn] - columnLower[iColumn] <= 0.0);
3400	// farkas is negative - relax lower bound all way
3401	thisRelax = farkas[iColumn] * (originalLower[iColumn] - columnLower[iColumn]);
3402	} else {
3403	assert(originalUpper[iColumn] - columnUpper[iColumn] >= 0.0);
3404	// farkas is positive - relax upper bound all way
3405	thisRelax = farkas[iColumn] * (originalUpper[iColumn] - columnUpper[iColumn]);
3406	}
3407	}
3408	}
3409	assert(thisRelax >= 0.0);
3410	relax += thisRelax;
3411	}
3412	if (relax + bSum > -1.0e-4 \|\| !nConflict) {
3413	if (relax + bSum > -1.0e-4) {
3414	#if PRINT_CONFLICT > 1 //ndef NDEBUG
3415	printf("General integers relax bSum to %g\n", relax + bSum);
3416	#endif
3417	} else {
3418	#if PRINT_CONFLICT
3419	printf("All variables relaxed and still infeasible - what does this mean?\n");
3420	#endif
3421	int nR = 0;
3422	for (int i = 0; i < numberRows; i++) {
3423	if (fabs(ray[i]) > 1.0e-10)
3424	nR++;
3425	else
3426	ray[i] = 0.0;
3427	}
3428	int nC = 0;
3429	for (int i = 0; i < numberColumns; i++) {
3430	if (fabs(farkas[i]) > 1.0e-10)
3431	nC++;
3432	else
3433	farkas[i] = 0.0;
3434	}
3435	if (nR < 3 && nC < 5) {
3436	printf("BAD %d nonzero rows, %d nonzero columns\n", nR, nC);
3437	}
3438	}
3439	} else if (nConflict < 1000) {
3440	#if PRINT_CONFLICT > 1 //ndef NDEBUG
3441	if (nConflict < 5)
3442	printf("BOUNDS violation bSum %g (relaxed %g) - %d at original bounds, %d fixed - %d in conflict\n", bSum,
3443	relax + bSum, nOriginal, nFixed, nConflict);
3444	#endif
3445	CoinSort_2(sort, sort + nConflict, conflict);
3446	if ((debugMode & 4) != 0) {
3447	double *temp = new double[numberColumns];
3448	int *tempC = new int[numberColumns];
3449	int n = 0;
3450	for (int j = 0; j < nConflict; j++) {
3451	int i = conflict[j];
3452	if (fabs(farkas[i]) > 1.0e-12) {
3453	temp[n] = farkas[i];
3454	tempC[n++] = i;
3455	}
3456	}
3457	debugModel.addRow(n, tempC, temp, bSum, bSum);
3458	delete[] tempC;
3459	delete[] temp;
3460	}
3461	int nC = nConflict;
3462	for (int i = 0; i < nConflict; i++) {
3463	int iColumn = conflict[i];
3464	if (fabs(sort[i]) != fabs(farkas[iColumn]) && originalUpper[iColumn] == 1.0)
3465	printf("odd %d %g %d %g\n", i, sort[i], iColumn, farkas[iColumn]);
3466	}
3467	bSum += relax;
3468	double saveBsum = bSum;
3469	// we can't use same values
3470	double totalChange = 0;
3471	while (nConflict) {
3472	double last = -sort[nConflict - 1];
3473	int kConflict = nConflict;
3474	while (kConflict) {
3475	double change = -sort[kConflict - 1];
3476	if (change > last + 1.0e-5)
3477	break;
3478	totalChange += change;
3479	kConflict--;
3480	}
3481	if (bSum + totalChange > -1.0e-4)
3482	break;
3483	#if 0
3484	//int iColumn=conflict[nConflict-1];
3485	double change=-sort[nConflict-1];
3486	if (bSum+change>-1.0e-4)
3487	break;
3488	nConflict--;
3489	bSum += change;
3490	#else
3491	nConflict = kConflict;
3492	bSum += totalChange;
3493	#endif
3494	}
3495	if (!nConflict) {
3496	int nR = 0;
3497	for (int i = 0; i < numberRows; i++) {
3498	if (fabs(ray[i]) > 1.0e-10)
3499	nR++;
3500	else
3501	ray[i] = 0.0;
3502	}
3503	int nC = 0;
3504	for (int i = 0; i < numberColumns; i++) {
3505	if (fabs(farkas[i]) > 1.0e-10)
3506	nC++;
3507	else
3508	farkas[i] = 0.0;
3509	}
3510	#if 1 //PRINT_CONFLICT>1 //ndef NDEBUG
3511	{
3512	printf("BAD2 - zero nConflict %d nonzero rows, %d nonzero columns\n", nR, nC);
3513	}
3514	#endif
3515	}
3516	// no point doing if no reduction (or big?) ?
3517	if (nConflict < nC + 1 && nConflict < 100 && nConflict) {
3518	cut = new OsiRowCut();
3519	cut->setUb(COIN_DBL_MAX);
3520	if (!typeCut) {
3521	double lo = 1.0;
3522	for (int i = 0; i < nConflict; i++) {
3523	int iColumn = conflict[i];
3524	if (originalLower[iColumn] == columnLower[iColumn]) {
3525	// must be at least one higher
3526	sort[i] = 1.0;
3527	lo += originalLower[iColumn];
3528	} else {
3529	// must be at least one lower
3530	sort[i] = -1.0;
3531	lo -= originalUpper[iColumn];
3532	}
3533	}
3534	cut->setLb(lo);
3535	cut->setRow(nConflict, conflict, sort);
3536	#if PRINT_CONFLICT
3537	printf("CUT has %d (started at %d) - final bSum %g\n", nConflict, nC, bSum);
3538	#endif
3539	if ((debugMode & 4) != 0) {
3540	debugModel.addRow(nConflict, conflict, sort, lo, COIN_DBL_MAX);
3541	debugModel.writeMps("bad.mps");
3542	}
3543	} else {
3544	// just save for use later
3545	// first take off small
3546	int nC2 = nC;
3547	while (nC2) {
3548	//int iColumn=conflict[nConflict-1];
3549	double change = -sort[nC2 - 1];
3550	if (saveBsum + change > -1.0e-4 \|\| change > 1.0e-4)
3551	break;
3552	nC2--;
3553	saveBsum += change;
3554	}
3555	cut->setLb(saveBsum);
3556	for (int i = 0; i < nC2; i++) {
3557	int iColumn = conflict[i];
3558	sort[i] = farkas[iColumn];
3559	}
3560	cut->setRow(nC2, conflict, sort);
3561	// mark as globally valid
3562	cut->setGloballyValid();
3563	#if PRINT_CONFLICT > 1 //ndef NDEBUG
3564	printf("Stem CUT has %d (greedy %d - with small %d) - saved bSum %g final greedy bSum %g\n",
3565	nC2, nConflict, nC, saveBsum, bSum);
3566	#endif
3567	}
3568	}
3569	}
3570	delete[] conflict;
3571	delete[] sort;
3572	}
3573	delete[] farkas;
3574	} else {
3575	#if PRINT_CONFLICT > 1 //ndef NDEBUG
3576	printf("Bad dual ray\n");
3577	#endif
3578	}
3579	modelPtr_->deleteRay();
3580	}
3581	return cut;
3582	}
3583	#endif
3584
3585	//#############################################################################
3586	// Problem information methods (original data)
3587	//#############################################################################
3588
3589	//------------------------------------------------------------------
3590	const char *OsiClpSolverInterface::getRowSense() const
3591	{
3592	extractSenseRhsRange();
3593	return rowsense_;
3594	}
3595	//------------------------------------------------------------------
3596	const double *OsiClpSolverInterface::getRightHandSide() const
3597	{
3598	extractSenseRhsRange();
3599	return rhs_;
3600	}
3601	//------------------------------------------------------------------
3602	const double *OsiClpSolverInterface::getRowRange() const
3603	{
3604	extractSenseRhsRange();
3605	return rowrange_;
3606	}
3607	//------------------------------------------------------------------
3608	// Return information on integrality
3609	//------------------------------------------------------------------
3610	bool OsiClpSolverInterface::isContinuous(int colNumber) const
3611	{
3612	if (integerInformation_ == NULL)
3613	return true;
3614	#ifndef NDEBUG
3615	int n = modelPtr_->numberColumns();
3616	if (colNumber < 0 \|\| colNumber >= n) {
3617	indexError(colNumber, "isContinuous");
3618	}
3619	#endif
3620	if (integerInformation_[colNumber] == 0)
3621	return true;
3622	return false;
3623	}
3624	bool OsiClpSolverInterface::isBinary(int colNumber) const
3625	{
3626	#ifndef NDEBUG
3627	int n = modelPtr_->numberColumns();
3628	if (colNumber < 0 \|\| colNumber >= n) {
3629	indexError(colNumber, "isBinary");
3630	}
3631	#endif
3632	if (integerInformation_ == NULL \|\| integerInformation_[colNumber] == 0) {
3633	return false;
3634	} else {
3635	const double *cu = getColUpper();
3636	const double *cl = getColLower();
3637	if ((cu[colNumber] == 1 \|\| cu[colNumber] == 0) && (cl[colNumber] == 0 \|\| cl[colNumber] == 1))
3638	return true;
3639	else
3640	return false;
3641	}
3642	}
3643	//-----------------------------------------------------------------------------
3644	bool OsiClpSolverInterface::isInteger(int colNumber) const
3645	{
3646	#ifndef NDEBUG
3647	int n = modelPtr_->numberColumns();
3648	if (colNumber < 0 \|\| colNumber >= n) {
3649	indexError(colNumber, "isInteger");
3650	}
3651	#endif
3652	if (integerInformation_ == NULL \|\| integerInformation_[colNumber] == 0)
3653	return false;
3654	else
3655	return true;
3656	}
3657	//-----------------------------------------------------------------------------
3658	bool OsiClpSolverInterface::isIntegerNonBinary(int colNumber) const
3659	{
3660	#ifndef NDEBUG
3661	int n = modelPtr_->numberColumns();
3662	if (colNumber < 0 \|\| colNumber >= n) {
3663	indexError(colNumber, "isIntegerNonBinary");
3664	}
3665	#endif
3666	if (integerInformation_ == NULL \|\| integerInformation_[colNumber] == 0) {
3667	return false;
3668	} else {
3669	return !isBinary(colNumber);
3670	}
3671	}
3672	//-----------------------------------------------------------------------------
3673	bool OsiClpSolverInterface::isFreeBinary(int colNumber) const
3674	{
3675	#ifndef NDEBUG
3676	int n = modelPtr_->numberColumns();
3677	if (colNumber < 0 \|\| colNumber >= n) {
3678	indexError(colNumber, "isFreeBinary");
3679	}
3680	#endif
3681	if (integerInformation_ == NULL \|\| integerInformation_[colNumber] == 0) {
3682	return false;
3683	} else {
3684	const double *cu = getColUpper();
3685	const double *cl = getColLower();
3686	if ((cu[colNumber] == 1) && (cl[colNumber] == 0))
3687	return true;
3688	else
3689	return false;
3690	}
3691	}
3692	/* Return array of column length
3693	0 - continuous
3694	1 - binary (may get fixed later)
3695	2 - general integer (may get fixed later)
3696	*/
3697	const char *
3698	OsiClpSolverInterface::getColType(bool refresh) const
3699	{
3700	if (!columnType_ \|\| refresh) {
3701	const int numCols = getNumCols();
3702	if (!columnType_)
3703	columnType_ = new char[numCols];
3704	if (integerInformation_ == NULL) {
3705	memset(columnType_, 0, numCols);
3706	} else {
3707	const double *cu = getColUpper();
3708	const double *cl = getColLower();
3709	for (int i = 0; i < numCols; ++i) {
3710	if (integerInformation_[i]) {
3711	if ((cu[i] == 1 \|\| cu[i] == 0) && (cl[i] == 0 \|\| cl[i] == 1))
3712	columnType_[i] = 1;
3713	else
3714	columnType_[i] = 2;
3715	} else {
3716	columnType_[i] = 0;
3717	}
3718	}
3719	}
3720	}
3721	return columnType_;
3722	}
3723
3724	/* Return true if column is integer but does not have to
3725	be declared as such.
3726	Note: This function returns true if the the column
3727	is binary or a general integer.
3728	*/
3729	bool OsiClpSolverInterface::isOptionalInteger(int colNumber) const
3730	{
3731	#ifndef NDEBUG
3732	int n = modelPtr_->numberColumns();
3733	if (colNumber < 0 \|\| colNumber >= n) {
3734	indexError(colNumber, "isInteger");
3735	}
3736	#endif
3737	if (integerInformation_ == NULL \|\| integerInformation_[colNumber] != 2)
3738	return false;
3739	else
3740	return true;
3741	}
3742	/* Set the index-th variable to be an optional integer variable */
3743	void OsiClpSolverInterface::setOptionalInteger(int index)
3744	{
3745	if (!integerInformation_) {
3746	integerInformation_ = new char[modelPtr_->numberColumns()];
3747	CoinFillN(integerInformation_, modelPtr_->numberColumns(), static_cast< char >(0));
3748	}
3749	#ifndef NDEBUG
3750	int n = modelPtr_->numberColumns();
3751	if (index < 0 \|\| index >= n) {
3752	indexError(index, "setInteger");
3753	}
3754	#endif
3755	integerInformation_[index] = 2;
3756	modelPtr_->setInteger(index);
3757	}
3758
3759	//------------------------------------------------------------------
3760	// Row and column copies of the matrix ...
3761	//------------------------------------------------------------------
3762	const CoinPackedMatrix *OsiClpSolverInterface::getMatrixByRow() const
3763	{
3764	if (matrixByRow_ == NULL \|\| matrixByRow_->getNumElements() != modelPtr_->clpMatrix()->getNumElements()) {
3765	delete matrixByRow_;
3766	matrixByRow_ = new CoinPackedMatrix();
3767	matrixByRow_->setExtraGap(0.0);
3768	matrixByRow_->setExtraMajor(0.0);
3769	matrixByRow_->reverseOrderedCopyOf(*modelPtr_->matrix());
3770	//matrixByRow_->removeGaps();
3771	#if 0
3772	CoinPackedMatrix back;
3773	std::cout<<"start check"<<std::endl;
3774	back.reverseOrderedCopyOf(*matrixByRow_);
3775	modelPtr_->matrix()->isEquivalent2(back);
3776	std::cout<<"stop check"<<std::endl;
3777	#endif
3778	}
3779	assert(matrixByRow_->getNumElements() == modelPtr_->clpMatrix()->getNumElements());
3780	return matrixByRow_;
3781	}
3782
3783	const CoinPackedMatrix *OsiClpSolverInterface::getMatrixByCol() const
3784	{
3785	return modelPtr_->matrix();
3786	}
3787
3788	// Get pointer to mutable column-wise copy of matrix (returns NULL if not meaningful)
3789	CoinPackedMatrix *
3790	OsiClpSolverInterface::getMutableMatrixByCol() const
3791	{
3792	ClpPackedMatrix matrix = dynamic_cast< ClpPackedMatrix >(modelPtr_->matrix_);
3793	if (matrix)
3794	return matrix->getPackedMatrix();
3795	else
3796	return NULL;
3797	}
3798
3799	//------------------------------------------------------------------
3800	std::vector< double * > OsiClpSolverInterface::getDualRays(int /maxNumRays/,
3801	bool fullRay) const
3802	{
3803	return std::vector< double * >(1, modelPtr_->infeasibilityRay(fullRay));
3804	}
3805	//------------------------------------------------------------------
3806	std::vector< double * > OsiClpSolverInterface::getPrimalRays(int /maxNumRays/) const
3807	{
3808	return std::vector< double * >(1, modelPtr_->unboundedRay());
3809	}
3810	//#############################################################################
3811	void OsiClpSolverInterface::setContinuous(int index)
3812	{
3813
3814	if (integerInformation_) {
3815	#ifndef NDEBUG
3816	int n = modelPtr_->numberColumns();
3817	if (index < 0 \|\| index >= n) {
3818	indexError(index, "setContinuous");
3819	}
3820	#endif
3821	integerInformation_[index] = 0;
3822	}
3823	modelPtr_->setContinuous(index);
3824	}
3825	//-----------------------------------------------------------------------------
3826	void OsiClpSolverInterface::setInteger(int index)
3827	{
3828	if (!integerInformation_) {
3829	integerInformation_ = new char[modelPtr_->numberColumns()];
3830	CoinFillN(integerInformation_, modelPtr_->numberColumns(), static_cast< char >(0));
3831	}
3832	#ifndef NDEBUG
3833	int n = modelPtr_->numberColumns();
3834	if (index < 0 \|\| index >= n) {
3835	indexError(index, "setInteger");
3836	}
3837	#endif
3838	integerInformation_[index] = 1;
3839	modelPtr_->setInteger(index);
3840	}
3841	//-----------------------------------------------------------------------------
3842	void OsiClpSolverInterface::setContinuous(const int *indices, int len)
3843	{
3844	if (integerInformation_) {
3845	#ifndef NDEBUG
3846	int n = modelPtr_->numberColumns();
3847	#endif
3848	int i;
3849	for (i = 0; i < len; i++) {
3850	int colNumber = indices[i];
3851	#ifndef NDEBUG
3852	if (colNumber < 0 \|\| colNumber >= n) {
3853	indexError(colNumber, "setContinuous");
3854	}
3855	#endif
3856	integerInformation_[colNumber] = 0;
3857	modelPtr_->setContinuous(colNumber);
3858	}
3859	}
3860	}
3861	//-----------------------------------------------------------------------------
3862	void OsiClpSolverInterface::setInteger(const int *indices, int len)
3863	{
3864	if (!integerInformation_) {
3865	integerInformation_ = new char[modelPtr_->numberColumns()];
3866	CoinFillN(integerInformation_, modelPtr_->numberColumns(), static_cast< char >(0));
3867	}
3868	#ifndef NDEBUG
3869	int n = modelPtr_->numberColumns();
3870	#endif
3871	int i;
3872	for (i = 0; i < len; i++) {
3873	int colNumber = indices[i];
3874	#ifndef NDEBUG
3875	if (colNumber < 0 \|\| colNumber >= n) {
3876	indexError(colNumber, "setInteger");
3877	}
3878	#endif
3879	integerInformation_[colNumber] = 1;
3880	modelPtr_->setInteger(colNumber);
3881	}
3882	}
3883	/* Set the objective coefficients for all columns
3884	array [getNumCols()] is an array of values for the objective.
3885	This defaults to a series of set operations and is here for speed.
3886	*/
3887	void OsiClpSolverInterface::setObjective(const double *array)
3888	{
3889	// Say can't gurantee optimal basis etc
3890	lastAlgorithm_ = 999;
3891	modelPtr_->whatsChanged_ &= (0xffff & ~64);
3892	int n = modelPtr_->numberColumns();
3893	if (fakeMinInSimplex_) {
3894	std::transform(array, array + n,
3895	modelPtr_->objective(), std::negate< double >());
3896	} else {
3897	CoinMemcpyN(array, n, modelPtr_->objective());
3898	}
3899	}
3900	/* Set the lower bounds for all columns
3901	array [getNumCols()] is an array of values for the objective.
3902	This defaults to a series of set operations and is here for speed.
3903	*/
3904	void OsiClpSolverInterface::setColLower(const double *array)
3905	{
3906	// Say can't gurantee optimal basis etc
3907	lastAlgorithm_ = 999;
3908	modelPtr_->whatsChanged_ &= (0x1ffff & 128);
3909	CoinMemcpyN(array, modelPtr_->numberColumns(),
3910	modelPtr_->columnLower());
3911	}
3912	/* Set the upper bounds for all columns
3913	array [getNumCols()] is an array of values for the objective.
3914	This defaults to a series of set operations and is here for speed.
3915	*/
3916	void OsiClpSolverInterface::setColUpper(const double *array)
3917	{
3918	// Say can't gurantee optimal basis etc
3919	lastAlgorithm_ = 999;
3920	modelPtr_->whatsChanged_ &= (0x1ffff & 256);
3921	CoinMemcpyN(array, modelPtr_->numberColumns(),
3922	modelPtr_->columnUpper());
3923	}
3924	//-----------------------------------------------------------------------------
3925	void OsiClpSolverInterface::setColSolution(const double *cs)
3926	{
3927	// Say can't gurantee optimal basis etc
3928	lastAlgorithm_ = 999;
3929	CoinDisjointCopyN(cs, modelPtr_->numberColumns(),
3930	modelPtr_->primalColumnSolution());
3931	if (modelPtr_->solveType() == 2) {
3932	// directly into code as well
3933	CoinDisjointCopyN(cs, modelPtr_->numberColumns(),
3934	modelPtr_->solutionRegion(1));
3935	}
3936	// compute row activity
3937	memset(modelPtr_->primalRowSolution(), 0, modelPtr_->numberRows() * sizeof(double));
3938	modelPtr_->times(1.0, modelPtr_->primalColumnSolution(), modelPtr_->primalRowSolution());
3939	}
3940	//-----------------------------------------------------------------------------
3941	void OsiClpSolverInterface::setRowPrice(const double *rs)
3942	{
3943	CoinDisjointCopyN(rs, modelPtr_->numberRows(),
3944	modelPtr_->dualRowSolution());
3945	if (modelPtr_->solveType() == 2) {
3946	// directly into code as well (? sign )
3947	CoinDisjointCopyN(rs, modelPtr_->numberRows(),
3948	modelPtr_->djRegion(0));
3949	}
3950	// compute reduced costs
3951	memcpy(modelPtr_->dualColumnSolution(), modelPtr_->objective(),
3952	modelPtr_->numberColumns() * sizeof(double));
3953	modelPtr_->transposeTimes(-1.0, modelPtr_->dualRowSolution(), modelPtr_->dualColumnSolution());
3954	}
3955
3956	//#############################################################################
3957	// Problem modifying methods (matrix)
3958	//#############################################################################
3959	void OsiClpSolverInterface::addCol(const CoinPackedVectorBase &vec,
3960	const double collb, const double colub,
3961	const double obj)
3962	{
3963	int numberColumns = modelPtr_->numberColumns();
3964	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 8 \| 64 \| 128 \| 256));
3965	modelPtr_->resize(modelPtr_->numberRows(), numberColumns + 1);
3966	linearObjective_ = modelPtr_->objective();
3967	basis_.resize(modelPtr_->numberRows(), numberColumns + 1);
3968	setColBounds(numberColumns, collb, colub);
3969	setObjCoeff(numberColumns, obj);
3970	if (!modelPtr_->clpMatrix())
3971	modelPtr_->createEmptyMatrix();
3972	modelPtr_->matrix()->appendCol(vec);
3973	if (integerInformation_) {
3974	char *temp = new char[numberColumns + 1];
3975	CoinMemcpyN(integerInformation_, numberColumns, temp);
3976	delete[] integerInformation_;
3977	integerInformation_ = temp;
3978	integerInformation_[numberColumns] = 0;
3979	}
3980	freeCachedResults();
3981	}
3982	//-----------------------------------------------------------------------------
3983	/* Add a column (primal variable) to the problem. */
3984	void OsiClpSolverInterface::addCol(int numberElements, const int rows, const double elements,
3985	const double collb, const double colub,
3986	const double obj)
3987	{
3988	CoinPackedVector column(numberElements, rows, elements);
3989	addCol(column, collb, colub, obj);
3990	}
3991	// Add a named column (primal variable) to the problem.
3992	void OsiClpSolverInterface::addCol(const CoinPackedVectorBase &vec,
3993	const double collb, const double colub,
3994	const double obj, std::string name)
3995	{
3996	int ndx = getNumCols();
3997	addCol(vec, collb, colub, obj);
3998	setColName(ndx, name);
3999	}
4000	//-----------------------------------------------------------------------------
4001	void OsiClpSolverInterface::addCols(const int numcols,
4002	const CoinPackedVectorBase const cols,
4003	const double collb, const double colub,
4004	const double *obj)
4005	{
4006	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 8 \| 64 \| 128 \| 256));
4007	int numberColumns = modelPtr_->numberColumns();
4008	modelPtr_->resize(modelPtr_->numberRows(), numberColumns + numcols);
4009	linearObjective_ = modelPtr_->objective();
4010	basis_.resize(modelPtr_->numberRows(), numberColumns + numcols);
4011	double *lower = modelPtr_->columnLower() + numberColumns;
4012	double *upper = modelPtr_->columnUpper() + numberColumns;
4013	double *objective = modelPtr_->objective() + numberColumns;
4014	int iCol;
4015	if (collb) {
4016	for (iCol = 0; iCol < numcols; iCol++) {
4017	lower[iCol] = forceIntoRange(collb[iCol], -OsiClpInfinity, OsiClpInfinity);
4018	if (lower[iCol] < -1.0e27)
4019	lower[iCol] = -COIN_DBL_MAX;
4020	}
4021	} else {
4022	CoinFillN(lower, numcols, 0.0);
4023	}
4024	if (colub) {
4025	for (iCol = 0; iCol < numcols; iCol++) {
4026	upper[iCol] = forceIntoRange(colub[iCol], -OsiClpInfinity, OsiClpInfinity);
4027	if (upper[iCol] > 1.0e27)
4028	upper[iCol] = COIN_DBL_MAX;
4029	}
4030	} else {
4031	CoinFillN(upper, numcols, COIN_DBL_MAX);
4032	}
4033	if (obj) {
4034	for (iCol = 0; iCol < numcols; iCol++) {
4035	objective[iCol] = obj[iCol];
4036	}
4037	} else {
4038	CoinFillN(objective, numcols, 0.0);
4039	}
4040	if (!modelPtr_->clpMatrix())
4041	modelPtr_->createEmptyMatrix();
4042	modelPtr_->matrix()->appendCols(numcols, cols);
4043	if (integerInformation_) {
4044	char *temp = new char[numberColumns + numcols];
4045	CoinMemcpyN(integerInformation_, numberColumns, temp);
4046	delete[] integerInformation_;
4047	integerInformation_ = temp;
4048	for (iCol = 0; iCol < numcols; iCol++)
4049	integerInformation_[numberColumns + iCol] = 0;
4050	}
4051	freeCachedResults();
4052	}
4053	void OsiClpSolverInterface::addCols(const int numcols,
4054	const CoinBigIndex columnStarts, const int rows, const double *elements,
4055	const double collb, const double colub,
4056	const double *obj)
4057	{
4058	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 8 \| 64 \| 128 \| 256));
4059	int numberColumns = modelPtr_->numberColumns();
4060	modelPtr_->resize(modelPtr_->numberRows(), numberColumns + numcols);
4061	linearObjective_ = modelPtr_->objective();
4062	basis_.resize(modelPtr_->numberRows(), numberColumns + numcols);
4063	double *lower = modelPtr_->columnLower() + numberColumns;
4064	double *upper = modelPtr_->columnUpper() + numberColumns;
4065	double *objective = modelPtr_->objective() + numberColumns;
4066	int iCol;
4067	if (collb) {
4068	for (iCol = 0; iCol < numcols; iCol++) {
4069	lower[iCol] = forceIntoRange(collb[iCol], -OsiClpInfinity, OsiClpInfinity);
4070	if (lower[iCol] < -1.0e27)
4071	lower[iCol] = -COIN_DBL_MAX;
4072	}
4073	} else {
4074	CoinFillN(lower, numcols, 0.0);
4075	}
4076	if (colub) {
4077	for (iCol = 0; iCol < numcols; iCol++) {
4078	upper[iCol] = forceIntoRange(colub[iCol], -OsiClpInfinity, OsiClpInfinity);
4079	if (upper[iCol] > 1.0e27)
4080	upper[iCol] = COIN_DBL_MAX;
4081	}
4082	} else {
4083	CoinFillN(upper, numcols, COIN_DBL_MAX);
4084	}
4085	if (obj) {
4086	for (iCol = 0; iCol < numcols; iCol++) {
4087	objective[iCol] = obj[iCol];
4088	}
4089	} else {
4090	CoinFillN(objective, numcols, 0.0);
4091	}
4092	if (!modelPtr_->clpMatrix())
4093	modelPtr_->createEmptyMatrix();
4094	modelPtr_->matrix()->appendCols(numcols, columnStarts, rows, elements);
4095	if (integerInformation_) {
4096	char *temp = new char[numberColumns + numcols];
4097	CoinMemcpyN(integerInformation_, numberColumns, temp);
4098	delete[] integerInformation_;
4099	integerInformation_ = temp;
4100	for (iCol = 0; iCol < numcols; iCol++)
4101	integerInformation_[numberColumns + iCol] = 0;
4102	}
4103	freeCachedResults();
4104	}
4105	//-----------------------------------------------------------------------------
4106	void OsiClpSolverInterface::deleteCols(const int num, const int *columnIndices)
4107	{
4108	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 8 \| 64 \| 128 \| 256));
4109	findIntegers(false);
4110	deleteBranchingInfo(num, columnIndices);
4111	modelPtr_->deleteColumns(num, columnIndices);
4112	int nameDiscipline;
4113	getIntParam(OsiNameDiscipline, nameDiscipline);
4114	if (num && nameDiscipline) {
4115	// Very clumsy (and inefficient) - need to sort and then go backwards in ? chunks
4116	int *indices = CoinCopyOfArray(columnIndices, num);
4117	std::sort(indices, indices + num);
4118	int num2 = num;
4119	while (num2) {
4120	int next = indices[num2 - 1];
4121	int firstDelete = num2 - 1;
4122	int i;
4123	for (i = num2 - 2; i >= 0; i--) {
4124	if (indices[i] + 1 == next) {
4125	next--;
4126	firstDelete = i;
4127	} else {
4128	break;
4129	}
4130	}
4131	OsiSolverInterface::deleteColNames(indices[firstDelete], num2 - firstDelete);
4132	num2 = firstDelete;
4133	assert(num2 >= 0);
4134	}
4135	delete[] indices;
4136	}
4137	// synchronize integers (again)
4138	if (integerInformation_) {
4139	int numberColumns = modelPtr_->numberColumns();
4140	for (int i = 0; i < numberColumns; i++) {
4141	if (modelPtr_->isInteger(i))
4142	integerInformation_[i] = 1;
4143	else
4144	integerInformation_[i] = 0;
4145	}
4146	}
4147	basis_.deleteColumns(num, columnIndices);
4148	linearObjective_ = modelPtr_->objective();
4149	freeCachedResults();
4150	}
4151	//-----------------------------------------------------------------------------
4152	void OsiClpSolverInterface::addRow(const CoinPackedVectorBase &vec,
4153	const double rowlb, const double rowub)
4154	{
4155	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4156	freeCachedResults0();
4157	int numberRows = modelPtr_->numberRows();
4158	modelPtr_->resize(numberRows + 1, modelPtr_->numberColumns());
4159	basis_.resize(numberRows + 1, modelPtr_->numberColumns());
4160	setRowBounds(numberRows, rowlb, rowub);
4161	if (!modelPtr_->clpMatrix())
4162	modelPtr_->createEmptyMatrix();
4163	modelPtr_->matrix()->appendRow(vec);
4164	freeCachedResults1();
4165	}
4166	//-----------------------------------------------------------------------------
4167	void OsiClpSolverInterface::addRow(const CoinPackedVectorBase &vec,
4168	const double rowlb, const double rowub,
4169	std::string name)
4170	{
4171	int ndx = getNumRows();
4172	addRow(vec, rowlb, rowub);
4173	setRowName(ndx, name);
4174	}
4175	//-----------------------------------------------------------------------------
4176	void OsiClpSolverInterface::addRow(const CoinPackedVectorBase &vec,
4177	const char rowsen, const double rowrhs,
4178	const double rowrng)
4179	{
4180	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4181	freeCachedResults0();
4182	int numberRows = modelPtr_->numberRows();
4183	modelPtr_->resize(numberRows + 1, modelPtr_->numberColumns());
4184	basis_.resize(numberRows + 1, modelPtr_->numberColumns());
4185	double rowlb = 0, rowub = 0;
4186	convertSenseToBound(rowsen, rowrhs, rowrng, rowlb, rowub);
4187	setRowBounds(numberRows, rowlb, rowub);
4188	if (!modelPtr_->clpMatrix())
4189	modelPtr_->createEmptyMatrix();
4190	modelPtr_->matrix()->appendRow(vec);
4191	freeCachedResults1();
4192	}
4193	//-----------------------------------------------------------------------------
4194	void OsiClpSolverInterface::addRow(int numberElements, const int columns, const double elements,
4195	const double rowlb, const double rowub)
4196	{
4197	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4198	freeCachedResults0();
4199	int numberRows = modelPtr_->numberRows();
4200	modelPtr_->resize(numberRows + 1, modelPtr_->numberColumns());
4201	basis_.resize(numberRows + 1, modelPtr_->numberColumns());
4202	setRowBounds(numberRows, rowlb, rowub);
4203	if (!modelPtr_->clpMatrix())
4204	modelPtr_->createEmptyMatrix();
4205	modelPtr_->matrix()->appendRow(numberElements, columns, elements);
4206	CoinBigIndex starts[2];
4207	starts[0] = 0;
4208	starts[1] = numberElements;
4209	redoScaleFactors(1, starts, columns, elements);
4210	freeCachedResults1();
4211	}
4212	//-----------------------------------------------------------------------------
4213	void OsiClpSolverInterface::addRows(const int numrows,
4214	const CoinPackedVectorBase const rows,
4215	const double rowlb, const double rowub)
4216	{
4217	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4218	freeCachedResults0();
4219	int numberRows = modelPtr_->numberRows();
4220	modelPtr_->resize(numberRows + numrows, modelPtr_->numberColumns());
4221	basis_.resize(numberRows + numrows, modelPtr_->numberColumns());
4222	double *lower = modelPtr_->rowLower() + numberRows;
4223	double *upper = modelPtr_->rowUpper() + numberRows;
4224	int iRow;
4225	for (iRow = 0; iRow < numrows; iRow++) {
4226	if (rowlb)
4227	lower[iRow] = forceIntoRange(rowlb[iRow], -OsiClpInfinity, OsiClpInfinity);
4228	else
4229	lower[iRow] = -OsiClpInfinity;
4230	if (rowub)
4231	upper[iRow] = forceIntoRange(rowub[iRow], -OsiClpInfinity, OsiClpInfinity);
4232	else
4233	upper[iRow] = OsiClpInfinity;
4234	if (lower[iRow] < -1.0e27)
4235	lower[iRow] = -COIN_DBL_MAX;
4236	if (upper[iRow] > 1.0e27)
4237	upper[iRow] = COIN_DBL_MAX;
4238	}
4239	if (!modelPtr_->clpMatrix())
4240	modelPtr_->createEmptyMatrix();
4241	modelPtr_->matrix()->appendRows(numrows, rows);
4242	freeCachedResults1();
4243	}
4244	//-----------------------------------------------------------------------------
4245	void OsiClpSolverInterface::addRows(const int numrows,
4246	const CoinPackedVectorBase const rows,
4247	const char rowsen, const double rowrhs,
4248	const double *rowrng)
4249	{
4250	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4251	freeCachedResults0();
4252	int numberRows = modelPtr_->numberRows();
4253	modelPtr_->resize(numberRows + numrows, modelPtr_->numberColumns());
4254	basis_.resize(numberRows + numrows, modelPtr_->numberColumns());
4255	double *lower = modelPtr_->rowLower() + numberRows;
4256	double *upper = modelPtr_->rowUpper() + numberRows;
4257	int iRow;
4258	for (iRow = 0; iRow < numrows; iRow++) {
4259	double rowlb = 0, rowub = 0;
4260	convertSenseToBound(rowsen[iRow], rowrhs[iRow], rowrng[iRow],
4261	rowlb, rowub);
4262	lower[iRow] = forceIntoRange(rowlb, -OsiClpInfinity, OsiClpInfinity);
4263	upper[iRow] = forceIntoRange(rowub, -OsiClpInfinity, OsiClpInfinity);
4264	if (lower[iRow] < -1.0e27)
4265	lower[iRow] = -COIN_DBL_MAX;
4266	if (upper[iRow] > 1.0e27)
4267	upper[iRow] = COIN_DBL_MAX;
4268	}
4269	if (!modelPtr_->clpMatrix())
4270	modelPtr_->createEmptyMatrix();
4271	modelPtr_->matrix()->appendRows(numrows, rows);
4272	freeCachedResults1();
4273	}
4274	void OsiClpSolverInterface::addRows(const int numrows,
4275	const CoinBigIndex rowStarts, const int columns, const double *element,
4276	const double rowlb, const double rowub)
4277	{
4278	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4279	freeCachedResults0();
4280	int numberRows = modelPtr_->numberRows();
4281	modelPtr_->resize(numberRows + numrows, modelPtr_->numberColumns());
4282	basis_.resize(numberRows + numrows, modelPtr_->numberColumns());
4283	double *lower = modelPtr_->rowLower() + numberRows;
4284	double *upper = modelPtr_->rowUpper() + numberRows;
4285	int iRow;
4286	for (iRow = 0; iRow < numrows; iRow++) {
4287	if (rowlb)
4288	lower[iRow] = forceIntoRange(rowlb[iRow], -OsiClpInfinity, OsiClpInfinity);
4289	else
4290	lower[iRow] = -OsiClpInfinity;
4291	if (rowub)
4292	upper[iRow] = forceIntoRange(rowub[iRow], -OsiClpInfinity, OsiClpInfinity);
4293	else
4294	upper[iRow] = OsiClpInfinity;
4295	if (lower[iRow] < -1.0e27)
4296	lower[iRow] = -COIN_DBL_MAX;
4297	if (upper[iRow] > 1.0e27)
4298	upper[iRow] = COIN_DBL_MAX;
4299	}
4300	if (!modelPtr_->clpMatrix())
4301	modelPtr_->createEmptyMatrix();
4302	modelPtr_->matrix()->appendRows(numrows, rowStarts, columns, element);
4303	redoScaleFactors(numrows, rowStarts, columns, element);
4304	freeCachedResults1();
4305	}
4306	//-----------------------------------------------------------------------------
4307	void OsiClpSolverInterface::deleteRows(const int num, const int *rowIndices)
4308	{
4309	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
4310	// will still be optimal if all rows basic
4311	bool allBasic = true;
4312	int numBasis = basis_.getNumArtificial();
4313	for (int i = 0; i < num; i++) {
4314	int iRow = rowIndices[i];
4315	if (iRow < numBasis) {
4316	if (basis_.getArtifStatus(iRow) != CoinWarmStartBasis::basic) {
4317	allBasic = false;
4318	break;
4319	}
4320	}
4321	}
4322	int saveAlgorithm = allBasic ? lastAlgorithm_ : 999;
4323	modelPtr_->deleteRows(num, rowIndices);
4324	int nameDiscipline;
4325	getIntParam(OsiNameDiscipline, nameDiscipline);
4326	if (num && nameDiscipline) {
4327	// Very clumsy (and inefficient) - need to sort and then go backwards in ? chunks
4328	int *indices = CoinCopyOfArray(rowIndices, num);
4329	std::sort(indices, indices + num);
4330	int num2 = num;
4331	while (num2) {
4332	int next = indices[num2 - 1];
4333	int firstDelete = num2 - 1;
4334	int i;
4335	for (i = num2 - 2; i >= 0; i--) {
4336	if (indices[i] + 1 == next) {
4337	next--;
4338	firstDelete = i;
4339	} else {
4340	break;
4341	}
4342	}
4343	OsiSolverInterface::deleteRowNames(indices[firstDelete], num2 - firstDelete);
4344	num2 = firstDelete;
4345	assert(num2 >= 0);
4346	}
4347	delete[] indices;
4348	}
4349	basis_.deleteRows(num, rowIndices);
4350	CoinPackedMatrix *saveRowCopy = matrixByRow_;
4351	matrixByRow_ = NULL;
4352	freeCachedResults();
4353	modelPtr_->setNewRowCopy(NULL);
4354	delete modelPtr_->scaledMatrix_;
4355	modelPtr_->scaledMatrix_ = NULL;
4356	if (saveRowCopy) {
4357	matrixByRow_ = saveRowCopy;
4358	matrixByRow_->deleteRows(num, rowIndices);
4359	if (matrixByRow_->getNumElements() != modelPtr_->clpMatrix()->getNumElements()) {
4360	delete matrixByRow_; // odd type matrix
4361	matrixByRow_ = NULL;
4362	}
4363	}
4364	lastAlgorithm_ = saveAlgorithm;
4365	if ((specialOptions_ & 131072) != 0)
4366	lastNumberRows_ = modelPtr_->numberRows();
4367	}
4368
4369	//#############################################################################
4370	// Methods to input a problem
4371	//#############################################################################
4372
4373	void OsiClpSolverInterface::loadProblem(const CoinPackedMatrix &matrix,
4374	const double collb, const double colub,
4375	const double *obj,
4376	const double rowlb, const double rowub)
4377	{
4378	modelPtr_->whatsChanged_ = 0;
4379	// Get rid of integer information (modelPtr will get rid of its copy)
4380	delete[] integerInformation_;
4381	integerInformation_ = NULL;
4382	modelPtr_->loadProblem(matrix, collb, colub, obj, rowlb, rowub);
4383	linearObjective_ = modelPtr_->objective();
4384	freeCachedResults();
4385	basis_ = CoinWarmStartBasis();
4386	if (ws_) {
4387	delete ws_;
4388	ws_ = 0;
4389	}
4390	}
4391
4392	//-----------------------------------------------------------------------------
4393
4394	/*
4395	Expose the method that takes ClpMatrixBase. User request. Can't hurt, given
4396	the number of non-OSI methods already here.
4397	*/
4398	void OsiClpSolverInterface::loadProblem(const ClpMatrixBase &matrix,
4399	const double *collb,
4400	const double *colub,
4401	const double *obj,
4402	const double *rowlb,
4403	const double *rowub)
4404	{
4405	modelPtr_->whatsChanged_ = 0;
4406	// Get rid of integer information (modelPtr will get rid of its copy)
4407	delete[] integerInformation_;
4408	integerInformation_ = NULL;
4409	modelPtr_->loadProblem(matrix, collb, colub, obj, rowlb, rowub);
4410	linearObjective_ = modelPtr_->objective();
4411	freeCachedResults();
4412	basis_ = CoinWarmStartBasis();
4413	if (ws_) {
4414	delete ws_;
4415	ws_ = 0;
4416	}
4417	}
4418
4419	//-----------------------------------------------------------------------------
4420
4421	void OsiClpSolverInterface::assignProblem(CoinPackedMatrix *&matrix,
4422	double &collb, double &colub,
4423	double *&obj,
4424	double &rowlb, double &rowub)
4425	{
4426	modelPtr_->whatsChanged_ = 0;
4427	// Get rid of integer information (modelPtr will get rid of its copy)
4428	loadProblem(*matrix, collb, colub, obj, rowlb, rowub);
4429	delete matrix;
4430	matrix = NULL;
4431	delete[] collb;
4432	collb = NULL;
4433	delete[] colub;
4434	colub = NULL;
4435	delete[] obj;
4436	obj = NULL;
4437	delete[] rowlb;
4438	rowlb = NULL;
4439	delete[] rowub;
4440	rowub = NULL;
4441	}
4442
4443	//-----------------------------------------------------------------------------
4444
4445	void OsiClpSolverInterface::loadProblem(const CoinPackedMatrix &matrix,
4446	const double collb, const double colub,
4447	const double *obj,
4448	const char rowsen, const double rowrhs,
4449	const double *rowrng)
4450	{
4451	modelPtr_->whatsChanged_ = 0;
4452	// Get rid of integer information (modelPtr will get rid of its copy)
4453	// assert( rowsen != NULL );
4454	// assert( rowrhs != NULL );
4455	// If any of Rhs NULLs then create arrays
4456	int numrows = matrix.getNumRows();
4457	const char *rowsenUse = rowsen;
4458	if (!rowsen) {
4459	char *rowsen = new char[numrows];
4460	for (int i = 0; i < numrows; i++)
4461	rowsen[i] = 'G';
4462	rowsenUse = rowsen;
4463	}
4464	const double *rowrhsUse = rowrhs;
4465	if (!rowrhs) {
4466	double *rowrhs = new double[numrows];
4467	for (int i = 0; i < numrows; i++)
4468	rowrhs[i] = 0.0;
4469	rowrhsUse = rowrhs;
4470	}
4471	const double *rowrngUse = rowrng;
4472	if (!rowrng) {
4473	double *rowrng = new double[numrows];
4474	for (int i = 0; i < numrows; i++)
4475	rowrng[i] = 0.0;
4476	rowrngUse = rowrng;
4477	}
4478	double *rowlb = new double[numrows];
4479	double *rowub = new double[numrows];
4480	for (int i = numrows - 1; i >= 0; --i) {
4481	convertSenseToBound(rowsenUse[i], rowrhsUse[i], rowrngUse[i], rowlb[i], rowub[i]);
4482	}
4483	if (rowsen != rowsenUse)
4484	delete[] rowsenUse;
4485	if (rowrhs != rowrhsUse)
4486	delete[] rowrhsUse;
4487	if (rowrng != rowrngUse)
4488	delete[] rowrngUse;
4489	loadProblem(matrix, collb, colub, obj, rowlb, rowub);
4490	delete[] rowlb;
4491	delete[] rowub;
4492	}
4493
4494	//-----------------------------------------------------------------------------
4495
4496	void OsiClpSolverInterface::assignProblem(CoinPackedMatrix *&matrix,
4497	double &collb, double &colub,
4498	double *&obj,
4499	char &rowsen, double &rowrhs,
4500	double *&rowrng)
4501	{
4502	modelPtr_->whatsChanged_ = 0;
4503	// Get rid of integer information (modelPtr will get rid of its copy)
4504	loadProblem(*matrix, collb, colub, obj, rowsen, rowrhs, rowrng);
4505	delete matrix;
4506	matrix = NULL;
4507	delete[] collb;
4508	collb = NULL;
4509	delete[] colub;
4510	colub = NULL;
4511	delete[] obj;
4512	obj = NULL;
4513	delete[] rowsen;
4514	rowsen = NULL;
4515	delete[] rowrhs;
4516	rowrhs = NULL;
4517	delete[] rowrng;
4518	rowrng = NULL;
4519	}
4520
4521	//-----------------------------------------------------------------------------
4522
4523	void OsiClpSolverInterface::loadProblem(const int numcols, const int numrows,
4524	const CoinBigIndex start, const int index,
4525	const double *value,
4526	const double collb, const double colub,
4527	const double *obj,
4528	const double rowlb, const double rowub)
4529	{
4530	modelPtr_->whatsChanged_ = 0;
4531	// Get rid of integer information (modelPtr will get rid of its copy)
4532	delete[] integerInformation_;
4533	integerInformation_ = NULL;
4534	modelPtr_->loadProblem(numcols, numrows, start, index,
4535	value, collb, colub, obj,
4536	rowlb, rowub);
4537	linearObjective_ = modelPtr_->objective();
4538	freeCachedResults();
4539	basis_ = CoinWarmStartBasis();
4540	if (ws_) {
4541	delete ws_;
4542	ws_ = 0;
4543	}
4544	}
4545	//-----------------------------------------------------------------------------
4546
4547	void OsiClpSolverInterface::loadProblem(const int numcols, const int numrows,
4548	const CoinBigIndex start, const int index,
4549	const double *value,
4550	const double collb, const double colub,
4551	const double *obj,
4552	const char rowsen, const double rowrhs,
4553	const double *rowrng)
4554	{
4555	modelPtr_->whatsChanged_ = 0;
4556	// Get rid of integer information (modelPtr will get rid of its copy)
4557	// If any of Rhs NULLs then create arrays
4558	const char *rowsenUse = rowsen;
4559	if (!rowsen) {
4560	char *rowsen = new char[numrows];
4561	for (int i = 0; i < numrows; i++)
4562	rowsen[i] = 'G';
4563	rowsenUse = rowsen;
4564	}
4565	const double *rowrhsUse = rowrhs;
4566	if (!rowrhs) {
4567	double *rowrhs = new double[numrows];
4568	for (int i = 0; i < numrows; i++)
4569	rowrhs[i] = 0.0;
4570	rowrhsUse = rowrhs;
4571	}
4572	const double *rowrngUse = rowrng;
4573	if (!rowrng) {
4574	double *rowrng = new double[numrows];
4575	for (int i = 0; i < numrows; i++)
4576	rowrng[i] = 0.0;
4577	rowrngUse = rowrng;
4578	}
4579	double *rowlb = new double[numrows];
4580	double *rowub = new double[numrows];
4581	for (int i = numrows - 1; i >= 0; --i) {
4582	convertSenseToBound(rowsenUse[i], rowrhsUse[i], rowrngUse[i], rowlb[i], rowub[i]);
4583	}
4584	if (rowsen != rowsenUse)
4585	delete[] rowsenUse;
4586	if (rowrhs != rowrhsUse)
4587	delete[] rowrhsUse;
4588	if (rowrng != rowrngUse)
4589	delete[] rowrngUse;
4590	loadProblem(numcols, numrows, start, index, value, collb, colub, obj,
4591	rowlb, rowub);
4592	delete[] rowlb;
4593	delete[] rowub;
4594	}
4595	// This loads a model from a coinModel object - returns number of errors
4596	int OsiClpSolverInterface::loadFromCoinModel(CoinModel &modelObject, bool keepSolution)
4597	{
4598	modelPtr_->whatsChanged_ = 0;
4599	int numberErrors = 0;
4600	// Set arrays for normal use
4601	double *rowLower = modelObject.rowLowerArray();
4602	double *rowUpper = modelObject.rowUpperArray();
4603	double *columnLower = modelObject.columnLowerArray();
4604	double *columnUpper = modelObject.columnUpperArray();
4605	double *objective = modelObject.objectiveArray();
4606	int *integerType = modelObject.integerTypeArray();
4607	double *associated = modelObject.associatedArray();
4608	// If strings then do copies
4609	if (modelObject.stringsExist()) {
4610	numberErrors = modelObject.createArrays(rowLower, rowUpper, columnLower, columnUpper,
4611	objective, integerType, associated);
4612	}
4613	CoinPackedMatrix matrix;
4614	modelObject.createPackedMatrix(matrix, associated);
4615	int numberRows = modelObject.numberRows();
4616	int numberColumns = modelObject.numberColumns();
4617	CoinWarmStart *ws = getWarmStart();
4618	bool restoreBasis = keepSolution && numberRows && numberRows == getNumRows() && numberColumns == getNumCols();
4619	loadProblem(matrix,
4620	columnLower, columnUpper, objective, rowLower, rowUpper);
4621	if (restoreBasis)
4622	setWarmStart(ws);
4623	delete ws;
4624	// Do names if wanted
4625	int numberItems;
4626	numberItems = modelObject.rowNames()->numberItems();
4627	if (numberItems) {
4628	const char const rowNames = modelObject.rowNames()->names();
4629	modelPtr_->copyRowNames(rowNames, 0, numberItems);
4630	}
4631	numberItems = modelObject.columnNames()->numberItems();
4632	if (numberItems) {
4633	const char const columnNames = modelObject.columnNames()->names();
4634	modelPtr_->copyColumnNames(columnNames, 0, numberItems);
4635	}
4636	// Do integers if wanted
4637	assert(integerType);
4638	for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
4639	if (integerType[iColumn])
4640	setInteger(iColumn);
4641	}
4642	if (rowLower != modelObject.rowLowerArray() \|\| columnLower != modelObject.columnLowerArray()) {
4643	delete[] rowLower;
4644	delete[] rowUpper;
4645	delete[] columnLower;
4646	delete[] columnUpper;
4647	delete[] objective;
4648	delete[] integerType;
4649	delete[] associated;
4650	//if (numberErrors)
4651	// handler_->message(CLP_BAD_STRING_VALUES,messages_)
4652	// <<numberErrors
4653	// <<CoinMessageEol;
4654	}
4655	modelPtr_->optimizationDirection_ = modelObject.optimizationDirection();
4656	return numberErrors;
4657	}
4658
4659	//-----------------------------------------------------------------------------
4660	// Write mps files
4661	//-----------------------------------------------------------------------------
4662
4663	void OsiClpSolverInterface::writeMps(const char *filename,
4664	const char *extension,
4665	double objSense) const
4666	{
4667	std::string f(filename);
4668	std::string e(extension);
4669	std::string fullname;
4670	if (e != "") {
4671	fullname = f + "." + e;
4672	} else {
4673	// no extension so no trailing period
4674	fullname = f;
4675	}
4676	// get names
4677	const char const const rowNames = modelPtr_->rowNamesAsChar();
4678	const char const const columnNames = modelPtr_->columnNamesAsChar();
4679	// Fall back on Osi version - possibly with names
4680	OsiSolverInterface::writeMpsNative(fullname.c_str(),
4681	const_cast< const char ** >(rowNames),
4682	const_cast< const char ** >(columnNames), 0, 2, objSense,
4683	numberSOS_, setInfo_);
4684	if (rowNames) {
4685	modelPtr_->deleteNamesAsChar(rowNames, modelPtr_->numberRows_ + 1);
4686	modelPtr_->deleteNamesAsChar(columnNames, modelPtr_->numberColumns_);
4687	}
4688	}
4689
4690	int OsiClpSolverInterface::writeMpsNative(const char *filename,
4691	const char rowNames, const char columnNames,
4692	int formatType, int numberAcross, double objSense) const
4693	{
4694	return OsiSolverInterface::writeMpsNative(filename, rowNames, columnNames,
4695	formatType, numberAcross, objSense,
4696	numberSOS_, setInfo_);
4697	}
4698
4699	//#############################################################################
4700	// CLP specific public interfaces
4701	//#############################################################################
4702
4703	ClpSimplex *OsiClpSolverInterface::getModelPtr() const
4704	{
4705	int saveAlgorithm = lastAlgorithm_;
4706	//freeCachedResults();
4707	lastAlgorithm_ = saveAlgorithm;
4708	//bool inCbcOrOther = (modelPtr_->specialOptions()&0x03000000)!=0;
4709	return modelPtr_;
4710	}
4711	//-------------------------------------------------------------------
4712
4713	//#############################################################################
4714	// Constructors, destructors clone and assignment
4715	//#############################################################################
4716	//-------------------------------------------------------------------
4717	// Default Constructor
4718	//-------------------------------------------------------------------
4719	OsiClpSolverInterface::OsiClpSolverInterface()
4720	: OsiSolverInterface()
4721	, rowsense_(NULL)
4722	, rhs_(NULL)
4723	, rowrange_(NULL)
4724	, ws_(NULL)
4725	, rowActivity_(NULL)
4726	, columnActivity_(NULL)
4727	, numberSOS_(0)
4728	, setInfo_(NULL)
4729	, smallModel_(NULL)
4730	, factorization_(NULL)
4731	, smallestElementInCut_(1.0e-15)
4732	, smallestChangeInCut_(1.0e-10)
4733	, largestAway_(-1.0)
4734	, spareArrays_(NULL)
4735	, matrixByRow_(NULL)
4736	, matrixByRowAtContinuous_(NULL)
4737	, integerInformation_(NULL)
4738	, whichRange_(NULL)
4739	, fakeMinInSimplex_(false)
4740	, linearObjective_(NULL)
4741	, cleanupScaling_(0)
4742	, specialOptions_(0x80000000)
4743	, baseModel_(NULL)
4744	, lastNumberRows_(0)
4745	, continuousModel_(NULL)
4746	, fakeObjective_(NULL)
4747	{
4748	//printf("%p %d null constructor\n",this,xxxxxx);xxxxxx++;
4749	modelPtr_ = NULL;
4750	notOwned_ = false;
4751	disasterHandler_ = new OsiClpDisasterHandler();
4752	reset();
4753	}
4754
4755	//-------------------------------------------------------------------
4756	// Clone
4757	//-------------------------------------------------------------------
4758	OsiSolverInterface *OsiClpSolverInterface::clone(bool CopyData) const
4759	{
4760	//printf("in clone %x\n",this);
4761	OsiClpSolverInterface *newSolver;
4762	if (CopyData) {
4763	newSolver = new OsiClpSolverInterface(*this);
4764	} else {
4765	newSolver = new OsiClpSolverInterface();
4766	}
4767	#if 0
4768	const double * obj = newSolver->getObjCoefficients();
4769	const double * oldObj = getObjCoefficients();
4770	if(newSolver->getNumCols()>3787)
4771	printf("%x - obj %x (from %x) val %g\n",newSolver,obj,oldObj,obj[3787]);
4772	#endif
4773	return newSolver;
4774	}
4775
4776	//-------------------------------------------------------------------
4777	// Copy constructor
4778	//-------------------------------------------------------------------
4779	OsiClpSolverInterface::OsiClpSolverInterface(
4780	const OsiClpSolverInterface &rhs)
4781	: OsiSolverInterface(rhs)
4782	, rowsense_(NULL)
4783	, rhs_(NULL)
4784	, rowrange_(NULL)
4785	, ws_(NULL)
4786	, rowActivity_(NULL)
4787	, columnActivity_(NULL)
4788	, stuff_(rhs.stuff_)
4789	, numberSOS_(rhs.numberSOS_)
4790	, setInfo_(NULL)
4791	, smallModel_(NULL)
4792	, factorization_(NULL)
4793	, smallestElementInCut_(rhs.smallestElementInCut_)
4794	, smallestChangeInCut_(rhs.smallestChangeInCut_)
4795	, largestAway_(-1.0)
4796	, spareArrays_(NULL)
4797	, basis_()
4798	, itlimOrig_(9999999)
4799	, lastAlgorithm_(0)
4800	, notOwned_(false)
4801	, matrixByRow_(NULL)
4802	, matrixByRowAtContinuous_(NULL)
4803	, integerInformation_(NULL)
4804	, whichRange_(NULL)
4805	, fakeMinInSimplex_(rhs.fakeMinInSimplex_)
4806	{
4807	//printf("%p %d copy constructor %p\n",this,xxxxxx,&rhs);xxxxxx++;
4808	if (rhs.modelPtr_)
4809	modelPtr_ = new ClpSimplex(*rhs.modelPtr_);
4810	else
4811	modelPtr_ = new ClpSimplex();
4812	if (rhs.baseModel_)
4813	baseModel_ = new ClpSimplex(*rhs.baseModel_);
4814	else
4815	baseModel_ = NULL;
4816	if (rhs.continuousModel_)
4817	continuousModel_ = new ClpSimplex(*rhs.continuousModel_);
4818	else
4819	continuousModel_ = NULL;
4820	if (rhs.matrixByRowAtContinuous_)
4821	matrixByRowAtContinuous_ = new CoinPackedMatrix(*rhs.matrixByRowAtContinuous_);
4822	if (rhs.disasterHandler_)
4823	disasterHandler_ = dynamic_cast< OsiClpDisasterHandler * >(rhs.disasterHandler_->clone());
4824	else
4825	disasterHandler_ = NULL;
4826	if (rhs.fakeObjective_)
4827	fakeObjective_ = new ClpLinearObjective(*rhs.fakeObjective_);
4828	else
4829	fakeObjective_ = NULL;
4830	linearObjective_ = modelPtr_->objective();
4831	if (rhs.ws_)
4832	ws_ = new CoinWarmStartBasis(*rhs.ws_);
4833	basis_ = rhs.basis_;
4834	if (rhs.integerInformation_) {
4835	int numberColumns = modelPtr_->numberColumns();
4836	integerInformation_ = new char[numberColumns];
4837	CoinMemcpyN(rhs.integerInformation_, numberColumns, integerInformation_);
4838	}
4839	saveData_ = rhs.saveData_;
4840	solveOptions_ = rhs.solveOptions_;
4841	cleanupScaling_ = rhs.cleanupScaling_;
4842	specialOptions_ = rhs.specialOptions_;
4843	lastNumberRows_ = rhs.lastNumberRows_;
4844	rowScale_ = rhs.rowScale_;
4845	columnScale_ = rhs.columnScale_;
4846	fillParamMaps();
4847	messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
4848	if (numberSOS_) {
4849	setInfo_ = new CoinSet[numberSOS_];
4850	for (int i = 0; i < numberSOS_; i++)
4851	setInfo_[i] = rhs.setInfo_[i];
4852	}
4853	}
4854
4855	// Borrow constructor - only delete one copy
4856	OsiClpSolverInterface::OsiClpSolverInterface(ClpSimplex *rhs,
4857	bool reallyOwn)
4858	: OsiSolverInterface()
4859	, rowsense_(NULL)
4860	, rhs_(NULL)
4861	, rowrange_(NULL)
4862	, ws_(NULL)
4863	, rowActivity_(NULL)
4864	, columnActivity_(NULL)
4865	, numberSOS_(0)
4866	, setInfo_(NULL)
4867	, smallModel_(NULL)
4868	, factorization_(NULL)
4869	, smallestElementInCut_(1.0e-15)
4870	, smallestChangeInCut_(1.0e-10)
4871	, largestAway_(-1.0)
4872	, spareArrays_(NULL)
4873	, basis_()
4874	, itlimOrig_(9999999)
4875	, lastAlgorithm_(0)
4876	, notOwned_(false)
4877	, matrixByRow_(NULL)
4878	, matrixByRowAtContinuous_(NULL)
4879	, integerInformation_(NULL)
4880	, whichRange_(NULL)
4881	, fakeMinInSimplex_(false)
4882	, cleanupScaling_(0)
4883	, specialOptions_(0x80000000)
4884	, baseModel_(NULL)
4885	, lastNumberRows_(0)
4886	, continuousModel_(NULL)
4887	, fakeObjective_(NULL)
4888	{
4889	disasterHandler_ = new OsiClpDisasterHandler();
4890	//printf("in borrow %x - > %x\n",&rhs,this);
4891	modelPtr_ = rhs;
4892	basis_.resize(modelPtr_->numberRows(), modelPtr_->numberColumns());
4893	linearObjective_ = modelPtr_->objective();
4894	if (rhs) {
4895	notOwned_ = !reallyOwn;
4896
4897	if (rhs->integerInformation()) {
4898	int numberColumns = modelPtr_->numberColumns();
4899	integerInformation_ = new char[numberColumns];
4900	CoinMemcpyN(rhs->integerInformation(), numberColumns, integerInformation_);
4901	}
4902	}
4903	fillParamMaps();
4904	}
4905
4906	// Releases so won't error
4907	void OsiClpSolverInterface::releaseClp()
4908	{
4909	modelPtr_ = NULL;
4910	notOwned_ = false;
4911	}
4912
4913	//-------------------------------------------------------------------
4914	// Destructor
4915	//-------------------------------------------------------------------
4916	OsiClpSolverInterface::~OsiClpSolverInterface()
4917	{
4918	//printf("%p destructor\n",this);
4919	freeCachedResults();
4920	if (!notOwned_)
4921	delete modelPtr_;
4922	delete baseModel_;
4923	delete continuousModel_;
4924	delete disasterHandler_;
4925	delete fakeObjective_;
4926	delete ws_;
4927	delete[] rowActivity_;
4928	delete[] columnActivity_;
4929	delete[] setInfo_;
4930	#ifdef KEEP_SMALL
4931	if (smallModel_) {
4932	delete[] spareArrays_;
4933	spareArrays_ = NULL;
4934	delete smallModel_;
4935	smallModel_ = NULL;
4936	}
4937	#endif
4938	assert(smallModel_ == NULL);
4939	assert(factorization_ == NULL);
4940	assert(spareArrays_ == NULL);
4941	delete[] integerInformation_;
4942	delete matrixByRowAtContinuous_;
4943	delete matrixByRow_;
4944	}
4945
4946	//-------------------------------------------------------------------
4947	// Assignment operator
4948	//-------------------------------------------------------------------
4949	OsiClpSolverInterface &
4950	OsiClpSolverInterface::operator=(const OsiClpSolverInterface &rhs)
4951	{
4952	if (this != &rhs) {
4953	//printf("in = %x - > %x\n",&rhs,this);
4954	OsiSolverInterface::operator=(rhs);
4955	freeCachedResults();
4956	if (!notOwned_)
4957	delete modelPtr_;
4958	delete ws_;
4959	if (rhs.modelPtr_)
4960	modelPtr_ = new ClpSimplex(*rhs.modelPtr_);
4961	delete baseModel_;
4962	if (rhs.baseModel_)
4963	baseModel_ = new ClpSimplex(*rhs.baseModel_);
4964	else
4965	baseModel_ = NULL;
4966	delete continuousModel_;
4967	if (rhs.continuousModel_)
4968	continuousModel_ = new ClpSimplex(*rhs.continuousModel_);
4969	else
4970	continuousModel_ = NULL;
4971	delete matrixByRowAtContinuous_;
4972	delete matrixByRow_;
4973	matrixByRow_ = NULL;
4974	if (rhs.matrixByRowAtContinuous_)
4975	matrixByRowAtContinuous_ = new CoinPackedMatrix(*rhs.matrixByRowAtContinuous_);
4976	else
4977	matrixByRowAtContinuous_ = NULL;
4978	delete disasterHandler_;
4979	if (rhs.disasterHandler_)
4980	disasterHandler_ = dynamic_cast< OsiClpDisasterHandler * >(rhs.disasterHandler_->clone());
4981	else
4982	disasterHandler_ = NULL;
4983	delete fakeObjective_;
4984	if (rhs.fakeObjective_)
4985	fakeObjective_ = new ClpLinearObjective(*rhs.fakeObjective_);
4986	else
4987	fakeObjective_ = NULL;
4988	notOwned_ = false;
4989	linearObjective_ = modelPtr_->objective();
4990	saveData_ = rhs.saveData_;
4991	solveOptions_ = rhs.solveOptions_;
4992	cleanupScaling_ = rhs.cleanupScaling_;
4993	specialOptions_ = rhs.specialOptions_;
4994	lastNumberRows_ = rhs.lastNumberRows_;
4995	rowScale_ = rhs.rowScale_;
4996	columnScale_ = rhs.columnScale_;
4997	basis_ = rhs.basis_;
4998	stuff_ = rhs.stuff_;
4999	delete[] integerInformation_;
5000	integerInformation_ = NULL;
5001	if (rhs.integerInformation_) {
5002	int numberColumns = modelPtr_->numberColumns();
5003	integerInformation_ = new char[numberColumns];
5004	CoinMemcpyN(rhs.integerInformation_, numberColumns, integerInformation_);
5005	}
5006	if (rhs.ws_)
5007	ws_ = new CoinWarmStartBasis(*rhs.ws_);
5008	else
5009	ws_ = NULL;
5010	delete[] rowActivity_;
5011	delete[] columnActivity_;
5012	rowActivity_ = NULL;
5013	columnActivity_ = NULL;
5014	delete[] setInfo_;
5015	numberSOS_ = rhs.numberSOS_;
5016	setInfo_ = NULL;
5017	if (numberSOS_) {
5018	setInfo_ = new CoinSet[numberSOS_];
5019	for (int i = 0; i < numberSOS_; i++)
5020	setInfo_[i] = rhs.setInfo_[i];
5021	}
5022	assert(smallModel_ == NULL);
5023	assert(factorization_ == NULL);
5024	smallestElementInCut_ = rhs.smallestElementInCut_;
5025	smallestChangeInCut_ = rhs.smallestChangeInCut_;
5026	largestAway_ = -1.0;
5027	assert(spareArrays_ == NULL);
5028	basis_ = rhs.basis_;
5029	fillParamMaps();
5030	messageHandler()->setLogLevel(rhs.messageHandler()->logLevel());
5031	}
5032	return *this;
5033	}
5034
5035	//#############################################################################
5036	// Applying cuts
5037	//#############################################################################
5038
5039	void OsiClpSolverInterface::applyRowCut(const OsiRowCut &rowCut)
5040	{
5041	applyRowCuts(1, &rowCut);
5042	}
5043	/* Apply a collection of row cuts which are all effective.
5044	applyCuts seems to do one at a time which seems inefficient.
5045	*/
5046	void OsiClpSolverInterface::applyRowCuts(int numberCuts, const OsiRowCut *cuts)
5047	{
5048	if (numberCuts) {
5049	// Say can't gurantee optimal basis etc
5050	lastAlgorithm_ = 999;
5051
5052	// Thanks to js
5053	const OsiRowCut *cutsp = new const OsiRowCut [numberCuts];
5054	for (int i = 0; i < numberCuts; i++)
5055	cutsp[i] = &cuts[i];
5056
5057	applyRowCuts(numberCuts, cutsp);
5058
5059	delete[] cutsp;
5060	}
5061	}
5062	/* Apply a collection of row cuts which are all effective.
5063	applyCuts seems to do one at a time which seems inefficient.
5064	*/
5065	void OsiClpSolverInterface::applyRowCuts(int numberCuts, const OsiRowCut **cuts)
5066	{
5067	int i;
5068	if (!numberCuts)
5069	return;
5070	modelPtr_->whatsChanged_ &= (0xffff & ~(1 \| 2 \| 4 \| 16 \| 32));
5071	CoinPackedMatrix *saveRowCopy = matrixByRow_;
5072	matrixByRow_ = NULL;
5073	#if 0 // was #ifndef NDEBUG
5074	int nameDiscipline;
5075	getIntParam(OsiNameDiscipline,nameDiscipline) ;
5076	assert (!nameDiscipline);
5077	#endif
5078	freeCachedResults0();
5079	// Say can't gurantee optimal basis etc
5080	lastAlgorithm_ = 999;
5081	int numberRows = modelPtr_->numberRows();
5082	modelPtr_->resize(numberRows + numberCuts, modelPtr_->numberColumns());
5083	basis_.resize(numberRows + numberCuts, modelPtr_->numberColumns());
5084	// redo as relaxed - use modelPtr_-> addRows with starts etc
5085	int size = 0;
5086	for (i = 0; i < numberCuts; i++)
5087	size += cuts[i]->row().getNumElements();
5088	CoinBigIndex *starts = new CoinBigIndex[numberCuts + 1];
5089	int *indices = new int[size];
5090	double *elements = new double[size];
5091	double *lower = modelPtr_->rowLower() + numberRows;
5092	double *upper = modelPtr_->rowUpper() + numberRows;
5093	const double *columnLower = modelPtr_->columnLower();
5094	const double *columnUpper = modelPtr_->columnUpper();
5095	size = 0;
5096	for (i = 0; i < numberCuts; i++) {
5097	double rowLb = cuts[i]->lb();
5098	double rowUb = cuts[i]->ub();
5099	int n = cuts[i]->row().getNumElements();
5100	const int *index = cuts[i]->row().getIndices();
5101	const double *elem = cuts[i]->row().getElements();
5102	starts[i] = size;
5103	for (int j = 0; j < n; j++) {
5104	double value = elem[j];
5105	int column = index[j];
5106	if (fabs(value) >= smallestChangeInCut_) {
5107	// always take
5108	indices[size] = column;
5109	elements[size++] = value;
5110	} else if (fabs(value) >= smallestElementInCut_) {
5111	double lowerValue = columnLower[column];
5112	double upperValue = columnUpper[column];
5113	double difference = upperValue - lowerValue;
5114	if (difference < 1.0e20 && difference * fabs(value) < smallestChangeInCut_ && (rowLb < -1.0e20 \|\| rowUb > 1.0e20)) {
5115	// Take out and adjust to relax
5116	//printf("small el %g adjusted\n",value);
5117	if (rowLb > -1.0e20) {
5118	// just lower bound on row
5119	if (value > 0.0) {
5120	// pretend at upper
5121	rowLb -= value * upperValue;
5122	} else {
5123	// pretend at lower
5124	rowLb -= value * lowerValue;
5125	}
5126	} else {
5127	// just upper bound on row
5128	if (value > 0.0) {
5129	// pretend at lower
5130	rowUb -= value * lowerValue;
5131	} else {
5132	// pretend at upper
5133	rowUb -= value * upperValue;
5134	}
5135	}
5136	} else {
5137	// take (unwillingly)
5138	indices[size] = column;
5139	elements[size++] = value;
5140	}
5141	} else {
5142	//printf("small el %g ignored\n",value);
5143	}
5144	}
5145	lower[i] = forceIntoRange(rowLb, -OsiClpInfinity, OsiClpInfinity);
5146	upper[i] = forceIntoRange(rowUb, -OsiClpInfinity, OsiClpInfinity);
5147	if (lower[i] < -1.0e27)
5148	lower[i] = -COIN_DBL_MAX;
5149	if (upper[i] > 1.0e27)
5150	upper[i] = COIN_DBL_MAX;
5151	}
5152	starts[numberCuts] = size;
5153	if (!modelPtr_->clpMatrix())
5154	modelPtr_->createEmptyMatrix();
5155	//modelPtr_->matrix()->appendRows(numberCuts,rows);
5156	modelPtr_->clpMatrix()->appendMatrix(numberCuts, 0, starts, indices, elements);
5157	modelPtr_->setNewRowCopy(NULL);
5158	modelPtr_->setClpScaledMatrix(NULL);
5159	freeCachedResults1();
5160	redoScaleFactors(numberCuts, starts, indices, elements);
5161	if (saveRowCopy) {
5162	#if 1
5163	matrixByRow_ = saveRowCopy;
5164	matrixByRow_->appendRows(numberCuts, starts, indices, elements, 0);
5165	if (matrixByRow_->getNumElements() != modelPtr_->clpMatrix()->getNumElements()) {
5166	delete matrixByRow_; // odd type matrix
5167	matrixByRow_ = NULL;
5168	}
5169	#else
5170	delete saveRowCopy;
5171	#endif
5172	}
5173	delete[] starts;
5174	delete[] indices;
5175	delete[] elements;
5176	}
5177	//#############################################################################
5178	// Apply Cuts
5179	//#############################################################################
5180
5181	OsiSolverInterface::ApplyCutsReturnCode
5182	OsiClpSolverInterface::applyCuts(const OsiCuts &cs, double effectivenessLb)
5183	{
5184	OsiSolverInterface::ApplyCutsReturnCode retVal;
5185	int i;
5186
5187	// Loop once for each column cut
5188	for (i = 0; i < cs.sizeColCuts(); i++) {
5189	if (cs.colCut(i).effectiveness() < effectivenessLb) {
5190	retVal.incrementIneffective();
5191	continue;
5192	}
5193	if (!cs.colCut(i).consistent()) {
5194	retVal.incrementInternallyInconsistent();
5195	continue;
5196	}
5197	if (!cs.colCut(i).consistent(*this)) {
5198	retVal.incrementExternallyInconsistent();
5199	continue;
5200	}
5201	if (cs.colCut(i).infeasible(*this)) {
5202	retVal.incrementInfeasible();
5203	continue;
5204	}
5205	applyColCut(cs.colCut(i));
5206	retVal.incrementApplied();
5207	}
5208
5209	// Loop once for each row cut
5210	const OsiRowCut *addCuts = new const OsiRowCut [cs.sizeRowCuts()];
5211	int nAdd = 0;
5212	for (i = 0; i < cs.sizeRowCuts(); i++) {
5213	if (cs.rowCut(i).effectiveness() < effectivenessLb) {
5214	retVal.incrementIneffective();
5215	continue;
5216	}
5217	if (!cs.rowCut(i).consistent()) {
5218	retVal.incrementInternallyInconsistent();
5219	continue;
5220	}
5221	if (!cs.rowCut(i).consistent(*this)) {
5222	retVal.incrementExternallyInconsistent();
5223	continue;
5224	}
5225	if (cs.rowCut(i).infeasible(*this)) {
5226	retVal.incrementInfeasible();
5227	continue;
5228	}
5229	addCuts[nAdd++] = cs.rowCutPtr(i);
5230	retVal.incrementApplied();
5231	}
5232	// now apply
5233	applyRowCuts(nAdd, addCuts);
5234	delete[] addCuts;
5235
5236	return retVal;
5237	}
5238	// Extend scale factors
5239	void OsiClpSolverInterface::redoScaleFactors(int numberAdd, const CoinBigIndex *starts,
5240	const int indices, const double elements)
5241	{
5242	if ((specialOptions_ & 131072) != 0) {
5243	int numberRows = modelPtr_->numberRows() - numberAdd;
5244	assert(lastNumberRows_ == numberRows); // ???
5245	int iRow;
5246	int newNumberRows = numberRows + numberAdd;
5247	rowScale_.extend(static_cast< int >(2 * newNumberRows * sizeof(double)));
5248	double *rowScale = rowScale_.array();
5249	double *oldInverseScale = rowScale + lastNumberRows_;
5250	double *inverseRowScale = rowScale + newNumberRows;
5251	for (iRow = lastNumberRows_ - 1; iRow >= 0; iRow--)
5252	inverseRowScale[iRow] = oldInverseScale[iRow];
5253	//int numberColumns = baseModel_->numberColumns();
5254	const double *columnScale = columnScale_.array();
5255	//const double * inverseColumnScale = columnScale + numberColumns;
5256	// Geometric mean on row scales
5257	// adjust arrays
5258	rowScale += lastNumberRows_;
5259	inverseRowScale += lastNumberRows_;
5260	for (iRow = 0; iRow < numberAdd; iRow++) {
5261	CoinBigIndex j;
5262	double largest = 1.0e-20;
5263	double smallest = 1.0e50;
5264	for (j = starts[iRow]; j < starts[iRow + 1]; j++) {
5265	int iColumn = indices[j];
5266	double value = fabs(elements[j]);
5267	// Don't bother with tiny elements
5268	if (value > 1.0e-20) {
5269	value *= columnScale[iColumn];
5270	largest = CoinMax(largest, value);
5271	smallest = CoinMin(smallest, value);
5272	}
5273	}
5274	double scale = sqrt(smallest * largest);
5275	scale = CoinMax(1.0e-10, CoinMin(1.0e10, scale));
5276	inverseRowScale[iRow] = scale;
5277	rowScale[iRow] = 1.0 / scale;
5278	}
5279	lastNumberRows_ = newNumberRows;
5280	}
5281	}
5282	// Delete all scale factor stuff and reset option
5283	void OsiClpSolverInterface::deleteScaleFactors()
5284	{
5285	delete baseModel_;
5286	baseModel_ = NULL;
5287	lastNumberRows_ = 0;
5288	specialOptions_ &= ~131072;
5289	}
5290	//-----------------------------------------------------------------------------
5291
5292	void OsiClpSolverInterface::applyColCut(const OsiColCut &cc)
5293	{
5294	modelPtr_->whatsChanged_ &= (0x1ffff & ~(128 \| 256));
5295	// Say can't gurantee optimal basis etc
5296	lastAlgorithm_ = 999;
5297	double *lower = modelPtr_->columnLower();
5298	double *upper = modelPtr_->columnUpper();
5299	const CoinPackedVector &lbs = cc.lbs();
5300	const CoinPackedVector &ubs = cc.ubs();
5301	int i;
5302
5303	for (i = 0; i < lbs.getNumElements(); i++) {
5304	int iCol = lbs.getIndices()[i];
5305	double value = lbs.getElements()[i];
5306	if (value > lower[iCol])
5307	lower[iCol] = value;
5308	}
5309	for (i = 0; i < ubs.getNumElements(); i++) {
5310	int iCol = ubs.getIndices()[i];
5311	double value = ubs.getElements()[i];
5312	if (value < upper[iCol])
5313	upper[iCol] = value;
5314	}
5315	}
5316	//#############################################################################
5317	// Private methods
5318	//#############################################################################
5319
5320	//-------------------------------------------------------------------
5321
5322	void OsiClpSolverInterface::freeCachedResults() const
5323	{
5324	// Say can't gurantee optimal basis etc
5325	lastAlgorithm_ = 999;
5326	delete[] rowsense_;
5327	delete[] rhs_;
5328	delete[] rowrange_;
5329	delete matrixByRow_;
5330	//delete ws_;
5331	rowsense_ = NULL;
5332	rhs_ = NULL;
5333	rowrange_ = NULL;
5334	matrixByRow_ = NULL;
5335	//ws_ = NULL;
5336	if (!notOwned_ && modelPtr_) {
5337	if (modelPtr_->scaledMatrix_) {
5338	delete modelPtr_->scaledMatrix_;
5339	modelPtr_->scaledMatrix_ = NULL;
5340	}
5341	if (modelPtr_->clpMatrix()) {
5342	modelPtr_->clpMatrix()->refresh(modelPtr_); // make sure all clean
5343	#ifndef NDEBUG
5344	ClpPackedMatrix clpMatrix = dynamic_cast< ClpPackedMatrix >(modelPtr_->clpMatrix());
5345	if (clpMatrix) {
5346	if (clpMatrix->getNumCols())
5347	assert(clpMatrix->getNumRows() == modelPtr_->getNumRows());
5348	if (clpMatrix->getNumRows())
5349	assert(clpMatrix->getNumCols() == modelPtr_->getNumCols());
5350	}
5351	#endif
5352	}
5353	}
5354	}
5355
5356	//-------------------------------------------------------------------
5357
5358	void OsiClpSolverInterface::freeCachedResults0() const
5359	{
5360	delete[] rowsense_;
5361	delete[] rhs_;
5362	delete[] rowrange_;
5363	rowsense_ = NULL;
5364	rhs_ = NULL;
5365	rowrange_ = NULL;
5366	}
5367
5368	//-------------------------------------------------------------------
5369
5370	void OsiClpSolverInterface::freeCachedResults1() const
5371	{
5372	// Say can't gurantee optimal basis etc
5373	lastAlgorithm_ = 999;
5374	delete matrixByRow_;
5375	matrixByRow_ = NULL;
5376	//ws_ = NULL;
5377	if (modelPtr_ && modelPtr_->clpMatrix()) {
5378	delete modelPtr_->scaledMatrix_;
5379	modelPtr_->scaledMatrix_ = NULL;
5380	modelPtr_->clpMatrix()->refresh(modelPtr_); // make sure all clean
5381	#ifndef NDEBUG
5382	ClpPackedMatrix clpMatrix = dynamic_cast< ClpPackedMatrix >(modelPtr_->clpMatrix());
5383	if (clpMatrix) {
5384	assert(clpMatrix->getNumRows() == modelPtr_->getNumRows());
5385	assert(clpMatrix->getNumCols() == modelPtr_->getNumCols());
5386	}
5387	#endif
5388	}
5389	}
5390
5391	//------------------------------------------------------------------
5392	void OsiClpSolverInterface::extractSenseRhsRange() const
5393	{
5394	if (rowsense_ == NULL) {
5395	// all three must be NULL
5396	assert((rhs_ == NULL) && (rowrange_ == NULL));
5397
5398	int nr = modelPtr_->numberRows();
5399	if (nr != 0) {
5400	rowsense_ = new char[nr];
5401	rhs_ = new double[nr];
5402	rowrange_ = new double[nr];
5403	std::fill(rowrange_, rowrange_ + nr, 0.0);
5404
5405	const double *lb = modelPtr_->rowLower();
5406	const double *ub = modelPtr_->rowUpper();
5407
5408	int i;
5409	for (i = 0; i < nr; i++) {
5410	convertBoundToSense(lb[i], ub[i], rowsense_[i], rhs_[i], rowrange_[i]);
5411	}
5412	}
5413	}
5414	}
5415	// Set language
5416	void OsiClpSolverInterface::newLanguage(CoinMessages::Language language)
5417	{
5418	modelPtr_->newLanguage(language);
5419	OsiSolverInterface::newLanguage(language);
5420	}
5421	//#############################################################################
5422
5423	void OsiClpSolverInterface::fillParamMaps()
5424	{
5425	assert(static_cast< int >(OsiMaxNumIteration) == static_cast< int >(ClpMaxNumIteration));
5426	assert(static_cast< int >(OsiMaxNumIterationHotStart) == static_cast< int >(ClpMaxNumIterationHotStart));
5427	//assert (static_cast<int> (OsiLastIntParam)== static_cast<int>(ClpLastIntParam));
5428
5429	assert(static_cast< int >(OsiDualObjectiveLimit) == static_cast< int >(ClpDualObjectiveLimit));
5430	assert(static_cast< int >(OsiPrimalObjectiveLimit) == static_cast< int >(ClpPrimalObjectiveLimit));
5431	assert(static_cast< int >(OsiDualTolerance) == static_cast< int >(ClpDualTolerance));
5432	assert(static_cast< int >(OsiPrimalTolerance) == static_cast< int >(ClpPrimalTolerance));
5433	assert(static_cast< int >(OsiObjOffset) == static_cast< int >(ClpObjOffset));
5434	//assert (static_cast<int> (OsiLastDblParam)== static_cast<int>(ClpLastDblParam));
5435
5436	assert(static_cast< int >(OsiProbName) == static_cast< int >(ClpProbName));
5437	//strParamMap_[OsiLastStrParam] = ClpLastStrParam;
5438	}
5439	// Sets up basis
5440	void OsiClpSolverInterface::setBasis(const CoinWarmStartBasis &basis)
5441	{
5442	setBasis(basis, modelPtr_);
5443	setWarmStart(&basis);
5444	}
5445	// Warm start
5446	CoinWarmStartBasis
5447	OsiClpSolverInterface::getBasis(ClpSimplex *model) const
5448	{
5449	int iRow, iColumn;
5450	int numberRows = model->numberRows();
5451	int numberColumns = model->numberColumns();
5452	CoinWarmStartBasis basis;
5453	basis.setSize(numberColumns, numberRows);
5454	if (model->statusExists()) {
5455	// Flip slacks
5456	int lookupA[] = { 0, 1, 3, 2, 0, 2 };
5457	for (iRow = 0; iRow < numberRows; iRow++) {
5458	int iStatus = model->getRowStatus(iRow);
5459	iStatus = lookupA[iStatus];
5460	basis.setArtifStatus(iRow, static_cast< CoinWarmStartBasis::Status >(iStatus));
5461	}
5462	int lookupS[] = { 0, 1, 2, 3, 0, 3 };
5463	for (iColumn = 0; iColumn < numberColumns; iColumn++) {
5464	int iStatus = model->getColumnStatus(iColumn);
5465	iStatus = lookupS[iStatus];
5466	basis.setStructStatus(iColumn, static_cast< CoinWarmStartBasis::Status >(iStatus));
5467	}
5468	}
5469	//basis.print();
5470	return basis;
5471	}
5472	// Warm start from statusArray
5473	CoinWarmStartBasis *
5474	OsiClpSolverInterface::getBasis(const unsigned char *statusArray) const
5475	{
5476	int iRow, iColumn;
5477	int numberRows = modelPtr_->numberRows();
5478	int numberColumns = modelPtr_->numberColumns();
5479	CoinWarmStartBasis *basis = new CoinWarmStartBasis();
5480	basis->setSize(numberColumns, numberRows);
5481	// Flip slacks
5482	int lookupA[] = { 0, 1, 3, 2, 0, 2 };
5483	for (iRow = 0; iRow < numberRows; iRow++) {
5484	int iStatus = statusArray[numberColumns + iRow] & 7;
5485	iStatus = lookupA[iStatus];
5486	basis->setArtifStatus(iRow, static_cast< CoinWarmStart