Context Navigation

source: trunk/Clp/src/ClpMain.cpp @ 1333

Last change on this file since 1333 was 1333, checked in by forrest, 12 years ago
version 1.9
Property svn:eol-style set to `native` Property svn:keywords set to `Author Date Id Revision`
File size: 91.3 KB

Line
1	// Copyright (C) 2002, International Business Machines
2	// Corporation and others. All Rights Reserved.
3
4	#include "CoinPragma.hpp"
5
6	#include <cassert>
7	#include <cstdio>
8	#include <cstdlib>
9	#include <cmath>
10	#include <cfloat>
11	#include <string>
12	#include <iostream>
13
14	int boundary_sort=1000;
15	int boundary_sort2=1000;
16	int boundary_sort3=10000;
17
18	#include "CoinPragma.hpp"
19	#include "CoinHelperFunctions.hpp"
20	#include "CoinSort.hpp"
21	// History since 1.0 at end
22	#define CLPVERSION "1.09.00"
23
24	#include "CoinMpsIO.hpp"
25	#include "CoinFileIO.hpp"
26
27	#include "ClpFactorization.hpp"
28	#include "CoinTime.hpp"
29	#include "ClpSimplex.hpp"
30	#include "ClpSimplexOther.hpp"
31	#include "ClpSolve.hpp"
32	#include "ClpPackedMatrix.hpp"
33	#include "ClpPlusMinusOneMatrix.hpp"
34	#include "ClpNetworkMatrix.hpp"
35	#include "ClpDualRowSteepest.hpp"
36	#include "ClpDualRowDantzig.hpp"
37	#include "ClpLinearObjective.hpp"
38	#include "ClpPrimalColumnSteepest.hpp"
39	#include "ClpPrimalColumnDantzig.hpp"
40	#include "ClpPresolve.hpp"
41	#include "CbcOrClpParam.hpp"
42	#include "CoinSignal.hpp"
43	#ifdef DMALLOC
44	#include "dmalloc.h"
45	#endif
46	#ifdef WSSMP_BARRIER
47	#define FOREIGN_BARRIER
48	#endif
49	#ifdef UFL_BARRIER
50	#define FOREIGN_BARRIER
51	#endif
52	#ifdef TAUCS_BARRIER
53	#define FOREIGN_BARRIER
54	#endif
55
56	static double totalTime=0.0;
57	static bool maskMatches(const int * starts, char ** masks,
58	std::string & check);
59	static ClpSimplex * currentModel = NULL;
60
61	extern "C" {
62	static void
63	#if defined(_MSC_VER)
64	__cdecl
65	#endif // _MSC_VER
66	signal_handler(int whichSignal)
67	{
68	if (currentModel!=NULL)
69	currentModel->setMaximumIterations(0); // stop at next iterations
70	return;
71	}
72	}
73
74	//#############################################################################
75
76	#ifdef NDEBUG
77	#undef NDEBUG
78	#endif
79
80	int mainTest (int argc, const char *argv[],int algorithm,
81	ClpSimplex empty, bool doPresolve,int switchOff,bool doVector);
82	static void statistics(ClpSimplex * originalModel, ClpSimplex * model);
83	static void generateCode(const char * fileName,int type);
84	// Returns next valid field
85	int CbcOrClpRead_mode=1;
86	FILE * CbcOrClpReadCommand=stdin;
87
88	int
89	#if defined(_MSC_VER)
90	__cdecl
91	#endif // _MSC_VER
92	main (int argc, const char *argv[])
93	{
94	// next {} is just to make sure all memory should be freed - for debug
95	{
96	double time1 = CoinCpuTime(),time2;
97	// Set up all non-standard stuff
98	//int numberModels=1;
99	ClpSimplex * models = new ClpSimplex[1];
100
101	// default action on import
102	int allowImportErrors=0;
103	int keepImportNames=1;
104	int doIdiot=-1;
105	int outputFormat=2;
106	int slpValue=-1;
107	int cppValue=-1;
108	int printOptions=0;
109	int printMode=0;
110	int presolveOptions=0;
111	int doCrash=0;
112	int doVector=0;
113	int doSprint=-1;
114	// set reasonable defaults
115	int preSolve=5;
116	bool preSolveFile=false;
117	models->setPerturbation(50);
118	models->messageHandler()->setPrefix(false);
119	const char dirsep = CoinFindDirSeparator();
120	std::string directory;
121	std::string dirSample;
122	std::string dirNetlib;
123	std::string dirMiplib;
124	if (dirsep == '/') {
125	directory = "./";
126	dirSample = "../../Data/Sample/";
127	dirNetlib = "../../Data/Netlib/";
128	dirMiplib = "../../Data/miplib3/";
129	} else {
130	directory = ".\\";
131	dirSample = "..\\..\\Data\\Sample\\";
132	dirNetlib = "..\\..\\Data\\Netlib\\";
133	dirMiplib = "..\\..\\Data\\miplib3\\";
134	}
135	std::string defaultDirectory = directory;
136	std::string importFile ="";
137	std::string exportFile ="default.mps";
138	std::string importBasisFile ="";
139	int basisHasValues=0;
140	int substitution=3;
141	int dualize=3; // dualize if looks promising
142	std::string exportBasisFile ="default.bas";
143	std::string saveFile ="default.prob";
144	std::string restoreFile ="default.prob";
145	std::string solutionFile ="stdout";
146	std::string solutionSaveFile ="solution.file";
147	std::string printMask="";
148	CbcOrClpParam parameters[CBCMAXPARAMETERS];
149	int numberParameters ;
150	establishParams(numberParameters,parameters) ;
151	parameters[whichParam(BASISIN,numberParameters,parameters)].setStringValue(importBasisFile);
152	parameters[whichParam(BASISOUT,numberParameters,parameters)].setStringValue(exportBasisFile);
153	parameters[whichParam(PRINTMASK,numberParameters,parameters)].setStringValue(printMask);
154	parameters[whichParam(DIRECTORY,numberParameters,parameters)].setStringValue(directory);
155	parameters[whichParam(DIRSAMPLE,numberParameters,parameters)].setStringValue(dirSample);
156	parameters[whichParam(DIRNETLIB,numberParameters,parameters)].setStringValue(dirNetlib);
157	parameters[whichParam(DIRMIPLIB,numberParameters,parameters)].setStringValue(dirMiplib);
158	parameters[whichParam(DUALBOUND,numberParameters,parameters)].setDoubleValue(models->dualBound());
159	parameters[whichParam(DUALTOLERANCE,numberParameters,parameters)].setDoubleValue(models->dualTolerance());
160	parameters[whichParam(EXPORT,numberParameters,parameters)].setStringValue(exportFile);
161	parameters[whichParam(IDIOT,numberParameters,parameters)].setIntValue(doIdiot);
162	parameters[whichParam(IMPORT,numberParameters,parameters)].setStringValue(importFile);
163	parameters[whichParam(SOLVERLOGLEVEL,numberParameters,parameters)].setIntValue(models->logLevel());
164	parameters[whichParam(MAXFACTOR,numberParameters,parameters)].setIntValue(models->factorizationFrequency());
165	parameters[whichParam(MAXITERATION,numberParameters,parameters)].setIntValue(models->maximumIterations());
166	parameters[whichParam(OUTPUTFORMAT,numberParameters,parameters)].setIntValue(outputFormat);
167	parameters[whichParam(PRESOLVEPASS,numberParameters,parameters)].setIntValue(preSolve);
168	parameters[whichParam(PERTVALUE,numberParameters,parameters)].setIntValue(models->perturbation());
169	parameters[whichParam(PRIMALTOLERANCE,numberParameters,parameters)].setDoubleValue(models->primalTolerance());
170	parameters[whichParam(PRIMALWEIGHT,numberParameters,parameters)].setDoubleValue(models->infeasibilityCost());
171	parameters[whichParam(RESTORE,numberParameters,parameters)].setStringValue(restoreFile);
172	parameters[whichParam(SAVE,numberParameters,parameters)].setStringValue(saveFile);
173	parameters[whichParam(TIMELIMIT,numberParameters,parameters)].setDoubleValue(models->maximumSeconds());
174	parameters[whichParam(SOLUTION,numberParameters,parameters)].setStringValue(solutionFile);
175	parameters[whichParam(SAVESOL,numberParameters,parameters)].setStringValue(solutionSaveFile);
176	parameters[whichParam(SPRINT,numberParameters,parameters)].setIntValue(doSprint);
177	parameters[whichParam(SUBSTITUTION,numberParameters,parameters)].setIntValue(substitution);
178	parameters[whichParam(DUALIZE,numberParameters,parameters)].setIntValue(dualize);
179	parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].setDoubleValue(1.0e-8);
180	int verbose=0;
181
182	// total number of commands read
183	int numberGoodCommands=0;
184	bool * goodModels = new bool[1];
185
186	// Hidden stuff for barrier
187	int choleskyType = 0;
188	int gamma=0;
189	int scaleBarrier=0;
190	int doKKT=0;
191	int crossover=2; // do crossover unless quadratic
192
193	int iModel=0;
194	goodModels[0]=false;
195	//models[0].scaling(1);
196	//models[0].setDualBound(1.0e6);
197	//models[0].setDualTolerance(1.0e-7);
198	//ClpDualRowSteepest steep;
199	//models[0].setDualRowPivotAlgorithm(steep);
200	//models[0].setPrimalTolerance(1.0e-7);
201	//ClpPrimalColumnSteepest steepP;
202	//models[0].setPrimalColumnPivotAlgorithm(steepP);
203	std::string field;
204	std::cout<<"Coin LP version "<<CLPVERSION
205	<<", build "<<__DATE__<<std::endl;
206	// Print command line
207	if (argc>1) {
208	printf("command line - ");
209	for (int i=0;i<argc;i++)
210	printf("%s ",argv[i]);
211	printf("\n");
212	}
213
214	while (1) {
215	// next command
216	field=CoinReadGetCommand(argc,argv);
217
218	// exit if null or similar
219	if (!field.length()) {
220	if (numberGoodCommands==1&&goodModels[0]) {
221	// we just had file name - do dual or primal
222	field="either";
223	} else if (!numberGoodCommands) {
224	// let's give the sucker a hint
225	std::cout
226	<<"Clp takes input from arguments ( - switches to stdin)"
227	<<std::endl
228	<<"Enter ? for list of commands or help"<<std::endl;
229	field="-";
230	} else {
231	break;
232	}
233	}
234
235	// see if ? at end
236	int numberQuery=0;
237	if (field!="?"&&field!="???") {
238	int length = field.length();
239	int i;
240	for (i=length-1;i>0;i--) {
241	if (field[i]=='?')
242	numberQuery++;
243	else
244	break;
245	}
246	field=field.substr(0,length-numberQuery);
247	}
248	// find out if valid command
249	int iParam;
250	int numberMatches=0;
251	int firstMatch=-1;
252	for ( iParam=0; iParam<numberParameters; iParam++ ) {
253	int match = parameters[iParam].matches(field);
254	if (match==1) {
255	numberMatches = 1;
256	firstMatch=iParam;
257	break;
258	} else {
259	if (match&&firstMatch<0)
260	firstMatch=iParam;
261	numberMatches += match>>1;
262	}
263	}
264	if (iParam<numberParameters&&!numberQuery) {
265	// found
266	CbcOrClpParam found = parameters[iParam];
267	CbcOrClpParameterType type = found.type();
268	int valid;
269	numberGoodCommands++;
270	if (type==GENERALQUERY) {
271	std::cout<<"In argument list keywords have leading - "
272	", -stdin or just - switches to stdin"<<std::endl;
273	std::cout<<"One command per line (and no -)"<<std::endl;
274	std::cout<<"abcd? gives list of possibilities, if only one + explanation"<<std::endl;
275	std::cout<<"abcd?? adds explanation, if only one fuller help"<<std::endl;
276	std::cout<<"abcd without value (where expected) gives current value"<<std::endl;
277	std::cout<<"abcd value sets value"<<std::endl;
278	std::cout<<"Commands are:"<<std::endl;
279	int maxAcross=5;
280	bool evenHidden=false;
281	if ((verbose&8)!=0) {
282	// even hidden
283	evenHidden = true;
284	verbose &= ~8;
285	}
286	if (verbose)
287	maxAcross=1;
288	int limits[]={1,101,201,301,401};
289	std::vector<std::string> types;
290	types.push_back("Double parameters:");
291	types.push_back("Int parameters:");
292	types.push_back("Keyword parameters:");
293	types.push_back("Actions or string parameters:");
294	int iType;
295	for (iType=0;iType<4;iType++) {
296	int across=0;
297	if ((verbose%4)!=0)
298	std::cout<<std::endl;
299	std::cout<<types[iType]<<std::endl;
300	if ((verbose&2)!=0)
301	std::cout<<std::endl;
302	for ( iParam=0; iParam<numberParameters; iParam++ ) {
303	int type = parameters[iParam].type();
304	if ((parameters[iParam].displayThis()\|\|evenHidden)&&
305	type>=limits[iType]
306	&&type<limits[iType+1]) {
307	if (!across) {
308	if ((verbose&2)==0)
309	std::cout<<" ";
310	else
311	std::cout<<"Command ";
312	}
313	std::cout<<parameters[iParam].matchName()<<" ";
314	across++;
315	if (across==maxAcross) {
316	across=0;
317	if (verbose) {
318	// put out description as well
319	if ((verbose&1)!=0)
320	std::cout<<parameters[iParam].shortHelp();
321	std::cout<<std::endl;
322	if ((verbose&2)!=0) {
323	std::cout<<"---- description"<<std::endl;
324	parameters[iParam].printLongHelp();
325	std::cout<<"----"<<std::endl<<std::endl;
326	}
327	} else {
328	std::cout<<std::endl;
329	}
330	}
331	}
332	}
333	if (across)
334	std::cout<<std::endl;
335	}
336	} else if (type==FULLGENERALQUERY) {
337	std::cout<<"Full list of commands is:"<<std::endl;
338	int maxAcross=5;
339	int limits[]={1,101,201,301,401};
340	std::vector<std::string> types;
341	types.push_back("Double parameters:");
342	types.push_back("Int parameters:");
343	types.push_back("Keyword parameters and others:");
344	types.push_back("Actions:");
345	int iType;
346	for (iType=0;iType<4;iType++) {
347	int across=0;
348	std::cout<<types[iType]<<std::endl;
349	for ( iParam=0; iParam<numberParameters; iParam++ ) {
350	int type = parameters[iParam].type();
351	if (type>=limits[iType]
352	&&type<limits[iType+1]) {
353	if (!across)
354	std::cout<<" ";
355	std::cout<<parameters[iParam].matchName()<<" ";
356	across++;
357	if (across==maxAcross) {
358	std::cout<<std::endl;
359	across=0;
360	}
361	}
362	}
363	if (across)
364	std::cout<<std::endl;
365	}
366	} else if (type<101) {
367	// get next field as double
368	double value = CoinReadGetDoubleField(argc,argv,&valid);
369	if (!valid) {
370	parameters[iParam].setDoubleParameter(models+iModel,value);
371	} else if (valid==1) {
372	std::cout<<" is illegal for double parameter "<<parameters[iParam].name()<<" value remains "<<
373	parameters[iParam].doubleValue()<<std::endl;
374	} else {
375	std::cout<<parameters[iParam].name()<<" has value "<<
376	parameters[iParam].doubleValue()<<std::endl;
377	}
378	} else if (type<201) {
379	// get next field as int
380	int value = CoinReadGetIntField(argc,argv,&valid);
381	if (!valid) {
382	if (parameters[iParam].type()==PRESOLVEPASS)
383	preSolve = value;
384	else if (parameters[iParam].type()==IDIOT)
385	doIdiot = value;
386	else if (parameters[iParam].type()==SPRINT)
387	doSprint = value;
388	else if (parameters[iParam].type()==OUTPUTFORMAT)
389	outputFormat = value;
390	else if (parameters[iParam].type()==SLPVALUE)
391	slpValue = value;
392	else if (parameters[iParam].type()==CPP)
393	cppValue = value;
394	else if (parameters[iParam].type()==PRESOLVEOPTIONS)
395	presolveOptions = value;
396	else if (parameters[iParam].type()==PRINTOPTIONS)
397	printOptions = value;
398	else if (parameters[iParam].type()==SUBSTITUTION)
399	substitution = value;
400	else if (parameters[iParam].type()==DUALIZE)
401	dualize = value;
402	else if (parameters[iParam].type()==VERBOSE)
403	verbose = value;
404	parameters[iParam].setIntParameter(models+iModel,value);
405	} else if (valid==1) {
406	std::cout<<" is illegal for integer parameter "<<parameters[iParam].name()<<" value remains "<<
407	parameters[iParam].intValue()<<std::endl;
408	} else {
409	std::cout<<parameters[iParam].name()<<" has value "<<
410	parameters[iParam].intValue()<<std::endl;
411	}
412	} else if (type<301) {
413	// one of several strings
414	std::string value = CoinReadGetString(argc,argv);
415	int action = parameters[iParam].parameterOption(value);
416	if (action<0) {
417	if (value!="EOL") {
418	// no match
419	parameters[iParam].printOptions();
420	} else {
421	// print current value
422	std::cout<<parameters[iParam].name()<<" has value "<<
423	parameters[iParam].currentOption()<<std::endl;
424	}
425	} else {
426	parameters[iParam].setCurrentOption(action);
427	// for now hard wired
428	switch (type) {
429	case DIRECTION:
430	if (action==0)
431	models[iModel].setOptimizationDirection(1);
432	else if (action==1)
433	models[iModel].setOptimizationDirection(-1);
434	else
435	models[iModel].setOptimizationDirection(0);
436	break;
437	case DUALPIVOT:
438	if (action==0) {
439	ClpDualRowSteepest steep(3);
440	models[iModel].setDualRowPivotAlgorithm(steep);
441	} else if (action==1) {
442	//ClpDualRowDantzig dantzig;
443	ClpDualRowSteepest dantzig(5);
444	models[iModel].setDualRowPivotAlgorithm(dantzig);
445	} else if (action==2) {
446	// partial steep
447	ClpDualRowSteepest steep(2);
448	models[iModel].setDualRowPivotAlgorithm(steep);
449	} else {
450	ClpDualRowSteepest steep;
451	models[iModel].setDualRowPivotAlgorithm(steep);
452	}
453	break;
454	case PRIMALPIVOT:
455	if (action==0) {
456	ClpPrimalColumnSteepest steep(3);
457	models[iModel].setPrimalColumnPivotAlgorithm(steep);
458	} else if (action==1) {
459	ClpPrimalColumnSteepest steep(0);
460	models[iModel].setPrimalColumnPivotAlgorithm(steep);
461	} else if (action==2) {
462	ClpPrimalColumnDantzig dantzig;
463	models[iModel].setPrimalColumnPivotAlgorithm(dantzig);
464	} else if (action==3) {
465	ClpPrimalColumnSteepest steep(4);
466	models[iModel].setPrimalColumnPivotAlgorithm(steep);
467	} else if (action==4) {
468	ClpPrimalColumnSteepest steep(1);
469	models[iModel].setPrimalColumnPivotAlgorithm(steep);
470	} else if (action==5) {
471	ClpPrimalColumnSteepest steep(2);
472	models[iModel].setPrimalColumnPivotAlgorithm(steep);
473	} else if (action==6) {
474	ClpPrimalColumnSteepest steep(10);
475	models[iModel].setPrimalColumnPivotAlgorithm(steep);
476	}
477	break;
478	case SCALING:
479	models[iModel].scaling(action);
480	break;
481	case AUTOSCALE:
482	models[iModel].setAutomaticScaling(action!=0);
483	break;
484	case SPARSEFACTOR:
485	models[iModel].setSparseFactorization((1-action)!=0);
486	break;
487	case BIASLU:
488	models[iModel].factorization()->setBiasLU(action);
489	break;
490	case PERTURBATION:
491	if (action==0)
492	models[iModel].setPerturbation(50);
493	else
494	models[iModel].setPerturbation(100);
495	break;
496	case ERRORSALLOWED:
497	allowImportErrors = action;
498	break;
499	case INTPRINT:
500	printMode=action;
501	break;
502	case KEEPNAMES:
503	keepImportNames = 1-action;
504	break;
505	case PRESOLVE:
506	if (action==0)
507	preSolve = 5;
508	else if (action==1)
509	preSolve=0;
510	else if (action==2)
511	preSolve=10;
512	else
513	preSolveFile=true;
514	break;
515	case PFI:
516	models[iModel].factorization()->setForrestTomlin(action==0);
517	break;
518	case CRASH:
519	doCrash=action;
520	break;
521	case VECTOR:
522	doVector=action;
523	break;
524	case MESSAGES:
525	models[iModel].messageHandler()->setPrefix(action!=0);
526	break;
527	case CHOLESKY:
528	choleskyType = action;
529	break;
530	case GAMMA:
531	gamma=action;
532	break;
533	case BARRIERSCALE:
534	scaleBarrier=action;
535	break;
536	case KKT:
537	doKKT=action;
538	break;
539	case CROSSOVER:
540	crossover=action;
541	break;
542	default:
543	abort();
544	}
545	}
546	} else {
547	// action
548	if (type==EXIT)
549	break; // stop all
550	switch (type) {
551	case DUALSIMPLEX:
552	case PRIMALSIMPLEX:
553	case EITHERSIMPLEX:
554	case BARRIER:
555	// synonym for dual
556	case BAB:
557	if (goodModels[iModel]) {
558	double objScale =
559	parameters[whichParam(OBJSCALE2,numberParameters,parameters)].doubleValue();
560	if (objScale!=1.0) {
561	int iColumn;
562	int numberColumns=models[iModel].numberColumns();
563	double * dualColumnSolution =
564	models[iModel].dualColumnSolution();
565	ClpObjective * obj = models[iModel].objectiveAsObject();
566	assert(dynamic_cast<ClpLinearObjective *> (obj));
567	double offset;
568	double * objective = obj->gradient(NULL,NULL,offset,true);
569	for (iColumn=0;iColumn<numberColumns;iColumn++) {
570	dualColumnSolution[iColumn] *= objScale;
571	objective[iColumn] *= objScale;;
572	}
573	int iRow;
574	int numberRows=models[iModel].numberRows();
575	double * dualRowSolution =
576	models[iModel].dualRowSolution();
577	for (iRow=0;iRow<numberRows;iRow++)
578	dualRowSolution[iRow] *= objScale;
579	models[iModel].setObjectiveOffset(objScale*models[iModel].objectiveOffset());
580	}
581	ClpSolve::SolveType method;
582	ClpSolve::PresolveType presolveType;
583	ClpSimplex * model2 = models+iModel;
584	if (dualize) {
585	bool tryIt=true;
586	double fractionColumn=1.0;
587	double fractionRow=1.0;
588	if (dualize==3) {
589	dualize=1;
590	int numberColumns=model2->numberColumns();
591	int numberRows=model2->numberRows();
592	if (numberRows<50000\|\|5*numberColumns>numberRows) {
593	tryIt=false;
594	} else {
595	fractionColumn=0.1;
596	fractionRow=0.1;
597	}
598	}
599	if (tryIt) {
600	model2 = static_cast<ClpSimplexOther *> (model2)->dualOfModel(fractionRow,fractionColumn);
601	if (model2) {
602	printf("Dual of model has %d rows and %d columns\n",
603	model2->numberRows(),model2->numberColumns());
604	model2->setOptimizationDirection(1.0);
605	} else {
606	model2 = models+iModel;
607	dualize=0;
608	}
609	} else {
610	dualize=0;
611	}
612	}
613	ClpSolve solveOptions;
614	solveOptions.setPresolveActions(presolveOptions);
615	solveOptions.setSubstitution(substitution);
616	if (preSolve!=5&&preSolve) {
617	presolveType=ClpSolve::presolveNumber;
618	if (preSolve<0) {
619	preSolve = - preSolve;
620	if (preSolve<=100) {
621	presolveType=ClpSolve::presolveNumber;
622	printf("Doing %d presolve passes - picking up non-costed slacks\n",
623	preSolve);
624	solveOptions.setDoSingletonColumn(true);
625	} else {
626	preSolve -=100;
627	presolveType=ClpSolve::presolveNumberCost;
628	printf("Doing %d presolve passes - picking up costed slacks\n",
629	preSolve);
630	}
631	}
632	} else if (preSolve) {
633	presolveType=ClpSolve::presolveOn;
634	} else {
635	presolveType=ClpSolve::presolveOff;
636	}
637	solveOptions.setPresolveType(presolveType,preSolve);
638	if (type==DUALSIMPLEX\|\|type==BAB) {
639	method=ClpSolve::useDual;
640	} else if (type==PRIMALSIMPLEX) {
641	method=ClpSolve::usePrimalorSprint;
642	} else if (type==EITHERSIMPLEX) {
643	method=ClpSolve::automatic;
644	} else {
645	method = ClpSolve::useBarrier;
646	if (crossover==1) {
647	method=ClpSolve::useBarrierNoCross;
648	} else if (crossover==2) {
649	ClpObjective * obj = models[iModel].objectiveAsObject();
650	if (obj->type()>1) {
651	method=ClpSolve::useBarrierNoCross;
652	presolveType=ClpSolve::presolveOff;
653	solveOptions.setPresolveType(presolveType,preSolve);
654	}
655	}
656	}
657	solveOptions.setSolveType(method);
658	if(preSolveFile)
659	presolveOptions \|= 0x40000000;
660	solveOptions.setSpecialOption(4,presolveOptions);
661	solveOptions.setSpecialOption(5,printOptions&1);
662	if (doVector) {
663	ClpMatrixBase * matrix = models[iModel].clpMatrix();
664	if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
665	ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
666	clpMatrix->makeSpecialColumnCopy();
667	}
668	}
669	if (method==ClpSolve::useDual) {
670	// dual
671	if (doCrash)
672	solveOptions.setSpecialOption(0,1,doCrash); // crash
673	else if (doIdiot)
674	solveOptions.setSpecialOption(0,2,doIdiot); // possible idiot
675	} else if (method==ClpSolve::usePrimalorSprint) {
676	// primal
677	// if slp turn everything off
678	if (slpValue>0) {
679	doCrash=false;
680	doSprint=0;
681	doIdiot=-1;
682	solveOptions.setSpecialOption(1,10,slpValue); // slp
683	method=ClpSolve::usePrimal;
684	}
685	if (doCrash) {
686	solveOptions.setSpecialOption(1,1,doCrash); // crash
687	} else if (doSprint>0) {
688	// sprint overrides idiot
689	solveOptions.setSpecialOption(1,3,doSprint); // sprint
690	} else if (doIdiot>0) {
691	solveOptions.setSpecialOption(1,2,doIdiot); // idiot
692	} else if (slpValue<=0) {
693	if (doIdiot==0) {
694	if (doSprint==0)
695	solveOptions.setSpecialOption(1,4); // all slack
696	else
697	solveOptions.setSpecialOption(1,9); // all slack or sprint
698	} else {
699	if (doSprint==0)
700	solveOptions.setSpecialOption(1,8); // all slack or idiot
701	else
702	solveOptions.setSpecialOption(1,7); // initiative
703	}
704	}
705	if (basisHasValues==-1)
706	solveOptions.setSpecialOption(1,11); // switch off values
707	} else if (method==ClpSolve::useBarrier\|\|method==ClpSolve::useBarrierNoCross) {
708	int barrierOptions = choleskyType;
709	if (scaleBarrier)
710	barrierOptions \|= 8;
711	if (doKKT)
712	barrierOptions \|= 16;
713	if (gamma)
714	barrierOptions \|= 32*gamma;
715	if (crossover==3)
716	barrierOptions \|= 256; // try presolve in crossover
717	solveOptions.setSpecialOption(4,barrierOptions);
718	}
719	int status;
720	if (cppValue>=0) {
721	// generate code
722	FILE * fp = fopen("user_driver.cpp","w");
723	if (fp) {
724	// generate enough to do solveOptions
725	model2->generateCpp(fp);
726	solveOptions.generateCpp(fp);
727	fclose(fp);
728	// now call generate code
729	generateCode("user_driver.cpp",cppValue);
730	} else {
731	std::cout<<"Unable to open file user_driver.cpp"<<std::endl;
732	}
733	}
734	#ifdef CLP_MULTIPLE_FACTORIZATIONS
735	int denseCode = parameters[whichParam(DENSE,numberParameters,parameters)].intValue();
736	model2->factorization()->setGoDenseThreshold(denseCode);
737	int smallCode = parameters[whichParam(SMALLFACT,numberParameters,parameters)].intValue();
738	model2->factorization()->setGoSmallThreshold(smallCode);
739	model2->factorization()->goDenseOrSmall(model2->numberRows());
740	#endif
741	try {
742	status=model2->initialSolve(solveOptions);
743	}
744	catch (CoinError e) {
745	e.print();
746	status=-1;
747	}
748	if (dualize) {
749	int returnCode=static_cast<ClpSimplexOther *> (models+iModel)->restoreFromDual(model2);
750	if (model2->status()==3)
751	returnCode=0;
752	delete model2;
753	if (returnCode&&dualize!=2) {
754	currentModel = models+iModel;
755	// register signal handler
756	signal(SIGINT,signal_handler);
757	models[iModel].primal(1);
758	currentModel=NULL;
759	}
760	}
761	if (status>=0)
762	basisHasValues=1;
763	} else {
764	std::cout<<"** Current model not valid"<<std::endl;
765	}
766	break;
767	case STATISTICS:
768	if (goodModels[iModel]) {
769	// If presolve on look at presolved
770	bool deleteModel2=false;
771	ClpSimplex * model2 = models+iModel;
772	if (preSolve) {
773	ClpPresolve pinfo;
774	int presolveOptions2 = presolveOptions&~0x40000000;
775	if ((presolveOptions2&0xffff)!=0)
776	pinfo.setPresolveActions(presolveOptions2);
777	pinfo.setSubstitution(substitution);
778	if ((printOptions&1)!=0)
779	pinfo.statistics();
780	double presolveTolerance =
781	parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
782	model2 =
783	pinfo.presolvedModel(models[iModel],presolveTolerance,
784	true,preSolve);
785	if (model2) {
786	printf("Statistics for presolved model\n");
787	deleteModel2=true;
788	} else {
789	printf("Presolved model looks infeasible - will use unpresolved\n");
790	model2 = models+iModel;
791	}
792	} else {
793	printf("Statistics for unpresolved model\n");
794	model2 = models+iModel;
795	}
796	statistics(models+iModel,model2);
797	if (deleteModel2)
798	delete model2;
799	} else {
800	std::cout<<"** Current model not valid"<<std::endl;
801	}
802	break;
803	case TIGHTEN:
804	if (goodModels[iModel]) {
805	int numberInfeasibilities = models[iModel].tightenPrimalBounds();
806	if (numberInfeasibilities)
807	std::cout<<"** Analysis indicates model infeasible"<<std::endl;
808	} else {
809	std::cout<<"** Current model not valid"<<std::endl;
810	}
811	break;
812	case PLUSMINUS:
813	if (goodModels[iModel]) {
814	ClpMatrixBase * saveMatrix = models[iModel].clpMatrix();
815	ClpPackedMatrix* clpMatrix =
816	dynamic_cast< ClpPackedMatrix*>(saveMatrix);
817	if (clpMatrix) {
818	ClpPlusMinusOneMatrix * newMatrix = new ClpPlusMinusOneMatrix(*(clpMatrix->matrix()));
819	if (newMatrix->getIndices()) {
820	models[iModel].replaceMatrix(newMatrix);
821	delete saveMatrix;
822	std::cout<<"Matrix converted to +- one matrix"<<std::endl;
823	} else {
824	std::cout<<"Matrix can not be converted to +- 1 matrix"<<std::endl;
825	}
826	} else {
827	std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
828	}
829	} else {
830	std::cout<<"** Current model not valid"<<std::endl;
831	}
832	break;
833	case NETWORK:
834	if (goodModels[iModel]) {
835	ClpMatrixBase * saveMatrix = models[iModel].clpMatrix();
836	ClpPackedMatrix* clpMatrix =
837	dynamic_cast< ClpPackedMatrix*>(saveMatrix);
838	if (clpMatrix) {
839	ClpNetworkMatrix * newMatrix = new ClpNetworkMatrix(*(clpMatrix->matrix()));
840	if (newMatrix->getIndices()) {
841	models[iModel].replaceMatrix(newMatrix);
842	delete saveMatrix;
843	std::cout<<"Matrix converted to network matrix"<<std::endl;
844	} else {
845	std::cout<<"Matrix can not be converted to network matrix"<<std::endl;
846	}
847	} else {
848	std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
849	}
850	} else {
851	std::cout<<"** Current model not valid"<<std::endl;
852	}
853	break;
854	case IMPORT:
855	{
856	// get next field
857	field = CoinReadGetString(argc,argv);
858	if (field=="$") {
859	field = parameters[iParam].stringValue();
860	} else if (field=="EOL") {
861	parameters[iParam].printString();
862	break;
863	} else {
864	parameters[iParam].setStringValue(field);
865	}
866	std::string fileName;
867	bool canOpen=false;
868	// See if gmpl file
869	int gmpl=0;
870	std::string gmplData;
871	if (field=="-") {
872	// stdin
873	canOpen=true;
874	fileName = "-";
875	} else {
876	// See if .lp
877	{
878	const char * c_name = field.c_str();
879	int length = strlen(c_name);
880	if (length>3&&!strncmp(c_name+length-3,".lp",3))
881	gmpl=-1; // .lp
882	}
883	bool absolutePath;
884	if (dirsep=='/') {
885	// non Windows (or cygwin)
886	absolutePath=(field[0]=='/');
887	} else {
888	//Windows (non cycgwin)
889	absolutePath=(field[0]=='\\');
890	// but allow for :
891	if (strchr(field.c_str(),':'))
892	absolutePath=true;
893	}
894	if (absolutePath) {
895	fileName = field;
896	} else if (field[0]=='~') {
897	char * environVar = getenv("HOME");
898	if (environVar) {
899	std::string home(environVar);
900	field=field.erase(0,1);
901	fileName = home+field;
902	} else {
903	fileName=field;
904	}
905	} else {
906	fileName = directory+field;
907	// See if gmpl (model & data) - or even lp file
908	int length = field.size();
909	int percent = field.find('%');
910	if (percent<length&&percent>0) {
911	gmpl=1;
912	fileName = directory+field.substr(0,percent);
913	gmplData = directory+field.substr(percent+1);
914	if (percent<length-1)
915	gmpl=2; // two files
916	printf("GMPL model file %s and data file %s\n",
917	fileName.c_str(),gmplData.c_str());
918	}
919	}
920	std::string name=fileName;
921	if (fileCoinReadable(name)) {
922	// can open - lets go for it
923	canOpen=true;
924	if (gmpl==2) {
925	FILE *fp;
926	fp=fopen(gmplData.c_str(),"r");
927	if (fp) {
928	fclose(fp);
929	} else {
930	canOpen=false;
931	std::cout<<"Unable to open file "<<gmplData<<std::endl;
932	}
933	}
934	} else {
935	std::cout<<"Unable to open file "<<fileName<<std::endl;
936	}
937	}
938	if (canOpen) {
939	int status;
940	if (!gmpl)
941	status =models[iModel].readMps(fileName.c_str(),
942	keepImportNames!=0,
943	allowImportErrors!=0);
944	else if (gmpl>0)
945	status= models[iModel].readGMPL(fileName.c_str(),
946	(gmpl==2) ? gmplData.c_str() : NULL,
947	keepImportNames!=0);
948	else
949	status= models[iModel].readLp(fileName.c_str(),1.0e-12);
950	if (!status\|\|(status>0&&allowImportErrors)) {
951	goodModels[iModel]=true;
952	// sets to all slack (not necessary?)
953	models[iModel].createStatus();
954	time2 = CoinCpuTime();
955	totalTime += time2-time1;
956	time1=time2;
957	// Go to canned file if just input file
958	if (CbcOrClpRead_mode==2&&argc==2) {
959	// only if ends .mps
960	char * find = const_cast<char *>(strstr(fileName.c_str(),".mps"));
961	if (find&&find[4]=='\0') {
962	find[1]='p'; find[2]='a';find[3]='r';
963	FILE *fp=fopen(fileName.c_str(),"r");
964	if (fp) {
965	CbcOrClpReadCommand=fp; // Read from that file
966	CbcOrClpRead_mode=-1;
967	}
968	}
969	}
970	} else {
971	// errors
972	std::cout<<"There were "<<status<<
973	" errors on input"<<std::endl;
974	}
975	}
976	}
977	break;
978	case EXPORT:
979	if (goodModels[iModel]) {
980	double objScale =
981	parameters[whichParam(OBJSCALE2,numberParameters,parameters)].doubleValue();
982	if (objScale!=1.0) {
983	int iColumn;
984	int numberColumns=models[iModel].numberColumns();
985	double * dualColumnSolution =
986	models[iModel].dualColumnSolution();
987	ClpObjective * obj = models[iModel].objectiveAsObject();
988	assert(dynamic_cast<ClpLinearObjective *> (obj));
989	double offset;
990	double * objective = obj->gradient(NULL,NULL,offset,true);
991	for (iColumn=0;iColumn<numberColumns;iColumn++) {
992	dualColumnSolution[iColumn] *= objScale;
993	objective[iColumn] *= objScale;;
994	}
995	int iRow;
996	int numberRows=models[iModel].numberRows();
997	double * dualRowSolution =
998	models[iModel].dualRowSolution();
999	for (iRow=0;iRow<numberRows;iRow++)
1000	dualRowSolution[iRow] *= objScale;
1001	models[iModel].setObjectiveOffset(objScale*models[iModel].objectiveOffset());
1002	}
1003	// get next field
1004	field = CoinReadGetString(argc,argv);
1005	if (field=="$") {
1006	field = parameters[iParam].stringValue();
1007	} else if (field=="EOL") {
1008	parameters[iParam].printString();
1009	break;
1010	} else {
1011	parameters[iParam].setStringValue(field);
1012	}
1013	std::string fileName;
1014	bool canOpen=false;
1015	if (field[0]=='/'\|\|field[0]=='\\') {
1016	fileName = field;
1017	} else if (field[0]=='~') {
1018	char * environVar = getenv("HOME");
1019	if (environVar) {
1020	std::string home(environVar);
1021	field=field.erase(0,1);
1022	fileName = home+field;
1023	} else {
1024	fileName=field;
1025	}
1026	} else {
1027	fileName = directory+field;
1028	}
1029	FILE *fp=fopen(fileName.c_str(),"w");
1030	if (fp) {
1031	// can open - lets go for it
1032	fclose(fp);
1033	canOpen=true;
1034	} else {
1035	std::cout<<"Unable to open file "<<fileName<<std::endl;
1036	}
1037	if (canOpen) {
1038	// If presolve on then save presolved
1039	bool deleteModel2=false;
1040	ClpSimplex * model2 = models+iModel;
1041	if (dualize&&dualize<3) {
1042	model2 = static_cast<ClpSimplexOther *> (model2)->dualOfModel();
1043	printf("Dual of model has %d rows and %d columns\n",
1044	model2->numberRows(),model2->numberColumns());
1045	model2->setOptimizationDirection(1.0);
1046	preSolve=0; // as picks up from model
1047	}
1048	if (preSolve) {
1049	ClpPresolve pinfo;
1050	int presolveOptions2 = presolveOptions&~0x40000000;
1051	if ((presolveOptions2&0xffff)!=0)
1052	pinfo.setPresolveActions(presolveOptions2);
1053	pinfo.setSubstitution(substitution);
1054	if ((printOptions&1)!=0)
1055	pinfo.statistics();
1056	double presolveTolerance =
1057	parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
1058	model2 =
1059	pinfo.presolvedModel(models[iModel],presolveTolerance,
1060	true,preSolve,false,false);
1061	if (model2) {
1062	printf("Saving presolved model on %s\n",
1063	fileName.c_str());
1064	deleteModel2=true;
1065	} else {
1066	printf("Presolved model looks infeasible - saving original on %s\n",
1067	fileName.c_str());
1068	deleteModel2=false;
1069	model2 = models+iModel;
1070
1071	}
1072	} else {
1073	printf("Saving model on %s\n",
1074	fileName.c_str());
1075	}
1076	#if 0
1077	// Convert names
1078	int iRow;
1079	int numberRows=model2->numberRows();
1080	int iColumn;
1081	int numberColumns=model2->numberColumns();
1082
1083	char ** rowNames = NULL;
1084	char ** columnNames = NULL;
1085	if (model2->lengthNames()) {
1086	rowNames = new char * [numberRows];
1087	for (iRow=0;iRow<numberRows;iRow++) {
1088	rowNames[iRow] =
1089	CoinStrdup(model2->rowName(iRow).c_str());
1090	#ifdef STRIPBLANKS
1091	char * xx = rowNames[iRow];
1092	int i;
1093	int length = strlen(xx);
1094	int n=0;
1095	for (i=0;i<length;i++) {
1096	if (xx[i]!=' ')
1097	xx[n++]=xx[i];
1098	}
1099	xx[n]='\0';
1100	#endif
1101	}
1102
1103	columnNames = new char * [numberColumns];
1104	for (iColumn=0;iColumn<numberColumns;iColumn++) {
1105	columnNames[iColumn] =
1106	CoinStrdup(model2->columnName(iColumn).c_str());
1107	#ifdef STRIPBLANKS
1108	char * xx = columnNames[iColumn];
1109	int i;
1110	int length = strlen(xx);
1111	int n=0;
1112	for (i=0;i<length;i++) {
1113	if (xx[i]!=' ')
1114	xx[n++]=xx[i];
1115	}
1116	xx[n]='\0';
1117	#endif
1118	}
1119	}
1120	CoinMpsIO writer;
1121	writer.setMpsData(*model2->matrix(), COIN_DBL_MAX,
1122	model2->getColLower(), model2->getColUpper(),
1123	model2->getObjCoefficients(),
1124	(const char) 0 /integrality*/,
1125	model2->getRowLower(), model2->getRowUpper(),
1126	columnNames, rowNames);
1127	// Pass in array saying if each variable integer
1128	writer.copyInIntegerInformation(model2->integerInformation());
1129	writer.setObjectiveOffset(model2->objectiveOffset());
1130	writer.writeMps(fileName.c_str(),0,1,1);
1131	if (rowNames) {
1132	for (iRow=0;iRow<numberRows;iRow++) {
1133	free(rowNames[iRow]);
1134	}
1135	delete [] rowNames;
1136	for (iColumn=0;iColumn<numberColumns;iColumn++) {
1137	free(columnNames[iColumn]);
1138	}
1139	delete [] columnNames;
1140	}
1141	#else
1142	model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
1143	#endif
1144	if (deleteModel2)
1145	delete model2;
1146	time2 = CoinCpuTime();
1147	totalTime += time2-time1;
1148	time1=time2;
1149	}
1150	} else {
1151	std::cout<<"** Current model not valid"<<std::endl;
1152	}
1153	break;
1154	case BASISIN:
1155	if (goodModels[iModel]) {
1156	// get next field
1157	field = CoinReadGetString(argc,argv);
1158	if (field=="$") {
1159	field = parameters[iParam].stringValue();
1160	} else if (field=="EOL") {
1161	parameters[iParam].printString();
1162	break;
1163	} else {
1164	parameters[iParam].setStringValue(field);
1165	}
1166	std::string fileName;
1167	bool canOpen=false;
1168	if (field=="-") {
1169	// stdin
1170	canOpen=true;
1171	fileName = "-";
1172	} else {
1173	if (field[0]=='/'\|\|field[0]=='\\') {
1174	fileName = field;
1175	} else if (field[0]=='~') {
1176	char * environVar = getenv("HOME");
1177	if (environVar) {
1178	std::string home(environVar);
1179	field=field.erase(0,1);
1180	fileName = home+field;
1181	} else {
1182	fileName=field;
1183	}
1184	} else {
1185	fileName = directory+field;
1186	}
1187	FILE *fp=fopen(fileName.c_str(),"r");
1188	if (fp) {
1189	// can open - lets go for it
1190	fclose(fp);
1191	canOpen=true;
1192	} else {
1193	std::cout<<"Unable to open file "<<fileName<<std::endl;
1194	}
1195	}
1196	if (canOpen) {
1197	int values = models[iModel].readBasis(fileName.c_str());
1198	if (values==0)
1199	basisHasValues=-1;
1200	else
1201	basisHasValues=1;
1202	}
1203	} else {
1204	std::cout<<"** Current model not valid"<<std::endl;
1205	}
1206	break;
1207	case PRINTMASK:
1208	// get next field
1209	{
1210	std::string name = CoinReadGetString(argc,argv);
1211	if (name!="EOL") {
1212	parameters[iParam].setStringValue(name);
1213	printMask = name;
1214	} else {
1215	parameters[iParam].printString();
1216	}
1217	}
1218	break;
1219	case BASISOUT:
1220	if (goodModels[iModel]) {
1221	// get next field
1222	field = CoinReadGetString(argc,argv);
1223	if (field=="$") {
1224	field = parameters[iParam].stringValue();
1225	} else if (field=="EOL") {
1226	parameters[iParam].printString();
1227	break;
1228	} else {
1229	parameters[iParam].setStringValue(field);
1230	}
1231	std::string fileName;
1232	bool canOpen=false;
1233	if (field[0]=='/'\|\|field[0]=='\\') {
1234	fileName = field;
1235	} else if (field[0]=='~') {
1236	char * environVar = getenv("HOME");
1237	if (environVar) {
1238	std::string home(environVar);
1239	field=field.erase(0,1);
1240	fileName = home+field;
1241	} else {
1242	fileName=field;
1243	}
1244	} else {
1245	fileName = directory+field;
1246	}
1247	FILE *fp=fopen(fileName.c_str(),"w");
1248	if (fp) {
1249	// can open - lets go for it
1250	fclose(fp);
1251	canOpen=true;
1252	} else {
1253	std::cout<<"Unable to open file "<<fileName<<std::endl;
1254	}
1255	if (canOpen) {
1256	ClpSimplex * model2 = models+iModel;
1257	model2->writeBasis(fileName.c_str(),outputFormat>1,outputFormat-2);
1258	time2 = CoinCpuTime();
1259	totalTime += time2-time1;
1260	time1=time2;
1261	}
1262	} else {
1263	std::cout<<"** Current model not valid"<<std::endl;
1264	}
1265	break;
1266	case SAVE:
1267	{
1268	// get next field
1269	field = CoinReadGetString(argc,argv);
1270	if (field=="$") {
1271	field = parameters[iParam].stringValue();
1272	} else if (field=="EOL") {
1273	parameters[iParam].printString();
1274	break;
1275	} else {
1276	parameters[iParam].setStringValue(field);
1277	}
1278	std::string fileName;
1279	bool canOpen=false;
1280	if (field[0]=='/'\|\|field[0]=='\\') {
1281	fileName = field;
1282	} else if (field[0]=='~') {
1283	char * environVar = getenv("HOME");
1284	if (environVar) {
1285	std::string home(environVar);
1286	field=field.erase(0,1);
1287	fileName = home+field;
1288	} else {
1289	fileName=field;
1290	}
1291	} else {
1292	fileName = directory+field;
1293	}
1294	FILE *fp=fopen(fileName.c_str(),"wb");
1295	if (fp) {
1296	// can open - lets go for it
1297	fclose(fp);
1298	canOpen=true;
1299	} else {
1300	std::cout<<"Unable to open file "<<fileName<<std::endl;
1301	}
1302	if (canOpen) {
1303	int status;
1304	// If presolve on then save presolved
1305	bool deleteModel2=false;
1306	ClpSimplex * model2 = models+iModel;
1307	if (preSolve) {
1308	ClpPresolve pinfo;
1309	double presolveTolerance =
1310	parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
1311	model2 =
1312	pinfo.presolvedModel(models[iModel],presolveTolerance,
1313	false,preSolve);
1314	if (model2) {
1315	printf("Saving presolved model on %s\n",
1316	fileName.c_str());
1317	deleteModel2=true;
1318	} else {
1319	printf("Presolved model looks infeasible - saving original on %s\n",
1320	fileName.c_str());
1321	deleteModel2=false;
1322	model2 = models+iModel;
1323
1324	}
1325	} else {
1326	printf("Saving model on %s\n",
1327	fileName.c_str());
1328	}
1329	status =model2->saveModel(fileName.c_str());
1330	if (deleteModel2)
1331	delete model2;
1332	if (!status) {
1333	goodModels[iModel]=true;
1334	time2 = CoinCpuTime();
1335	totalTime += time2-time1;
1336	time1=time2;
1337	} else {
1338	// errors
1339	std::cout<<"There were errors on output"<<std::endl;
1340	}
1341	}
1342	}
1343	break;
1344	case RESTORE:
1345	{
1346	// get next field
1347	field = CoinReadGetString(argc,argv);
1348	if (field=="$") {
1349	field = parameters[iParam].stringValue();
1350	} else if (field=="EOL") {
1351	parameters[iParam].printString();
1352	break;
1353	} else {
1354	parameters[iParam].setStringValue(field);
1355	}
1356	std::string fileName;
1357	bool canOpen=false;
1358	if (field[0]=='/'\|\|field[0]=='\\') {
1359	fileName = field;
1360	} else if (field[0]=='~') {
1361	char * environVar = getenv("HOME");
1362	if (environVar) {
1363	std::string home(environVar);
1364	field=field.erase(0,1);
1365	fileName = home+field;
1366	} else {
1367	fileName=field;
1368	}
1369	} else {
1370	fileName = directory+field;
1371	}
1372	FILE *fp=fopen(fileName.c_str(),"rb");
1373	if (fp) {
1374	// can open - lets go for it
1375	fclose(fp);
1376	canOpen=true;
1377	} else {
1378	std::cout<<"Unable to open file "<<fileName<<std::endl;
1379	}
1380	if (canOpen) {
1381	int status =models[iModel].restoreModel(fileName.c_str());
1382	if (!status) {
1383	goodModels[iModel]=true;
1384	time2 = CoinCpuTime();
1385	totalTime += time2-time1;
1386	time1=time2;
1387	} else {
1388	// errors
1389	std::cout<<"There were errors on input"<<std::endl;
1390	}
1391	}
1392	}
1393	break;
1394	case MAXIMIZE:
1395	models[iModel].setOptimizationDirection(-1);
1396	break;
1397	case MINIMIZE:
1398	models[iModel].setOptimizationDirection(1);
1399	break;
1400	case ALLSLACK:
1401	models[iModel].allSlackBasis(true);
1402	break;
1403	case REVERSE:
1404	if (goodModels[iModel]) {
1405	int iColumn;
1406	int numberColumns=models[iModel].numberColumns();
1407	double * dualColumnSolution =
1408	models[iModel].dualColumnSolution();
1409	ClpObjective * obj = models[iModel].objectiveAsObject();
1410	assert(dynamic_cast<ClpLinearObjective *> (obj));
1411	double offset;
1412	double * objective = obj->gradient(NULL,NULL,offset,true);
1413	for (iColumn=0;iColumn<numberColumns;iColumn++) {
1414	dualColumnSolution[iColumn] = -dualColumnSolution[iColumn];
1415	objective[iColumn] = -objective[iColumn];
1416	}
1417	int iRow;
1418	int numberRows=models[iModel].numberRows();
1419	double * dualRowSolution =
1420	models[iModel].dualRowSolution();
1421	for (iRow=0;iRow<numberRows;iRow++) {
1422	dualRowSolution[iRow] = -dualRowSolution[iRow];
1423	}
1424	models[iModel].setObjectiveOffset(-models[iModel].objectiveOffset());
1425	}
1426	break;
1427	case DIRECTORY:
1428	{
1429	std::string name = CoinReadGetString(argc,argv);
1430	if (name!="EOL") {
1431	int length=name.length();
1432	if (name[length-1]==dirsep) {
1433	directory = name;
1434	} else {
1435	directory = name+dirsep;
1436	}
1437	parameters[iParam].setStringValue(directory);
1438	} else {
1439	parameters[iParam].printString();
1440	}
1441	}
1442	break;
1443	case DIRSAMPLE:
1444	{
1445	std::string name = CoinReadGetString(argc,argv);
1446	if (name!="EOL") {
1447	int length=name.length();
1448	if (name[length-1]==dirsep) {
1449	dirSample = name;
1450	} else {
1451	dirSample = name+dirsep;
1452	}
1453	parameters[iParam].setStringValue(dirSample);
1454	} else {
1455	parameters[iParam].printString();
1456	}
1457	}
1458	break;
1459	case DIRNETLIB:
1460	{
1461	std::string name = CoinReadGetString(argc,argv);
1462	if (name!="EOL") {
1463	int length=name.length();
1464	if (name[length-1]==dirsep) {
1465	dirNetlib = name;
1466	} else {
1467	dirNetlib = name+dirsep;
1468	}
1469	parameters[iParam].setStringValue(dirNetlib);
1470	} else {
1471	parameters[iParam].printString();
1472	}
1473	}
1474	break;
1475	case DIRMIPLIB:
1476	{
1477	std::string name = CoinReadGetString(argc,argv);
1478	if (name!="EOL") {
1479	int length=name.length();
1480	if (name[length-1]==dirsep) {
1481	dirMiplib = name;
1482	} else {
1483	dirMiplib = name+dirsep;
1484	}
1485	parameters[iParam].setStringValue(dirMiplib);
1486	} else {
1487	parameters[iParam].printString();
1488	}
1489	}
1490	break;
1491	case STDIN:
1492	CbcOrClpRead_mode=-1;
1493	break;
1494	case NETLIB_DUAL:
1495	case NETLIB_EITHER:
1496	case NETLIB_BARRIER:
1497	case NETLIB_PRIMAL:
1498	case NETLIB_TUNE:
1499	{
1500	// create fields for unitTest
1501	const char * fields[4];
1502	int nFields=4;
1503	fields[0]="fake main from unitTest";
1504	std::string mpsfield = "-dirSample=";
1505	mpsfield += dirSample.c_str();
1506	fields[1]=mpsfield.c_str();
1507	std::string netfield = "-dirNetlib=";
1508	netfield += dirNetlib.c_str();
1509	fields[2]=netfield.c_str();
1510	fields[3]="-netlib";
1511	int algorithm;
1512	if (type==NETLIB_DUAL) {
1513	std::cerr<<"Doing netlib with dual algorithm"<<std::endl;
1514	algorithm =0;
1515	} else if (type==NETLIB_BARRIER) {
1516	std::cerr<<"Doing netlib with barrier algorithm"<<std::endl;
1517	algorithm =2;
1518	} else if (type==NETLIB_EITHER) {
1519	std::cerr<<"Doing netlib with dual or primal algorithm"<<std::endl;
1520	algorithm =3;
1521	} else if (type==NETLIB_TUNE) {
1522	std::cerr<<"Doing netlib with best algorithm!"<<std::endl;
1523	algorithm =5;
1524	// uncomment next to get active tuning
1525	// algorithm=6;
1526	} else {
1527	std::cerr<<"Doing netlib with primal agorithm"<<std::endl;
1528	algorithm=1;
1529	}
1530	int specialOptions = models[iModel].specialOptions();
1531	models[iModel].setSpecialOptions(0);
1532	mainTest(nFields,fields,algorithm,models[iModel],
1533	(preSolve!=0),specialOptions,doVector!=0);
1534	}
1535	break;
1536	case UNITTEST:
1537	{
1538	// create fields for unitTest
1539	const char * fields[2];
1540	int nFields=2;
1541	fields[0]="fake main from unitTest";
1542	std::string dirfield = "-dirSample=";
1543	dirfield += dirSample.c_str();
1544	fields[1]=dirfield.c_str();
1545	int specialOptions = models[iModel].specialOptions();
1546	models[iModel].setSpecialOptions(0);
1547	int algorithm=-1;
1548	if (models[iModel].numberRows())
1549	algorithm=7;
1550	mainTest(nFields,fields,algorithm,models[iModel],(preSolve!=0),specialOptions,doVector!=0);
1551	}
1552	break;
1553	case FAKEBOUND:
1554	if (goodModels[iModel]) {
1555	// get bound
1556	double value = CoinReadGetDoubleField(argc,argv,&valid);
1557	if (!valid) {
1558	std::cout<<"Setting "<<parameters[iParam].name()<<
1559	" to DEBUG "<<value<<std::endl;
1560	int iRow;
1561	int numberRows=models[iModel].numberRows();
1562	double * rowLower = models[iModel].rowLower();
1563	double * rowUpper = models[iModel].rowUpper();
1564	for (iRow=0;iRow<numberRows;iRow++) {
1565	// leave free ones for now
1566	if (rowLower[iRow]>-1.0e20\|\|rowUpper[iRow]<1.0e20) {
1567	rowLower[iRow]=CoinMax(rowLower[iRow],-value);
1568	rowUpper[iRow]=CoinMin(rowUpper[iRow],value);
1569	}
1570	}
1571	int iColumn;
1572	int numberColumns=models[iModel].numberColumns();
1573	double * columnLower = models[iModel].columnLower();
1574	double * columnUpper = models[iModel].columnUpper();
1575	for (iColumn=0;iColumn<numberColumns;iColumn++) {
1576	// leave free ones for now
1577	if (columnLower[iColumn]>-1.0e20\|\|
1578	columnUpper[iColumn]<1.0e20) {
1579	columnLower[iColumn]=CoinMax(columnLower[iColumn],-value);
1580	columnUpper[iColumn]=CoinMin(columnUpper[iColumn],value);
1581	}
1582	}
1583	} else if (valid==1) {
1584	abort();
1585	} else {
1586	std::cout<<"enter value for "<<parameters[iParam].name()<<
1587	std::endl;
1588	}
1589	}
1590	break;
1591	case REALLY_SCALE:
1592	if (goodModels[iModel]) {
1593	ClpSimplex newModel(models[iModel],
1594	models[iModel].scalingFlag());
1595	printf("model really really scaled\n");
1596	models[iModel]=newModel;
1597	}
1598	break;
1599	case USERCLP:
1600	// Replace the sample code by whatever you want
1601	if (goodModels[iModel]) {
1602	ClpSimplex * thisModel = &models[iModel];
1603	printf("Dummy user code - model has %d rows and %d columns\n",
1604	thisModel->numberRows(),thisModel->numberColumns());
1605	}
1606	break;
1607	case HELP:
1608	std::cout<<"Coin LP version "<<CLPVERSION
1609	<<", build "<<__DATE__<<std::endl;
1610	std::cout<<"Non default values:-"<<std::endl;
1611	std::cout<<"Perturbation "<<models[0].perturbation()<<" (default 100)"
1612	<<std::endl;
1613	CoinReadPrintit(
1614	"Presolve being done with 5 passes\n\
1615	Dual steepest edge steep/partial on matrix shape and factorization density\n\
1616	Clpnnnn taken out of messages\n\
1617	If Factorization frequency default then done on size of matrix\n\n\
1618	(-)unitTest, (-)netlib or (-)netlibp will do standard tests\n\n\
1619	You can switch to interactive mode at any time so\n\
1620	clp watson.mps -scaling off -primalsimplex\nis the same as\n\
1621	clp watson.mps -\nscaling off\nprimalsimplex"
1622	);
1623	break;
1624	case SOLUTION:
1625	if (goodModels[iModel]) {
1626	// get next field
1627	field = CoinReadGetString(argc,argv);
1628	if (field=="$") {
1629	field = parameters[iParam].stringValue();
1630	} else if (field=="EOL") {
1631	parameters[iParam].printString();
1632	break;
1633	} else {
1634	parameters[iParam].setStringValue(field);
1635	}
1636	std::string fileName;
1637	FILE *fp=NULL;
1638	if (field=="-"\|\|field=="EOL"\|\|field=="stdout") {
1639	// stdout
1640	fp=stdout;
1641	fprintf(fp,"\n");
1642	} else if (field=="stderr") {
1643	// stderr
1644	fp=stderr;
1645	fprintf(fp,"\n");
1646	} else {
1647	if (field[0]=='/'\|\|field[0]=='\\') {
1648	fileName = field;
1649	} else if (field[0]=='~') {
1650	char * environVar = getenv("HOME");
1651	if (environVar) {
1652	std::string home(environVar);
1653	field=field.erase(0,1);
1654	fileName = home+field;
1655	} else {
1656	fileName=field;
1657	}
1658	} else {
1659	fileName = directory+field;
1660	}
1661	fp=fopen(fileName.c_str(),"w");
1662	}
1663	if (fp) {
1664	// Write solution header (suggested by Luigi Poderico)
1665	double objValue = models[iModel].getObjValue()*models[iModel].getObjSense();
1666	int iStat = models[iModel].status();
1667	if (iStat==0) {
1668	fprintf(fp, "optimal\n" );
1669	} else if (iStat==1) {
1670	// infeasible
1671	fprintf(fp, "infeasible\n" );
1672	} else if (iStat==2) {
1673	// unbounded
1674	fprintf(fp, "unbounded\n" );
1675	} else if (iStat==3) {
1676	fprintf(fp, "stopped on iterations or time\n" );
1677	} else if (iStat==4) {
1678	fprintf(fp, "stopped on difficulties\n" );
1679	} else if (iStat==5) {
1680	fprintf(fp, "stopped on ctrl-c\n" );
1681	} else {
1682	fprintf(fp, "status unknown\n" );
1683	}
1684	fprintf(fp, "Objective value %15.8g\n", objValue);
1685	// make fancy later on
1686	int iRow;
1687	int numberRows=models[iModel].numberRows();
1688	int lengthName = models[iModel].lengthNames(); // 0 if no names
1689	// in general I don't want to pass around massive
1690	// amounts of data but seems simpler here
1691	std::vector<std::string> rowNames =
1692	*(models[iModel].rowNames());
1693	std::vector<std::string> columnNames =
1694	*(models[iModel].columnNames());
1695
1696	double * dualRowSolution = models[iModel].dualRowSolution();
1697	double * primalRowSolution =
1698	models[iModel].primalRowSolution();
1699	double * rowLower = models[iModel].rowLower();
1700	double * rowUpper = models[iModel].rowUpper();
1701	double primalTolerance = models[iModel].primalTolerance();
1702	char format[6];
1703	sprintf(format,"%%-%ds",CoinMax(lengthName,8));
1704	bool doMask = (printMask!=""&&lengthName);
1705	int * maskStarts=NULL;
1706	int maxMasks=0;
1707	char ** masks =NULL;
1708	if (doMask) {
1709	int nAst =0;
1710	const char * pMask2 = printMask.c_str();
1711	char pMask[100];
1712	int iChar;
1713	int lengthMask = strlen(pMask2);
1714	assert (lengthMask<100);
1715	if (*pMask2=='"') {
1716	if (pMask2[lengthMask-1]!='"') {
1717	printf("mismatched \" in mask %s\n",pMask2);
1718	break;
1719	} else {
1720	strcpy(pMask,pMask2+1);
1721	*strchr(pMask,'"')='\0';
1722	}
1723	} else if (*pMask2=='\'') {
1724	if (pMask2[lengthMask-1]!='\'') {
1725	printf("mismatched ' in mask %s\n",pMask2);
1726	break;
1727	} else {
1728	strcpy(pMask,pMask2+1);
1729	*strchr(pMask,'\'')='\0';
1730	}
1731	} else {
1732	strcpy(pMask,pMask2);
1733	}
1734	if (lengthMask>lengthName) {
1735	printf("mask %s too long - skipping\n",pMask);
1736	break;
1737	}
1738	maxMasks = 1;
1739	for (iChar=0;iChar<lengthMask;iChar++) {
1740	if (pMask[iChar]=='*') {
1741	nAst++;
1742	maxMasks *= (lengthName+1);
1743	}
1744	}
1745	int nEntries = 1;
1746	maskStarts = new int[lengthName+2];
1747	masks = new char * [maxMasks];
1748	char ** newMasks = new char * [maxMasks];
1749	int i;
1750	for (i=0;i<maxMasks;i++) {
1751	masks[i] = new char[lengthName+1];
1752	newMasks[i] = new char[lengthName+1];
1753	}
1754	strcpy(masks[0],pMask);
1755	for (int iAst=0;iAst<nAst;iAst++) {
1756	int nOldEntries = nEntries;
1757	nEntries=0;
1758	for (int iEntry = 0;iEntry<nOldEntries;iEntry++) {
1759	char * oldMask = masks[iEntry];
1760	char * ast = strchr(oldMask,'*');
1761	assert (ast);
1762	int length = strlen(oldMask)-1;
1763	int nBefore = ast-oldMask;
1764	int nAfter = length-nBefore;
1765	// and add null
1766	nAfter++;
1767	for (int i=0;i<=lengthName-length;i++) {
1768	char * maskOut = newMasks[nEntries];
1769	CoinMemcpyN(oldMask,nBefore,maskOut);
1770	for (int k=0;k<i;k++)
1771	maskOut[k+nBefore]='?';
1772	CoinMemcpyN(ast+1,nAfter,maskOut+nBefore+i);
1773	nEntries++;
1774	assert (nEntries<=maxMasks);
1775	}
1776	}
1777	char ** temp = masks;
1778	masks = newMasks;
1779	newMasks = temp;
1780	}
1781	// Now extend and sort
1782	int * sort = new int[nEntries];
1783	for (i=0;i<nEntries;i++) {
1784	char * maskThis = masks[i];
1785	int length = strlen(maskThis);
1786	while (maskThis[length-1]==' ')
1787	length--;
1788	maskThis[length]='\0';
1789	sort[i]=length;
1790	}
1791	CoinSort_2(sort,sort+nEntries,masks);
1792	int lastLength=-1;
1793	for (i=0;i<nEntries;i++) {
1794	int length = sort[i];
1795	while (length>lastLength)
1796	maskStarts[++lastLength] = i;
1797	}
1798	maskStarts[++lastLength]=nEntries;
1799	delete [] sort;
1800	for (i=0;i<maxMasks;i++)
1801	delete [] newMasks[i];
1802	delete [] newMasks;
1803	}
1804	if (printMode>2) {
1805	for (iRow=0;iRow<numberRows;iRow++) {
1806	int type=printMode-3;
1807	if (primalRowSolution[iRow]>rowUpper[iRow]+primalTolerance\|\|
1808	primalRowSolution[iRow]<rowLower[iRow]-primalTolerance) {
1809	fprintf(fp,"** ");
1810	type=2;
1811	} else if (fabs(primalRowSolution[iRow])>1.0e-8) {
1812	type=1;
1813	} else if (numberRows<50) {
1814	type=3;
1815	}
1816	if (doMask&&!maskMatches(maskStarts,masks,rowNames[iRow]))
1817	type=0;
1818	if (type) {
1819	fprintf(fp,"%7d ",iRow);
1820	if (lengthName)
1821	fprintf(fp,format,rowNames[iRow].c_str());
1822	fprintf(fp,"%15.8g %15.8g\n",primalRowSolution[iRow],
1823	dualRowSolution[iRow]);
1824	}
1825	}
1826	}
1827	int iColumn;
1828	int numberColumns=models[iModel].numberColumns();
1829	double * dualColumnSolution =
1830	models[iModel].dualColumnSolution();
1831	double * primalColumnSolution =
1832	models[iModel].primalColumnSolution();
1833	double * columnLower = models[iModel].columnLower();
1834	double * columnUpper = models[iModel].columnUpper();
1835	for (iColumn=0;iColumn<numberColumns;iColumn++) {
1836	int type=(printMode>3) ? 1 : 0;
1837	if (primalColumnSolution[iColumn]>columnUpper[iColumn]+primalTolerance\|\|
1838	primalColumnSolution[iColumn]<columnLower[iColumn]-primalTolerance) {
1839	fprintf(fp,"** ");
1840	type=2;
1841	} else if (fabs(primalColumnSolution[iColumn])>1.0e-8) {
1842	type=1;
1843	} else if (numberColumns<50) {
1844	type=3;
1845	}
1846	if (doMask&&!maskMatches(maskStarts,masks,
1847	columnNames[iColumn]))
1848	type =0;
1849	if (type) {
1850	fprintf(fp,"%7d ",iColumn);
1851	if (lengthName)
1852	fprintf(fp,format,columnNames[iColumn].c_str());
1853	fprintf(fp,"%15.8g %15.8g\n",
1854	primalColumnSolution[iColumn],
1855	dualColumnSolution[iColumn]);
1856	}
1857	}
1858	if (fp!=stdout)
1859	fclose(fp);
1860	if (masks) {
1861	delete [] maskStarts;
1862	for (int i=0;i<maxMasks;i++)
1863	delete [] masks[i];
1864	delete [] masks;
1865	}
1866	} else {
1867	std::cout<<"Unable to open file "<<fileName<<std::endl;
1868	}
1869	} else {
1870	std::cout<<"** Current model not valid"<<std::endl;
1871
1872	}
1873
1874	break;
1875	case SAVESOL:
1876	if (goodModels[iModel]) {
1877	// get next field
1878	field = CoinReadGetString(argc,argv);
1879	if (field=="$") {
1880	field = parameters[iParam].stringValue();
1881	} else if (field=="EOL") {
1882	parameters[iParam].printString();
1883	break;
1884	} else {
1885	parameters[iParam].setStringValue(field);
1886	}
1887	std::string fileName;
1888	if (field[0]=='/'\|\|field[0]=='\\') {
1889	fileName = field;
1890	} else if (field[0]=='~') {
1891	char * environVar = getenv("HOME");
1892	if (environVar) {
1893	std::string home(environVar);
1894	field=field.erase(0,1);
1895	fileName = home+field;
1896	} else {
1897	fileName=field;
1898	}
1899	} else {
1900	fileName = directory+field;
1901	}
1902	saveSolution(models+iModel,fileName);
1903	} else {
1904	std::cout<<"** Current model not valid"<<std::endl;
1905
1906	}
1907	break;
1908	default:
1909	abort();
1910	}
1911	}
1912	} else if (!numberMatches) {
1913	std::cout<<"No match for "<<field<<" - ? for list of commands"
1914	<<std::endl;
1915	} else if (numberMatches==1) {
1916	if (!numberQuery) {
1917	std::cout<<"Short match for "<<field<<" - completion: ";
1918	std::cout<<parameters[firstMatch].matchName()<<std::endl;
1919	} else if (numberQuery) {
1920	std::cout<<parameters[firstMatch].matchName()<<" : ";
1921	std::cout<<parameters[firstMatch].shortHelp()<<std::endl;
1922	if (numberQuery>=2)
1923	parameters[firstMatch].printLongHelp();
1924	}
1925	} else {
1926	if (!numberQuery)
1927	std::cout<<"Multiple matches for "<<field<<" - possible completions:"
1928	<<std::endl;
1929	else
1930	std::cout<<"Completions of "<<field<<":"<<std::endl;
1931	for ( iParam=0; iParam<numberParameters; iParam++ ) {
1932	int match = parameters[iParam].matches(field);
1933	if (match&&parameters[iParam].displayThis()) {
1934	std::cout<<parameters[iParam].matchName();
1935	if (numberQuery>=2)
1936	std::cout<<" : "<<parameters[iParam].shortHelp();
1937	std::cout<<std::endl;
1938	}
1939	}
1940	}
1941	}
1942	delete [] models;
1943	delete [] goodModels;
1944	}
1945	// By now all memory should be freed
1946	#ifdef DMALLOC
1947	dmalloc_log_unfreed();
1948	dmalloc_shutdown();
1949	#endif
1950	return 0;
1951	}
1952	static void breakdown(const char * name, int numberLook, const double * region)
1953	{
1954	double range[] = {
1955	-COIN_DBL_MAX,
1956	-1.0e15,-1.0e11,-1.0e8,-1.0e5,-1.0e4,-1.0e3,-1.0e2,-1.0e1,
1957	-1.0,
1958	-1.0e-1,-1.0e-2,-1.0e-3,-1.0e-4,-1.0e-5,-1.0e-8,-1.0e-11,-1.0e-15,
1959	0.0,
1960	1.0e-15,1.0e-11,1.0e-8,1.0e-5,1.0e-4,1.0e-3,1.0e-2,1.0e-1,
1961	1.0,
1962	1.0e1,1.0e2,1.0e3,1.0e4,1.0e5,1.0e8,1.0e11,1.0e15,
1963	COIN_DBL_MAX};
1964	int nRanges = static_cast<int> (sizeof(range)/sizeof(double));
1965	int * number = new int[nRanges];
1966	memset(number,0,nRanges*sizeof(int));
1967	int * numberExact = new int[nRanges];
1968	memset(numberExact,0,nRanges*sizeof(int));
1969	int i;
1970	for ( i=0;i<numberLook;i++) {
1971	double value = region[i];
1972	for (int j=0;j<nRanges;j++) {
1973	if (value==range[j]) {
1974	numberExact[j]++;
1975	break;
1976	} else if (value<range[j]) {
1977	number[j]++;
1978	break;
1979	}
1980	}
1981	}
1982	printf("\n%s has %d entries\n",name,numberLook);
1983	for (i=0;i<nRanges;i++) {
1984	if (number[i])
1985	printf("%d between %g and %g",number[i],range[i-1],range[i]);
1986	if (numberExact[i]) {
1987	if (number[i])
1988	printf(", ");
1989	printf("%d exactly at %g",numberExact[i],range[i]);
1990	}
1991	if (number[i]+numberExact[i])
1992	printf("\n");
1993	}
1994	delete [] number;
1995	delete [] numberExact;
1996	}
1997	void sortOnOther(int * column,
1998	const CoinBigIndex * rowStart,
1999	int * order,
2000	int * other,
2001	int nRow,
2002	int nInRow,
2003	int where)
2004	{
2005	if (nRow<2\|\|where>=nInRow)
2006	return;
2007	// do initial sort
2008	int kRow;
2009	int iRow;
2010	for ( kRow=0;kRow<nRow;kRow++) {
2011	iRow = order[kRow];
2012	other[kRow]=column[rowStart[iRow]+where];
2013	}
2014	CoinSort_2(other,other+nRow,order);
2015	int first=0;
2016	iRow=order[0];
2017	int firstC=column[rowStart[iRow]+where];
2018	kRow=1;
2019	while (kRow<nRow) {
2020	int lastC=9999999;;
2021	for (;kRow<nRow+1;kRow++) {
2022	if (kRow<nRow) {
2023	iRow=order[kRow];
2024	lastC=column[rowStart[iRow]+where];
2025	} else {
2026	lastC=9999999;
2027	}
2028	if (lastC>firstC)
2029	break;
2030	}
2031	// sort
2032	sortOnOther(column,rowStart,order+first,other,kRow-first,
2033	nInRow,where+1);
2034	firstC=lastC;
2035	first=kRow;
2036	}
2037	}
2038	static void statistics(ClpSimplex * originalModel, ClpSimplex * model)
2039	{
2040	int numberColumns = originalModel->numberColumns();
2041	const char * integerInformation = originalModel->integerInformation();
2042	const double * columnLower = originalModel->columnLower();
2043	const double * columnUpper = originalModel->columnUpper();
2044	int numberIntegers=0;
2045	int numberBinary=0;
2046	int iRow,iColumn;
2047	if (integerInformation) {
2048	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2049	if (integerInformation[iColumn]) {
2050	if (columnUpper[iColumn]>columnLower[iColumn]) {
2051	numberIntegers++;
2052	if (columnUpper[iColumn]==0.0&&columnLower[iColumn]==1)
2053	numberBinary++;
2054	}
2055	}
2056	}
2057	}
2058	numberColumns = model->numberColumns();
2059	int numberRows = model->numberRows();
2060	columnLower = model->columnLower();
2061	columnUpper = model->columnUpper();
2062	const double * rowLower = model->rowLower();
2063	const double * rowUpper = model->rowUpper();
2064	const double * objective = model->objective();
2065	CoinPackedMatrix * matrix = model->matrix();
2066	CoinBigIndex numberElements = matrix->getNumElements();
2067	const int * columnLength = matrix->getVectorLengths();
2068	//const CoinBigIndex * columnStart = matrix->getVectorStarts();
2069	const double * elementByColumn = matrix->getElements();
2070	int * number = new int[numberRows+1];
2071	memset(number,0,(numberRows+1)*sizeof(int));
2072	int numberObjSingletons=0;
2073	/* cType
2074	0 0/inf, 1 0/up, 2 lo/inf, 3 lo/up, 4 free, 5 fix, 6 -inf/0, 7 -inf/up,
2075	8 0/1
2076	*/
2077	int cType[9];
2078	std::string cName[]={"0.0->inf,","0.0->up,","lo->inf,","lo->up,","free,","fixed,","-inf->0.0,",
2079	"-inf->up,","0.0->1.0"};
2080	int nObjective=0;
2081	memset(cType,0,sizeof(cType));
2082	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2083	int length=columnLength[iColumn];
2084	if (length==1&&objective[iColumn])
2085	numberObjSingletons++;
2086	number[length]++;
2087	if (objective[iColumn])
2088	nObjective++;
2089	if (columnLower[iColumn]>-1.0e20) {
2090	if (columnLower[iColumn]==0.0) {
2091	if (columnUpper[iColumn]>1.0e20)
2092	cType[0]++;
2093	else if (columnUpper[iColumn]==1.0)
2094	cType[8]++;
2095	else if (columnUpper[iColumn]==0.0)
2096	cType[5]++;
2097	else
2098	cType[1]++;
2099	} else {
2100	if (columnUpper[iColumn]>1.0e20)
2101	cType[2]++;
2102	else if (columnUpper[iColumn]==columnLower[iColumn])
2103	cType[5]++;
2104	else
2105	cType[3]++;
2106	}
2107	} else {
2108	if (columnUpper[iColumn]>1.0e20)
2109	cType[4]++;
2110	else if (columnUpper[iColumn]==0.0)
2111	cType[6]++;
2112	else
2113	cType[7]++;
2114	}
2115	}
2116	/* rType
2117	0 E 0, 1 E 1, 2 E -1, 3 E other, 4 G 0, 5 G 1, 6 G other,
2118	7 L 0, 8 L 1, 9 L other, 10 Range 0/1, 11 Range other, 12 free
2119	*/
2120	int rType[13];
2121	std::string rName[]={"E 0.0,","E 1.0,","E -1.0,","E other,","G 0.0,","G 1.0,","G other,",
2122	"L 0.0,","L 1.0,","L other,","Range 0.0->1.0,","Range other,","Free"};
2123	memset(rType,0,sizeof(rType));
2124	for (iRow=0;iRow<numberRows;iRow++) {
2125	if (rowLower[iRow]>-1.0e20) {
2126	if (rowLower[iRow]==0.0) {
2127	if (rowUpper[iRow]>1.0e20)
2128	rType[4]++;
2129	else if (rowUpper[iRow]==1.0)
2130	rType[10]++;
2131	else if (rowUpper[iRow]==0.0)
2132	rType[0]++;
2133	else
2134	rType[11]++;
2135	} else if (rowLower[iRow]==1.0) {
2136	if (rowUpper[iRow]>1.0e20)
2137	rType[5]++;
2138	else if (rowUpper[iRow]==rowLower[iRow])
2139	rType[1]++;
2140	else
2141	rType[11]++;
2142	} else if (rowLower[iRow]==-1.0) {
2143	if (rowUpper[iRow]>1.0e20)
2144	rType[6]++;
2145	else if (rowUpper[iRow]==rowLower[iRow])
2146	rType[2]++;
2147	else
2148	rType[11]++;
2149	} else {
2150	if (rowUpper[iRow]>1.0e20)
2151	rType[6]++;
2152	else if (rowUpper[iRow]==rowLower[iRow])
2153	rType[3]++;
2154	else
2155	rType[11]++;
2156	}
2157	} else {
2158	if (rowUpper[iRow]>1.0e20)
2159	rType[12]++;
2160	else if (rowUpper[iRow]==0.0)
2161	rType[7]++;
2162	else if (rowUpper[iRow]==1.0)
2163	rType[8]++;
2164	else
2165	rType[9]++;
2166	}
2167	}
2168	// Basic statistics
2169	printf("\n\nProblem has %d rows, %d columns (%d with objective) and %d elements\n",
2170	numberRows,numberColumns,nObjective,numberElements);
2171	if (number[0]+number[1]) {
2172	printf("There are ");
2173	if (numberObjSingletons)
2174	printf("%d singletons with objective ",numberObjSingletons);
2175	int numberNoObj = number[1]-numberObjSingletons;
2176	if (numberNoObj)
2177	printf("%d singletons with no objective ",numberNoObj);
2178	if (number[0])
2179	printf("** %d columns have no entries",number[0]);
2180	printf("\n");
2181	}
2182	printf("Column breakdown:\n");
2183	int k;
2184	for (k=0;k<static_cast<int> (sizeof(cType)/sizeof(int));k++) {
2185	printf("%d of type %s ",cType[k],cName[k].c_str());
2186	if (((k+1)%3)==0)
2187	printf("\n");
2188	}
2189	if ((k%3)!=0)
2190	printf("\n");
2191	printf("Row breakdown:\n");
2192	for (k=0;k<static_cast<int> (sizeof(rType)/sizeof(int));k++) {
2193	printf("%d of type %s ",rType[k],rName[k].c_str());
2194	if (((k+1)%3)==0)
2195	printf("\n");
2196	}
2197	if ((k%3)!=0)
2198	printf("\n");
2199	#define SYM
2200	#ifndef SYM
2201	if (model->logLevel()<2)
2202	return ;
2203	#endif
2204	int kMax = model->logLevel()>3 ? 1000000 : 10;
2205	k=0;
2206	for (iRow=1;iRow<=numberRows;iRow++) {
2207	if (number[iRow]) {
2208	k++;
2209	printf("%d columns have %d entries\n",number[iRow],iRow);
2210	if (k==kMax)
2211	break;
2212	}
2213	}
2214	if (k<numberRows) {
2215	int kk=k;
2216	k=0;
2217	for (iRow=numberRows;iRow>=1;iRow--) {
2218	if (number[iRow]) {
2219	k++;
2220	if (k==kMax)
2221	break;
2222	}
2223	}
2224	if (k>kk) {
2225	printf("\n .........\n\n");
2226	iRow=k;
2227	k=0;
2228	for (;iRow<numberRows;iRow++) {
2229	if (number[iRow]) {
2230	k++;
2231	printf("%d columns have %d entries\n",number[iRow],iRow);
2232	if (k==kMax)
2233	break;
2234	}
2235	}
2236	}
2237	}
2238	delete [] number;
2239	printf("\n\n");
2240	if (model->logLevel()==63
2241	#ifdef SYM
2242	\|\|true
2243	#endif
2244	) {
2245	// get column copy
2246	CoinPackedMatrix columnCopy = *matrix;
2247	const int * columnLength = columnCopy.getVectorLengths();
2248	number = new int[numberRows+1];
2249	memset(number,0,(numberRows+1)*sizeof(int));
2250	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2251	int length=columnLength[iColumn];
2252	number[length]++;
2253	}
2254	k=0;
2255	for (iRow=1;iRow<=numberRows;iRow++) {
2256	if (number[iRow]) {
2257	k++;
2258	}
2259	}
2260	int * row = columnCopy.getMutableIndices();
2261	const CoinBigIndex * columnStart = columnCopy.getVectorStarts();
2262	double * element = columnCopy.getMutableElements();
2263	int * order = new int[numberColumns];
2264	int * other = new int[numberColumns];
2265	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2266	int length=columnLength[iColumn];
2267	order[iColumn]=iColumn;
2268	other[iColumn]=length;
2269	CoinBigIndex start = columnStart[iColumn];
2270	CoinSort_2(row+start,row+start+length,element+start);
2271	}
2272	CoinSort_2(other,other+numberColumns,order);
2273	int jColumn=number[0]+number[1];
2274	for (iRow=2;iRow<=numberRows;iRow++) {
2275	if (number[iRow]) {
2276	printf("XX %d columns have %d entries\n",number[iRow],iRow);
2277	int kColumn = jColumn+number[iRow];
2278	sortOnOther(row,columnStart,
2279	order+jColumn,other,number[iRow],iRow,0);
2280	// Now print etc
2281	if (iRow<500000) {
2282	for (int lColumn =jColumn;lColumn<kColumn;lColumn++) {
2283	iColumn = order[lColumn];
2284	CoinBigIndex start = columnStart[iColumn];
2285	if (model->logLevel()==63) {
2286	printf("column %d %g <= ",iColumn,columnLower[iColumn]);
2287	for (CoinBigIndex i=start;i<start+iRow;i++)
2288	printf("( %d, %g) ",row[i],element[i]);
2289	printf("<= %g\n",columnUpper[iColumn]);
2290	}
2291	}
2292	}
2293	jColumn =kColumn;
2294	}
2295	}
2296	delete [] order;
2297	delete [] other;
2298	delete [] number;
2299	}
2300	// get row copy
2301	CoinPackedMatrix rowCopy = *matrix;
2302	rowCopy.reverseOrdering();
2303	const int * rowLength = rowCopy.getVectorLengths();
2304	number = new int[numberColumns+1];
2305	memset(number,0,(numberColumns+1)*sizeof(int));
2306	for (iRow=0;iRow<numberRows;iRow++) {
2307	int length=rowLength[iRow];
2308	number[length]++;
2309	}
2310	if (number[0])
2311	printf("** %d rows have no entries\n",number[0]);
2312	k=0;
2313	for (iColumn=1;iColumn<=numberColumns;iColumn++) {
2314	if (number[iColumn]) {
2315	k++;
2316	printf("%d rows have %d entries\n",number[iColumn],iColumn);
2317	if (k==kMax)
2318	break;
2319	}
2320	}
2321	if (k<numberColumns) {
2322	int kk=k;
2323	k=0;
2324	for (iColumn=numberColumns;iColumn>=1;iColumn--) {
2325	if (number[iColumn]) {
2326	k++;
2327	if (k==kMax)
2328	break;
2329	}
2330	}
2331	if (k>kk) {
2332	printf("\n .........\n\n");
2333	iColumn=k;
2334	k=0;
2335	for (;iColumn<numberColumns;iColumn++) {
2336	if (number[iColumn]) {
2337	k++;
2338	printf("%d rows have %d entries\n",number[iColumn],iColumn);
2339	if (k==kMax)
2340	break;
2341	}
2342	}
2343	}
2344	}
2345	if (model->logLevel()==63
2346	#ifdef SYM
2347	\|\|true
2348	#endif
2349	) {
2350	int * column = rowCopy.getMutableIndices();
2351	const CoinBigIndex * rowStart = rowCopy.getVectorStarts();
2352	double * element = rowCopy.getMutableElements();
2353	int * order = new int[numberRows];
2354	int * other = new int[numberRows];
2355	for (iRow=0;iRow<numberRows;iRow++) {
2356	int length=rowLength[iRow];
2357	order[iRow]=iRow;
2358	other[iRow]=length;
2359	CoinBigIndex start = rowStart[iRow];
2360	CoinSort_2(column+start,column+start+length,element+start);
2361	}
2362	CoinSort_2(other,other+numberRows,order);
2363	int jRow=number[0]+number[1];
2364	double * weight = new double[numberRows];
2365	double * randomColumn = new double[numberColumns+1];
2366	double * randomRow = new double [numberRows+1];
2367	int * sortRow = new int [numberRows];
2368	int * possibleRow = new int [numberRows];
2369	int * backRow = new int [numberRows];
2370	int * stackRow = new int [numberRows];
2371	int * sortColumn = new int [numberColumns];
2372	int * possibleColumn = new int [numberColumns];
2373	int * backColumn = new int [numberColumns];
2374	int * backColumn2 = new int [numberColumns];
2375	int * mapRow = new int [numberRows];
2376	int * mapColumn = new int [numberColumns];
2377	int * stackColumn = new int [numberColumns];
2378	double randomLower = CoinDrand48();
2379	double randomUpper = CoinDrand48();
2380	double randomInteger = CoinDrand48();
2381	int * startAdd = new int[numberRows+1];
2382	int * columnAdd = new int [2*numberElements];
2383	double * elementAdd = new double[2*numberElements];
2384	int nAddRows=0;
2385	startAdd[0]=0;
2386	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2387	randomColumn[iColumn] = CoinDrand48();
2388	backColumn2[iColumn]=-1;
2389	}
2390	for (iColumn=2;iColumn<=numberColumns;iColumn++) {
2391	if (number[iColumn]) {
2392	printf("XX %d rows have %d entries\n",number[iColumn],iColumn);
2393	int kRow = jRow+number[iColumn];
2394	sortOnOther(column,rowStart,
2395	order+jRow,other,number[iColumn],iColumn,0);
2396	// Now print etc
2397	if (iColumn<500000) {
2398	int nLook=0;
2399	for (int lRow =jRow;lRow<kRow;lRow++) {
2400	iRow = order[lRow];
2401	CoinBigIndex start = rowStart[iRow];
2402	if (model->logLevel()==63) {
2403	printf("row %d %g <= ",iRow,rowLower[iRow]);
2404	for (CoinBigIndex i=start;i<start+iColumn;i++)
2405	printf("( %d, %g) ",column[i],element[i]);
2406	printf("<= %g\n",rowUpper[iRow]);
2407	}
2408	int first = column[start];
2409	double sum=0.0;
2410	for (CoinBigIndex i=start;i<start+iColumn;i++) {
2411	int jColumn = column[i];
2412	double value = element[i];
2413	jColumn -= first;
2414	assert (jColumn>=0);
2415	sum += value*randomColumn[jColumn];
2416	}
2417	if (rowLower[iRow]>-1.0e30&&rowLower[iRow])
2418	sum += rowLower[iRow]*randomLower;
2419	else if (!rowLower[iRow])
2420	sum += 1.234567e-7*randomLower;
2421	if (rowUpper[iRow]<1.0e30&&rowUpper[iRow])
2422	sum += rowUpper[iRow]*randomUpper;
2423	else if (!rowUpper[iRow])
2424	sum += 1.234567e-7*randomUpper;
2425	sortRow[nLook]=iRow;
2426	randomRow[nLook++]=sum;
2427	// best way is to number unique elements and bounds and use
2428	if (fabs(sum)>1.0e4)
2429	sum *= 1.0e-6;
2430	weight[iRow]=sum;
2431	}
2432	assert (nLook<=numberRows);
2433	CoinSort_2(randomRow,randomRow+nLook,sortRow);
2434	randomRow[nLook]=COIN_DBL_MAX;
2435	double last=-COIN_DBL_MAX;
2436	int iLast=-1;
2437	for (int iLook=0;iLook<nLook+1;iLook++) {
2438	if (randomRow[iLook]>last) {
2439	if (iLast>=0) {
2440	int n=iLook-iLast;
2441	if (n>1) {
2442	//printf("%d rows possible?\n",n);
2443	}
2444	}
2445	iLast=iLook;
2446	last=randomRow[iLook];
2447	}
2448	}
2449	}
2450	jRow =kRow;
2451	}
2452	}
2453	CoinPackedMatrix columnCopy = *matrix;
2454	const int * columnLength = columnCopy.getVectorLengths();
2455	const int * row = columnCopy.getIndices();
2456	const CoinBigIndex * columnStart = columnCopy.getVectorStarts();
2457	const double * elementByColumn = columnCopy.getElements();
2458	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2459	int length=columnLength[iColumn];
2460	CoinBigIndex start = columnStart[iColumn];
2461	double sum = objective[iColumn];
2462	if (columnLower[iColumn]>-1.0e30&&columnLower[iColumn])
2463	sum += columnLower[iColumn]*randomLower;
2464	else if (!columnLower[iColumn])
2465	sum += 1.234567e-7*randomLower;
2466	if (columnUpper[iColumn]<1.0e30&&columnUpper[iColumn])
2467	sum += columnUpper[iColumn]*randomUpper;
2468	else if (!columnUpper[iColumn])
2469	sum += 1.234567e-7*randomUpper;
2470	if (model->isInteger(iColumn))
2471	sum += 9.87654321e-6*randomInteger;
2472	for (CoinBigIndex i=start;i<start+length;i++) {
2473	int iRow = row[i];
2474	sum += elementByColumn[i]*weight[iRow];
2475	}
2476	sortColumn[iColumn]=iColumn;
2477	randomColumn[iColumn]=sum;
2478	}
2479	{
2480	CoinSort_2(randomColumn,randomColumn+numberColumns,sortColumn);
2481	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2482	int i=sortColumn[iColumn];
2483	backColumn[i]=iColumn;
2484	}
2485	randomColumn[numberColumns]=COIN_DBL_MAX;
2486	double last=-COIN_DBL_MAX;
2487	int iLast=-1;
2488	for (int iLook=0;iLook<numberColumns+1;iLook++) {
2489	if (randomColumn[iLook]>last) {
2490	if (iLast>=0) {
2491	int n=iLook-iLast;
2492	if (n>1) {
2493	//printf("%d columns possible?\n",n);
2494	}
2495	for (int i=iLast;i<iLook;i++) {
2496	possibleColumn[sortColumn[i]]=n;
2497	}
2498	}
2499	iLast=iLook;
2500	last=randomColumn[iLook];
2501	}
2502	}
2503	for (iRow =0;iRow<numberRows;iRow++) {
2504	CoinBigIndex start = rowStart[iRow];
2505	double sum=0.0;
2506	int length=rowLength[iRow];
2507	for (CoinBigIndex i=start;i<start+length;i++) {
2508	int jColumn = column[i];
2509	double value = element[i];
2510	jColumn=backColumn[jColumn];
2511	sum += value*randomColumn[jColumn];
2512	//if (iColumn==23089\|\|iRow==23729)
2513	//printf("row %d cola %d colb %d value %g rand %g sum %g\n",
2514	// iRow,jColumn,column[i],value,randomColumn[jColumn],sum);
2515	}
2516	sortRow[iRow]=iRow;
2517	randomRow[iRow]=weight[iRow];
2518	randomRow[iRow]=sum;
2519	}
2520	CoinSort_2(randomRow,randomRow+numberRows,sortRow);
2521	for (iRow=0;iRow<numberRows;iRow++) {
2522	int i=sortRow[iRow];
2523	backRow[i]=iRow;
2524	}
2525	randomRow[numberRows]=COIN_DBL_MAX;
2526	last=-COIN_DBL_MAX;
2527	iLast=-1;
2528	// Do backward indices from order
2529	for (iRow=0;iRow<numberRows;iRow++) {
2530	other[order[iRow]]=iRow;
2531	}
2532	for (int iLook=0;iLook<numberRows+1;iLook++) {
2533	if (randomRow[iLook]>last) {
2534	if (iLast>=0) {
2535	int n=iLook-iLast;
2536	if (n>1) {
2537	//printf("%d rows possible?\n",n);
2538	// Within group sort as for original "order"
2539	for (int i=iLast;i<iLook;i++) {
2540	int jRow=sortRow[i];
2541	order[i]=other[jRow];
2542	}
2543	CoinSort_2(order+iLast,order+iLook,sortRow+iLast);
2544	}
2545	for (int i=iLast;i<iLook;i++) {
2546	possibleRow[sortRow[i]]=n;
2547	}
2548	}
2549	iLast=iLook;
2550	last=randomRow[iLook];
2551	}
2552	}
2553	// Temp out
2554	for (int iLook=0;iLook<numberRows-1000000;iLook++) {
2555	iRow=sortRow[iLook];
2556	CoinBigIndex start = rowStart[iRow];
2557	int length=rowLength[iRow];
2558	int numberPossible = possibleRow[iRow];
2559	for (CoinBigIndex i=start;i<start+length;i++) {
2560	int jColumn = column[i];
2561	if (possibleColumn[jColumn]!=numberPossible)
2562	numberPossible=-1;
2563	}
2564	int n=numberPossible;
2565	if (numberPossible>1) {
2566	//printf("pppppossible %d\n",numberPossible);
2567	for (int jLook=iLook+1;jLook<iLook+numberPossible;jLook++) {
2568	int jRow=sortRow[jLook];
2569	CoinBigIndex start2 = rowStart[jRow];
2570	assert (numberPossible==possibleRow[jRow]);
2571	assert(length==rowLength[jRow]);
2572	for (CoinBigIndex i=start2;i<start2+length;i++) {
2573	int jColumn = column[i];
2574	if (possibleColumn[jColumn]!=numberPossible)
2575	numberPossible=-1;
2576	}
2577	}
2578	if (numberPossible<2) {
2579	// switch off
2580	for (int jLook=iLook;jLook<iLook+n;jLook++)
2581	possibleRow[sortRow[jLook]]=-1;
2582	}
2583	// skip rest
2584	iLook += n-1;
2585	} else {
2586	possibleRow[iRow]=-1;
2587	}
2588	}
2589	for (int iLook=0;iLook<numberRows;iLook++) {
2590	iRow=sortRow[iLook];
2591	int numberPossible = possibleRow[iRow];
2592	// Only if any integers
2593	int numberIntegers=0;
2594	CoinBigIndex start = rowStart[iRow];
2595	int length=rowLength[iRow];
2596	for (CoinBigIndex i=start;i<start+length;i++) {
2597	int jColumn = column[i];
2598	if (model->isInteger(jColumn))
2599	numberIntegers++;
2600	}
2601	if (numberPossible>1&&!numberIntegers) {
2602	//printf("possible %d - but no integers\n",numberPossible);
2603	}
2604	if (numberPossible>1&&(numberIntegers\|\|true)) {
2605	//
2606	printf("possible %d - %d integers\n",numberPossible,numberIntegers);
2607	int lastLook=iLook;
2608	int nMapRow=-1;
2609	for (int jLook=iLook+1;jLook<iLook+numberPossible;jLook++) {
2610	// stop if too many failures
2611	if (jLook>iLook+10&&nMapRow<0)
2612	break;
2613	// Create identity mapping
2614	int i;
2615	for (i=0;i<numberRows;i++)
2616	mapRow[i]=i;
2617	for (i=0;i<numberColumns;i++)
2618	mapColumn[i]=i;
2619	int offset=jLook-iLook;
2620	int nStackC=0;
2621	// build up row and column mapping
2622	int nStackR=1;
2623	stackRow[0]=iLook;
2624	bool good=true;
2625	while (nStackR) {
2626	nStackR--;
2627	int look1 = stackRow[nStackR];
2628	int look2 = look1 + offset;
2629	assert (randomRow[look1]==randomRow[look2]);
2630	int row1=sortRow[look1];
2631	int row2=sortRow[look2];
2632	assert (mapRow[row1]==row1);
2633	assert (mapRow[row2]==row2);
2634	mapRow[row1]=row2;
2635	mapRow[row2]=row1;
2636	CoinBigIndex start1 = rowStart[row1];
2637	CoinBigIndex offset2 = rowStart[row2]-start1;
2638	int length=rowLength[row1];
2639	assert( length==rowLength[row2]);
2640	for (CoinBigIndex i=start1;i<start1+length;i++) {
2641	int jColumn1 = column[i];
2642	int jColumn2 = column[i+offset2];
2643	if (randomColumn[backColumn[jColumn1]]!=
2644	randomColumn[backColumn[jColumn2]]) {
2645	good=false;
2646	break;
2647	}
2648	if (mapColumn[jColumn1]==jColumn1) {
2649	// not touched
2650	assert (mapColumn[jColumn2]==jColumn2);
2651	if (jColumn1!=jColumn2) {
2652	// Put on stack
2653	mapColumn[jColumn1]=jColumn2;
2654	mapColumn[jColumn2]=jColumn1;
2655	stackColumn[nStackC++]=jColumn1;
2656	}
2657	} else {
2658	if(mapColumn[jColumn1]!=jColumn2\|\|
2659	mapColumn[jColumn2]!=jColumn1) {
2660	// bad
2661	good=false;
2662	printf("bad col\n");
2663	break;
2664	}
2665	}
2666	}
2667	if (!good)
2668	break;
2669	while (nStackC) {
2670	nStackC--;
2671	int iColumn = stackColumn[nStackC];
2672	int iColumn2 = mapColumn[iColumn];
2673	assert (iColumn!=iColumn2);
2674	int length=columnLength[iColumn];
2675	assert (length==columnLength[iColumn2]);
2676	CoinBigIndex start = columnStart[iColumn];
2677	CoinBigIndex offset2 = columnStart[iColumn2]-start;
2678	for (CoinBigIndex i=start;i<start+length;i++) {
2679	int iRow = row[i];
2680	int iRow2 = row[i+offset2];
2681	if (mapRow[iRow]==iRow) {
2682	// First (but be careful)
2683	if (iRow!=iRow2) {
2684	//mapRow[iRow]=iRow2;
2685	//mapRow[iRow2]=iRow;
2686	int iBack=backRow[iRow];
2687	int iBack2=backRow[iRow2];
2688	if (randomRow[iBack]==randomRow[iBack2]&&
2689	iBack2-iBack==offset) {
2690	stackRow[nStackR++]=iBack;
2691	} else {
2692	//printf("randomRow diff - weights %g %g\n",
2693	// weight[iRow],weight[iRow2]);
2694	// bad
2695	good=false;
2696	break;
2697	}
2698	}
2699	} else {
2700	if(mapRow[iRow]!=iRow2\|\|
2701	mapRow[iRow2]!=iRow) {
2702	// bad
2703	good=false;
2704	printf("bad row\n");
2705	break;
2706	}
2707	}
2708	}
2709	if (!good)
2710	break;
2711	}
2712	}
2713	// then check OK
2714	if (good) {
2715	for (iRow =0;iRow<numberRows;iRow++) {
2716	CoinBigIndex start = rowStart[iRow];
2717	int length=rowLength[iRow];
2718	if (mapRow[iRow]==iRow) {
2719	for (CoinBigIndex i=start;i<start+length;i++) {
2720	int jColumn = column[i];
2721	backColumn2[jColumn]=i-start;
2722	}
2723	for (CoinBigIndex i=start;i<start+length;i++) {
2724	int jColumn = column[i];
2725	if (mapColumn[jColumn]!=jColumn) {
2726	int jColumn2 =mapColumn[jColumn];
2727	CoinBigIndex i2 = backColumn2[jColumn2];
2728	if (i2<0) {
2729	good=false;
2730	} else if (element[i]!=element[i2+start]) {
2731	good=false;
2732	}
2733	}
2734	}
2735	for (CoinBigIndex i=start;i<start+length;i++) {
2736	int jColumn = column[i];
2737	backColumn2[jColumn]=-1;
2738	}
2739	} else {
2740	int row2 = mapRow[iRow];
2741	assert (iRow = mapRow[row2]);
2742	if (rowLower[iRow]!=rowLower[row2]\|\|
2743	rowLower[row2]!=rowLower[iRow])
2744	good=false;
2745	CoinBigIndex offset2 = rowStart[row2]-start;
2746	for (CoinBigIndex i=start;i<start+length;i++) {
2747	int jColumn = column[i];
2748	double value = element[i];
2749	int jColumn2 = column[i+offset2];
2750	double value2 = element[i+offset2];
2751	if (value!=value2\|\|mapColumn[jColumn]!=jColumn2\|\|
2752	mapColumn[jColumn2]!=jColumn)
2753	good=false;
2754	}
2755	}
2756	}
2757	if (good) {
2758	// check rim
2759	for (iColumn=0;iColumn<numberColumns;iColumn++) {
2760	if (mapColumn[iColumn]!=iColumn) {
2761	int iColumn2 = mapColumn[iColumn];
2762	if(objective[iColumn]!=objective[iColumn2])
2763	good=false;
2764	if (columnLower[iColumn]!=columnLower[iColumn2])
2765	good=false;
2766	if (columnUpper[iColumn]!=columnUpper[iColumn2])
2767	good=false;
2768	if (model->isInteger(iColumn)!=model->isInteger(iColumn2))
2769	good=false;
2770	}
2771	}
2772	}
2773	if (good) {
2774	// temp
2775	if (nMapRow<0) {
2776	//const double * solution = model->primalColumnSolution();
2777	// find mapped
2778	int nMapColumn=0;
2779	for (int i=0;i<numberColumns;i++) {
2780	if (mapColumn[i]>i)
2781	nMapColumn++;
2782	}
2783	nMapRow=0;
2784	int kRow=-1;
2785	for (int i=0;i<numberRows;i++) {
2786	if (mapRow[i]>i) {
2787	nMapRow++;
2788	kRow=i;
2789	}
2790	}
2791	printf("%d columns, %d rows\n",nMapColumn,nMapRow);
2792	if (nMapRow==1) {
2793	CoinBigIndex start = rowStart[kRow];
2794	int length=rowLength[kRow];
2795	printf("%g <= ",rowLower[kRow]);
2796	for (CoinBigIndex i=start;i<start+length;i++) {
2797	int jColumn = column[i];
2798	if (mapColumn[jColumn]!=jColumn)
2799	printf("* ");
2800	printf("%d,%g ",jColumn,element[i]);
2801	}
2802	printf("<= %g\n",rowUpper[kRow]);
2803	}
2804	}
2805	// temp
2806	int row1=sortRow[lastLook];
2807	int row2=sortRow[jLook];
2808	lastLook=jLook;
2809	CoinBigIndex start1 = rowStart[row1];
2810	CoinBigIndex offset2 = rowStart[row2]-start1;
2811	int length=rowLength[row1];
2812	assert( length==rowLength[row2]);
2813	CoinBigIndex put=startAdd[nAddRows];
2814	double multiplier=length<11 ? 2.0 : 1.125;
2815	double value=1.0;
2816	for (CoinBigIndex i=start1;i<start1+length;i++) {
2817	int jColumn1 = column[i];
2818	int jColumn2 = column[i+offset2];
2819	columnAdd[put]=jColumn1;
2820	elementAdd[put++]=value;
2821	columnAdd[put]=jColumn2;
2822	elementAdd[put++]=-value;
2823	value *= multiplier;
2824	}
2825	nAddRows++;
2826	startAdd[nAddRows]=put;
2827	} else {
2828	printf("ouch - did not check out as good\n");
2829	}
2830	}
2831	}
2832	// skip rest
2833	iLook += numberPossible-1;
2834	}
2835	}
2836	}
2837	if (nAddRows) {
2838	double * lower = new double [nAddRows];
2839	double * upper = new double[nAddRows];
2840	int i;
2841	//const double * solution = model->primalColumnSolution();
2842	for (i=0;i<nAddRows;i++) {
2843	lower[i]=0.0;
2844	upper[i]=COIN_DBL_MAX;
2845	}
2846	printf("Adding %d rows with %d elements\n",nAddRows,
2847	startAdd[nAddRows]);
2848	//ClpSimplex newModel(*model);
2849	//newModel.addRows(nAddRows,lower,upper,startAdd,columnAdd,elementAdd);
2850	//newModel.writeMps("modified.mps");
2851	delete [] lower;
2852	delete [] upper;
2853	}
2854	delete [] startAdd;
2855	delete [] columnAdd;
2856	delete [] elementAdd;
2857	delete [] order;
2858	delete [] other;
2859	delete [] randomColumn;
2860	delete [] weight;
2861	delete [] randomRow;
2862	delete [] sortRow;
2863	delete [] backRow;
2864	delete [] possibleRow;
2865	delete [] sortColumn;
2866	delete [] backColumn;
2867	delete [] backColumn2;
2868	delete [] possibleColumn;
2869	delete [] mapRow;
2870	delete [] mapColumn;
2871	delete [] stackRow;
2872	delete [] stackColumn;
2873	}
2874	delete [] number;
2875	// Now do breakdown of ranges
2876	breakdown("Elements",numberElements,elementByColumn);
2877	breakdown("RowLower",numberRows,rowLower);
2878	breakdown("RowUpper",numberRows,rowUpper);
2879	breakdown("ColumnLower",numberColumns,columnLower);
2880	breakdown("ColumnUpper",numberColumns,columnUpper);
2881	breakdown("Objective",numberColumns,objective);
2882	}
2883	static bool maskMatches(const int * starts, char ** masks,
2884	std::string & check)
2885	{
2886	// back to char as I am old fashioned
2887	const char * checkC = check.c_str();
2888	int length = strlen(checkC);
2889	while (checkC[length-1]==' ')
2890	length--;
2891	for (int i=starts[length];i<starts[length+1];i++) {
2892	char * thisMask = masks[i];
2893	int k;
2894	for ( k=0;k<length;k++) {
2895	if (thisMask[k]!='?'&&thisMask[k]!=checkC[k])
2896	break;
2897	}
2898	if (k==length)
2899	return true;
2900	}
2901	return false;
2902	}
2903	static void clean(char * temp)
2904	{
2905	char * put = temp;
2906	while (*put>=' ')
2907	put++;
2908	*put='\0';
2909	}
2910	static void generateCode(const char * fileName,int type)
2911	{
2912	FILE * fp = fopen(fileName,"r");
2913	assert (fp);
2914	int numberLines=0;
2915	#define MAXLINES 500
2916	#define MAXONELINE 200
2917	char line[MAXLINES][MAXONELINE];
2918	while (fgets(line[numberLines],MAXONELINE,fp)) {
2919	assert (numberLines<MAXLINES);
2920	clean(line[numberLines]);
2921	numberLines++;
2922	}
2923	fclose(fp);
2924	// add in actual solve
2925	strcpy(line[numberLines],"5 clpModel->initialSolve(clpSolve);");
2926	numberLines++;
2927	fp = fopen(fileName,"w");
2928	assert (fp);
2929	char apo='"';
2930	char backslash = '\\';
2931
2932	fprintf(fp,"#include %cClpSimplex.hpp%c\n",apo,apo);
2933	fprintf(fp,"#include %cClpSolve.hpp%c\n",apo,apo);
2934
2935	fprintf(fp,"\nint main (int argc, const char *argv[])\n{\n");
2936	fprintf(fp," ClpSimplex model;\n");
2937	fprintf(fp," int status=1;\n");
2938	fprintf(fp," if (argc<2)\n");
2939	fprintf(fp," fprintf(stderr,%cPlease give file name%cn%c);\n",
2940	apo,backslash,apo);
2941	fprintf(fp," else\n");
2942	fprintf(fp," status=model.readMps(argv[1],true);\n");
2943	fprintf(fp," if (status) {\n");
2944	fprintf(fp," fprintf(stderr,%cBad readMps %%s%cn%c,argv[1]);\n",
2945	apo,backslash,apo);
2946	fprintf(fp," exit(1);\n");
2947	fprintf(fp," }\n\n");
2948	fprintf(fp," // Now do requested saves and modifications\n");
2949	fprintf(fp," ClpSimplex * clpModel = & model;\n");
2950	int wanted[9];
2951	memset(wanted,0,sizeof(wanted));
2952	wanted[0]=wanted[3]=wanted[5]=wanted[8]=1;
2953	if (type>0)
2954	wanted[1]=wanted[6]=1;
2955	if (type>1)
2956	wanted[2]=wanted[4]=wanted[7]=1;
2957	std::string header[9]=
2958	{ "","Save values","Redundant save of default values","Set changed values",
2959	"Redundant set default values","Solve","Restore values","Redundant restore values","Add to model"};
2960	for (int iType=0;iType<9;iType++) {
2961	if (!wanted[iType])
2962	continue;
2963	int n=0;
2964	int iLine;
2965	for (iLine=0;iLine<numberLines;iLine++) {
2966	if (line[iLine][0]=='0'+iType) {
2967	if (!n)
2968	fprintf(fp,"\n // %s\n\n",header[iType].c_str());
2969	n++;
2970	fprintf(fp,"%s\n",line[iLine]+1);
2971	}
2972	}
2973	}
2974	fprintf(fp,"\n // Now you would use solution etc etc\n\n");
2975	fprintf(fp," return 0;\n}\n");
2976	fclose(fp);
2977	printf("C++ file written to %s\n",fileName);
2978	}
2979	/*
2980	Version 1.00.00 October 13 2004.
2981	1.00.01 October 18. Added basis handling helped/prodded by Thorsten Koch.
2982	Also modifications to make faster with sbb (I hope I haven't broken anything).
2983	1.00.02 March 21 2005. Redid ClpNonLinearCost to save memory also redid
2984	createRim to try and improve cache characteristics.
2985	1.00.03 April 8 2005. Added Volume algorithm as crash and made code more
2986	robust on testing. Also added "either" and "tune" algorithm.
2987	1.01.01 April 12 2005. Decided to go to different numbering. Backups will
2988	be last 2 digits while middle 2 are for improvements. Still take a long
2989	time to get to 2.00.01
2990	1.01.02 May 4 2005. Will be putting in many changes - so saving stable version
2991	1.02.01 May 6 2005. Lots of changes to try and make faster and more stable in
2992	branch and cut.
2993	1.02.02 May 19 2005. Stuff for strong branching and some improvements to simplex
2994	1.03.01 May 24 2006. Lots done but I can't remember what!
2995	1.03.03 June 13 2006. For clean up after dual perturbation
2996	1.04.01 June 26 2007. Lots of changes but I got lazy
2997	1.05.00 June 27 2007. This is trunk so when gets to stable will be 1.5
2998	*/

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

Download in other formats: