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

Last change on this file since 1682 was 1682, checked in by forrest, 9 years ago

stupid me - I did it again - glpk fix

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