Changeset 2030 for trunk/Clp/src


Ignore:
Timestamp:
Apr 15, 2014 11:54:11 AM (5 years ago)
Author:
forrest
Message:

fix some ampl stuff, add ClpSolver? and a few fixes

Location:
trunk/Clp/src
Files:
1 added
18 edited

Legend:

Unmodified
Added
Removed
  • trunk/Clp/src/AbcSimplex.cpp

    r2025 r2030  
    40444044#endif
    40454045  numberIterations_=totalNumberIterations;
    4046   return (problemStatus_<0||problemStatus_>3) ? 1 : 0;
     4046  return (problemStatus_<0||problemStatus_==4) ? 1 : 0;
    40474047}
    40484048int AbcSimplex::dual ()
  • trunk/Clp/src/AbcSimplexPrimal.cpp

    r2015 r2030  
    11941194    lastObj3 += infeasibilityCost_ * 2.0 * lastInf3;
    11951195    if (lastObj < thisObj - 1.0e-5 * CoinMax(fabs(thisObj), fabs(lastObj)) - 1.0e-7
    1196         && firstFree_ < 0) {
     1196        && firstFree_ < 0 && thisInf >= lastInf) {
    11971197#if ABC_NORMAL_DEBUG>0
    11981198      if (handler_->logLevel() == 63)
     
    12101210      }
    12111211    } else if (lastObj3 < thisObj - 1.0e-5 * CoinMax(fabs(thisObj), fabs(lastObj3)) - 1.0e-7
    1212                && firstFree_ < 0) {
     1212               && firstFree_ < 0 && thisInf >= lastInf) {
    12131213#if ABC_NORMAL_DEBUG>0
    12141214      if (handler_->logLevel() == 63)
     
    30153015  double tolerance = 100.0 * primalTolerance_;
    30163016  int numberChanged=0;
     3017  // Set bit if fixed
     3018  for (int i=0;i<numberRows_;i++) {
     3019    if (rowLower_[i]!=rowUpper_[i])
     3020      internalStatus_[i] &= ~128;
     3021    else
     3022      internalStatus_[i] |= 128;
     3023  }
     3024  for (int i=0;i<numberColumns_;i++) {
     3025    if (columnLower_[i]!=columnUpper_[i])
     3026      internalStatus_[i+numberRows_] &= ~128;
     3027    else
     3028      internalStatus_[i+numberRows_] |= 128;
     3029  }
    30173030  //double multiplier = perturbation*maximumFraction;
    30183031  for (iSequence = 0; iSequence < numberRows_ + numberColumns_; iSequence++) {
    30193032    if (getInternalStatus(iSequence) == basic) {
     3033      if ((internalStatus_[i] &128)!=0)
     3034        continue;
    30203035      double lowerValue = abcLower_[iSequence];
    30213036      double upperValue = abcUpper_[iSequence];
     
    30793094    for (iSequence = 0; iSequence < maximumAbcNumberRows_ + numberColumns_; iSequence++) {
    30803095      if (getInternalStatus(iSequence) != basic) {
     3096        if ((internalStatus_[i] &128)!=0)
     3097          continue;
    30813098        double lowerValue = abcLower_[iSequence];
    30823099        double upperValue = abcUpper_[iSequence];
     
    31263143  // Clean up
    31273144  for (i = 0; i < numberColumns_ + numberRows_; i++) {
     3145    internalStatus_[i] &= ~128;
    31283146    switch(getInternalStatus(i)) {
    31293147     
  • trunk/Clp/src/CbcOrClpParam.cpp

    r2016 r2030  
    36913691     parameters[numberParameters++] =
    36923692          CbcOrClpParam("tune!PreProcess", "Dubious tuning parameters",
    3693                         0, 20000000, CLP_PARAM_INT_PROCESSTUNE, 1);
    3694      parameters[numberParameters-1].setLonghelp
    3695      (
    3696           "For making equality cliques this is minimumsize.  Also for adding \
    3697 integer slacks.  May be used for more later \
    3698 If <10000 that is what it does.  If <1000000 - numberPasses is (value/10000)-1 and tune is tune %10000. \
    3699 If >= 1000000! - numberPasses is (value/1000000)-1 and tune is tune %1000000.  In this case if tune is now still >=10000 \
    3700 numberPassesPerInnerLoop is changed from 10 to (tune-10000)-1 and tune becomes tune % 10000!!!!! - happy? - \
    3701 so to keep normal limit on cliques of 5, do 3 major passes (include presolves) but only doing one tightening pass per major pass - \
    3702 you would use 3010005 (I think)"
     3693                        0, 2000000000, CLP_PARAM_INT_PROCESSTUNE, 1);
     3694     parameters[numberParameters-1].setLonghelp
     3695     (
     3696          "Format aabbcccc - \n If aa then this is number of major passes (i.e. with presolve) \n \
     3697If bb and bb>0 then this is number of minor passes (if unset or 0 then 10) \n \
     3698cccc is bit set \n 0 - 1 Heavy probing \n 1 - 2 Make variables integer if possible (if obj value)\n \
     36992 - 4 As above but even if zero objective value\n \
     37007 - 128 Try and create cliques\n 8 - 256 If all +1 try hard for dominated rows\n \
     370110 - 1024 Use a larger feasibility tolerance in presolve\n \
     370211 - 2048 Try probing before creating cliques"
    37033703     );
    37043704     parameters[numberParameters++] =
  • trunk/Clp/src/ClpCholeskyBase.cpp

    r1974 r2030  
    196196          workInteger_ = ClpCopyOfArray(rhs.workInteger_, numberRows_);
    197197          clique_ = ClpCopyOfArray(rhs.clique_, numberRows_);
    198           delete rowCopy_;
    199198          rowCopy_ = rhs.rowCopy_->clone();
    200199          whichDense_ = NULL;
  • trunk/Clp/src/ClpMain.cpp

    r1975 r2030  
    66#include "CoinPragma.hpp"
    77
    8 #include <cassert>
    9 #include <cstdio>
    10 #include <cstdlib>
    11 #include <cmath>
    12 #include <cfloat>
    13 #include <string>
    14 #include <iostream>
    15 
    16 int boundary_sort = 1000;
    17 int boundary_sort2 = 1000;
    18 int boundary_sort3 = 10000;
    19 // for printing
    20 #ifndef CLP_OUTPUT_FORMAT
    21 #define CLP_OUTPUT_FORMAT %15.8g
    22 #endif
    23 #define CLP_QUOTE(s) CLP_STRING(s)
    24 #define CLP_STRING(s) #s
    25 #include "CoinPragma.hpp"
    26 #include "CoinHelperFunctions.hpp"
    27 #include "CoinSort.hpp"
    28 // History since 1.0 at end
    29 #include "ClpConfig.h"
    30 #include "CoinMpsIO.hpp"
    31 #include "CoinFileIO.hpp"
    32 #ifdef COIN_HAS_GLPK
    33 #include "glpk.h"
    34 extern glp_tran* cbc_glp_tran;
    35 extern glp_prob* cbc_glp_prob;
    36 #else
    37 #define GLP_UNDEF 1
    38 #define GLP_FEAS 2
    39 #define GLP_INFEAS 3
    40 #define GLP_NOFEAS 4
    41 #define GLP_OPT 5
    42 #endif
    43 
    448#include "AbcCommon.hpp"
    45 #include "ClpFactorization.hpp"
    46 #include "CoinTime.hpp"
    47 #include "CoinWarmStartBasis.hpp"
    489#include "ClpSimplex.hpp"
    49 #include "ClpSimplexOther.hpp"
    50 #include "ClpSolve.hpp"
    51 #include "ClpMessage.hpp"
    52 #include "ClpPackedMatrix.hpp"
    53 #include "ClpPlusMinusOneMatrix.hpp"
    54 #include "ClpNetworkMatrix.hpp"
    55 #include "ClpDualRowSteepest.hpp"
    56 #include "ClpDualRowDantzig.hpp"
    57 #include "ClpLinearObjective.hpp"
    58 #include "ClpPrimalColumnSteepest.hpp"
    59 #include "ClpPrimalColumnDantzig.hpp"
    60 #include "ClpPresolve.hpp"
    61 #include "CbcOrClpParam.hpp"
    62 #include "CoinSignal.hpp"
    6310#ifdef ABC_INHERIT
    6411#include "AbcSimplex.hpp"
    65 #include "AbcSimplexFactorization.hpp"
    66 #include "AbcDualRowSteepest.hpp"
    67 #include "AbcDualRowDantzig.hpp"
    6812#endif
    69 #ifdef COIN_HAS_ASL
    70 #include "Clp_ampl.h"
    71 #endif
    72 #ifdef DMALLOC
    73 #include "dmalloc.h"
    74 #endif
    75 #if defined(COIN_HAS_WSMP) || defined(COIN_HAS_AMD) || defined(COIN_HAS_CHOLMOD) || defined(TAUCS_BARRIER) || defined(COIN_HAS_MUMPS)
    76 #define FOREIGN_BARRIER
    77 #endif
    78 
    79 static double totalTime = 0.0;
    80 static bool maskMatches(const int * starts, char ** masks,
    81                         std::string & check);
    8213#ifndef ABC_INHERIT
    83 static ClpSimplex * currentModel = NULL;
     14void ClpMain0(ClpSimplex * models);
     15int ClpMain1(int argc, const char *argv[],ClpSimplex * model);
    8416#else
    85 static AbcSimplex * currentModel = NULL;
    86 #endif
    87 
    88 extern "C" {
    89      static void
    90 #if defined(_MSC_VER)
    91      __cdecl
    92 #endif // _MSC_VER
    93      signal_handler(int /*whichSignal*/)
    94      {
    95           if (currentModel != NULL)
    96                currentModel->setMaximumIterations(0); // stop at next iterations
    97           return;
    98      }
    99 }
    100 
    101 //#############################################################################
    102 
    103 #ifdef NDEBUG
    104 #undef NDEBUG
    105 #endif
    106 
    107 #ifndef ABC_INHERIT
    108 int mainTest (int argc, const char *argv[], int algorithm,
    109               ClpSimplex empty, ClpSolve solveOptions, int switchOff, bool doVector);
    110 #else
    111 int mainTest (int argc, const char *argv[], int algorithm,
    112               AbcSimplex empty, ClpSolve solveOptions, int switchOff, bool doVector);
    113 #endif
    114 static void statistics(ClpSimplex * originalModel, ClpSimplex * model);
    115 static void generateCode(const char * fileName, int type);
    116 // Returns next valid field
    117 int CbcOrClpRead_mode = 1;
    118 FILE * CbcOrClpReadCommand = stdin;
    119 extern int CbcOrClpEnvironmentIndex;
    120 #ifdef CLP_USER_DRIVEN1
    121 /* Returns true if variable sequenceOut can leave basis when
    122    model->sequenceIn() enters.
    123    This function may be entered several times for each sequenceOut. 
    124    The first time realAlpha will be positive if going to lower bound
    125    and negative if going to upper bound (scaled bounds in lower,upper) - then will be zero.
    126    currentValue is distance to bound.
    127    currentTheta is current theta.
    128    alpha is fabs(pivot element).
    129    Variable will change theta if currentValue - currentTheta*alpha < 0.0
    130 */
    131 bool userChoiceValid1(const ClpSimplex * model,
    132                       int sequenceOut,
    133                       double currentValue,
    134                       double currentTheta,
    135                       double alpha,
    136                       double realAlpha)
    137 {
    138   return true;
    139 }
    140 /* This returns true if chosen in/out pair valid.
    141    The main thing to check would be variable flipping bounds may be
    142    OK.  This would be signaled by reasonable theta_ and valueOut_.
    143    If you return false sequenceIn_ will be flagged as ineligible.
    144 */
    145 bool userChoiceValid2(const ClpSimplex * model)
    146 {
    147   return true;
    148 }
    149 /* If a good pivot then you may wish to unflag some variables.
    150  */
    151 void userChoiceWasGood(ClpSimplex * model)
    152 {
    153 }
     17void ClpMain0(AbcSimplex * models);
     18int ClpMain1(int argc, const char *argv[],AbcSimplex * model);
    15419#endif
    15520//#define CILK_TEST
     
    16631  cilkTest();
    16732#endif
    168      // next {} is just to make sure all memory should be freed - for debug
    169      {
    170           double time1 = CoinCpuTime(), time2;
    171           // Set up all non-standard stuff
    172           //int numberModels=1;
    17333#ifndef ABC_INHERIT
    174           ClpSimplex * models = new ClpSimplex[1];
     34  ClpSimplex * models = new ClpSimplex[1];
    17535#else
    176           AbcSimplex * models = new AbcSimplex[1];
     36  AbcSimplex * models = new AbcSimplex[1];
    17737#endif
    178 
    179           // default action on import
    180           int allowImportErrors = 0;
    181           int keepImportNames = 1;
    182           int doIdiot = -1;
    183           int outputFormat = 2;
    184           int slpValue = -1;
    185           int cppValue = -1;
    186           int printOptions = 0;
    187           int printMode = 0;
    188           int presolveOptions = 0;
    189           int doCrash = 0;
    190           int doVector = 0;
    191           int doSprint = -1;
    192           // set reasonable defaults
    193 #ifdef ABC_INHERIT
    194 #define DEFAULT_PRESOLVE_PASSES 20
    195 #else
    196 #define DEFAULT_PRESOLVE_PASSES 5
    197 #endif
    198           int preSolve = DEFAULT_PRESOLVE_PASSES;
    199           bool preSolveFile = false;
    200           models->setPerturbation(50);
    201           models->messageHandler()->setPrefix(false);
    202           const char dirsep =  CoinFindDirSeparator();
    203           std::string directory;
    204           std::string dirSample;
    205           std::string dirNetlib;
    206           std::string dirMiplib;
    207           if (dirsep == '/') {
    208                directory = "./";
    209                dirSample = "../../Data/Sample/";
    210                dirNetlib = "../../Data/Netlib/";
    211                dirMiplib = "../../Data/miplib3/";
    212           } else {
    213                directory = ".\\";
    214 #              ifdef COIN_MSVS
    215                // Visual Studio builds are deeper
    216                dirSample = "..\\..\\..\\..\\Data\\Sample\\";
    217                dirNetlib = "..\\..\\..\\..\\Data\\Netlib\\";
    218                dirMiplib = "..\\..\\..\\..\\Data\\miplib3\\";
    219 #              else
    220                dirSample = "..\\..\\Data\\Sample\\";
    221                dirNetlib = "..\\..\\Data\\Netlib\\";
    222                dirMiplib = "..\\..\\Data\\miplib3\\";
    223 #              endif
    224           }
    225           std::string defaultDirectory = directory;
    226           std::string importFile = "";
    227           std::string exportFile = "default.mps";
    228           std::string importBasisFile = "";
    229           int basisHasValues = 0;
    230           int substitution = 3;
    231           int dualize = 3;  // dualize if looks promising
    232           std::string exportBasisFile = "default.bas";
    233           std::string saveFile = "default.prob";
    234           std::string restoreFile = "default.prob";
    235           std::string solutionFile = "stdout";
    236           std::string solutionSaveFile = "solution.file";
    237           std::string printMask = "";
    238           CbcOrClpParam parameters[CBCMAXPARAMETERS];
    239           int numberParameters ;
    240           establishParams(numberParameters, parameters) ;
    241           parameters[whichParam(CLP_PARAM_ACTION_BASISIN, numberParameters, parameters)].setStringValue(importBasisFile);
    242           parameters[whichParam(CLP_PARAM_ACTION_BASISOUT, numberParameters, parameters)].setStringValue(exportBasisFile);
    243           parameters[whichParam(CLP_PARAM_ACTION_PRINTMASK, numberParameters, parameters)].setStringValue(printMask);
    244           parameters[whichParam(CLP_PARAM_ACTION_DIRECTORY, numberParameters, parameters)].setStringValue(directory);
    245           parameters[whichParam(CLP_PARAM_ACTION_DIRSAMPLE, numberParameters, parameters)].setStringValue(dirSample);
    246           parameters[whichParam(CLP_PARAM_ACTION_DIRNETLIB, numberParameters, parameters)].setStringValue(dirNetlib);
    247           parameters[whichParam(CBC_PARAM_ACTION_DIRMIPLIB, numberParameters, parameters)].setStringValue(dirMiplib);
    248           parameters[whichParam(CLP_PARAM_DBL_DUALBOUND, numberParameters, parameters)].setDoubleValue(models->dualBound());
    249           parameters[whichParam(CLP_PARAM_DBL_DUALTOLERANCE, numberParameters, parameters)].setDoubleValue(models->dualTolerance());
    250           parameters[whichParam(CLP_PARAM_ACTION_EXPORT, numberParameters, parameters)].setStringValue(exportFile);
    251           parameters[whichParam(CLP_PARAM_INT_IDIOT, numberParameters, parameters)].setIntValue(doIdiot);
    252           parameters[whichParam(CLP_PARAM_ACTION_IMPORT, numberParameters, parameters)].setStringValue(importFile);
    253           parameters[whichParam(CLP_PARAM_INT_SOLVERLOGLEVEL, numberParameters, parameters)].setIntValue(models->logLevel());
    254           parameters[whichParam(CLP_PARAM_INT_MAXFACTOR, numberParameters, parameters)].setIntValue(models->factorizationFrequency());
    255           parameters[whichParam(CLP_PARAM_INT_MAXITERATION, numberParameters, parameters)].setIntValue(models->maximumIterations());
    256           parameters[whichParam(CLP_PARAM_INT_OUTPUTFORMAT, numberParameters, parameters)].setIntValue(outputFormat);
    257           parameters[whichParam(CLP_PARAM_INT_PRESOLVEPASS, numberParameters, parameters)].setIntValue(preSolve);
    258           parameters[whichParam(CLP_PARAM_INT_PERTVALUE, numberParameters, parameters)].setIntValue(models->perturbation());
    259           parameters[whichParam(CLP_PARAM_DBL_PRIMALTOLERANCE, numberParameters, parameters)].setDoubleValue(models->primalTolerance());
    260           parameters[whichParam(CLP_PARAM_DBL_PRIMALWEIGHT, numberParameters, parameters)].setDoubleValue(models->infeasibilityCost());
    261           parameters[whichParam(CLP_PARAM_ACTION_RESTORE, numberParameters, parameters)].setStringValue(restoreFile);
    262           parameters[whichParam(CLP_PARAM_ACTION_SAVE, numberParameters, parameters)].setStringValue(saveFile);
    263           parameters[whichParam(CLP_PARAM_DBL_TIMELIMIT, numberParameters, parameters)].setDoubleValue(models->maximumSeconds());
    264           parameters[whichParam(CLP_PARAM_ACTION_SOLUTION, numberParameters, parameters)].setStringValue(solutionFile);
    265           parameters[whichParam(CLP_PARAM_ACTION_SAVESOL, numberParameters, parameters)].setStringValue(solutionSaveFile);
    266           parameters[whichParam(CLP_PARAM_INT_SPRINT, numberParameters, parameters)].setIntValue(doSprint);
    267           parameters[whichParam(CLP_PARAM_INT_SUBSTITUTION, numberParameters, parameters)].setIntValue(substitution);
    268           parameters[whichParam(CLP_PARAM_INT_DUALIZE, numberParameters, parameters)].setIntValue(dualize);
    269           parameters[whichParam(CLP_PARAM_DBL_PRESOLVETOLERANCE, numberParameters, parameters)].setDoubleValue(1.0e-8);
    270           int verbose = 0;
    271 
    272           // total number of commands read
    273           int numberGoodCommands = 0;
    274           bool * goodModels = new bool[1];
    275           goodModels[0] = false;
    276 #ifdef COIN_HAS_ASL
    277           ampl_info info;
    278           int usingAmpl=0;
    279           CoinMessageHandler * generalMessageHandler = models->messageHandler();
    280           generalMessageHandler->setPrefix(false);
    281           CoinMessages generalMessages = models->messages();
    282           char generalPrint[10000];
    283           {
    284             bool noPrinting_=false;
    285             memset(&info, 0, sizeof(info));
    286             if (argc > 2 && !strcmp(argv[2], "-AMPL")) {
    287               usingAmpl = 1;
    288               // see if log in list
    289               noPrinting_ = true;
    290               for (int i = 1; i < argc; i++) {
    291                     if (!strncmp(argv[i], "log", 3)) {
    292                         const char * equals = strchr(argv[i], '=');
    293                         if (equals && atoi(equals + 1) > 0) {
    294                             noPrinting_ = false;
    295                             info.logLevel = atoi(equals + 1);
    296                             int log = whichParam(CLP_PARAM_INT_LOGLEVEL, numberParameters, parameters);
    297                             parameters[log].setIntValue(info.logLevel);
    298                             // mark so won't be overWritten
    299                             info.numberRows = -1234567;
    300                             break;
    301                         }
    302                     }
    303                 }
    304 
    305                 union {
    306                     void * voidModel;
    307                     CoinModel * model;
    308                 } coinModelStart;
    309                 coinModelStart.model = NULL;
    310                 int returnCode = readAmpl(&info, argc, const_cast<char **>(argv), & coinModelStart.voidModel);
    311                 if (returnCode)
    312                     return returnCode;
    313                 if (info.numberBinary+info.numberIntegers+info.numberSos
    314                     &&!info.starts) {
    315                   printf("Unable to handle integer problems\n");
    316                   return 1;
    317                 }
    318                 CbcOrClpRead_mode = 2; // so will start with parameters
    319                 // see if log in list (including environment)
    320                 for (int i = 1; i < info.numberArguments; i++) {
    321                     if (!strcmp(info.arguments[i], "log")) {
    322                         if (i < info.numberArguments - 1 && atoi(info.arguments[i+1]) > 0)
    323                             noPrinting_ = false;
    324                         break;
    325                     }
    326                 }
    327                 if (noPrinting_) {
    328                     models->messageHandler()->setLogLevel(0);
    329                     setCbcOrClpPrinting(false);
    330                 }
    331                 if (!noPrinting_)
    332                     printf("%d rows, %d columns and %d elements\n",
    333                            info.numberRows, info.numberColumns, info.numberElements);
    334                 if (!coinModelStart.model) {
    335                   // linear
    336                     models->loadProblem(info.numberColumns, info.numberRows, info.starts,
    337                                         info.rows, info.elements,
    338                                         info.columnLower, info.columnUpper, info.objective,
    339                                         info.rowLower, info.rowUpper);
    340                 } else {
    341                   // QP
    342                   models->loadProblem(*(coinModelStart.model));
    343                 }
    344                 // If we had a solution use it
    345                 if (info.primalSolution) {
    346                     models->setColSolution(info.primalSolution);
    347                 }
    348                 // status
    349                 if (info.rowStatus) {
    350                     unsigned char * statusArray = models->statusArray();
    351                     int i;
    352                     for (i = 0; i < info.numberColumns; i++)
    353                         statusArray[i] = static_cast<unsigned char>(info.columnStatus[i]);
    354                     statusArray += info.numberColumns;
    355                     for (i = 0; i < info.numberRows; i++)
    356                         statusArray[i] = static_cast<unsigned char>(info.rowStatus[i]);
    357                 }
    358                 freeArrays1(&info);
    359                 // modify objective if necessary
    360                 models->setOptimizationDirection(info.direction);
    361                 models->setObjectiveOffset(info.offset);
    362                 if (info.offset) {
    363                     sprintf(generalPrint, "Ampl objective offset is %g",
    364                             info.offset);
    365                     generalMessageHandler->message(CLP_GENERAL, generalMessages)
    366                     << generalPrint
    367                     << CoinMessageEol;
    368                 }
    369                 goodModels[0] = true;
    370                 // change argc etc
    371                 argc = info.numberArguments;
    372                 argv = const_cast<const char **>(info.arguments);
    373             }
    374         }
    375 #endif
    376 
    377           // Hidden stuff for barrier
    378           int choleskyType = 0;
    379           int gamma = 0;
    380           parameters[whichParam(CLP_PARAM_STR_BARRIERSCALE, numberParameters, parameters)].setCurrentOption(2);
    381           int scaleBarrier = 2;
    382           int doKKT = 0;
    383           int crossover = 2; // do crossover unless quadratic
    384 
    385           int iModel = 0;
    386           //models[0].scaling(1);
    387           //models[0].setDualBound(1.0e6);
    388           //models[0].setDualTolerance(1.0e-7);
    389           //ClpDualRowSteepest steep;
    390           //models[0].setDualRowPivotAlgorithm(steep);
    391 #ifdef ABC_INHERIT
    392           models[0].setDualTolerance(1.0e-6);
    393           models[0].setPrimalTolerance(1.0e-6);
    394 #endif
    395           //ClpPrimalColumnSteepest steepP;
    396           //models[0].setPrimalColumnPivotAlgorithm(steepP);
    397           std::string field;
    398           std::cout << "Coin LP version " << CLP_VERSION
    399                     << ", build " << __DATE__ << std::endl;
    400           // Print command line
    401           if (argc > 1) {
    402                printf("command line - ");
    403                for (int i = 0; i < argc; i++)
    404                     printf("%s ", argv[i]);
    405                printf("\n");
    406           }
    407 
    408           while (1) {
    409                // next command
    410                field = CoinReadGetCommand(argc, argv);
    411 
    412                // exit if null or similar
    413                if (!field.length()) {
    414                     if (numberGoodCommands == 1 && goodModels[0]) {
    415                          // we just had file name - do dual or primal
    416                          field = "either";
    417                     } else if (!numberGoodCommands) {
    418                          // let's give the sucker a hint
    419                          std::cout
    420                                    << "Clp takes input from arguments ( - switches to stdin)"
    421                                    << std::endl
    422                                    << "Enter ? for list of commands or help" << std::endl;
    423                          field = "-";
    424                     } else {
    425                          break;
    426                     }
    427                }
    428 
    429                // see if ? at end
    430                int numberQuery = 0;
    431                if (field != "?" && field != "???") {
    432                     size_t length = field.length();
    433                     size_t i;
    434                     for (i = length - 1; i > 0; i--) {
    435                          if (field[i] == '?')
    436                               numberQuery++;
    437                          else
    438                               break;
    439                     }
    440                     field = field.substr(0, length - numberQuery);
    441                }
    442                // find out if valid command
    443                int iParam;
    444                int numberMatches = 0;
    445                int firstMatch = -1;
    446                for ( iParam = 0; iParam < numberParameters; iParam++ ) {
    447                     int match = parameters[iParam].matches(field);
    448                     if (match == 1) {
    449                          numberMatches = 1;
    450                          firstMatch = iParam;
    451                          break;
    452                     } else {
    453                          if (match && firstMatch < 0)
    454                               firstMatch = iParam;
    455                          numberMatches += match >> 1;
    456                     }
    457                }
    458                ClpSimplex * thisModel=models+iModel;
    459                if (iParam < numberParameters && !numberQuery) {
    460                     // found
    461                     CbcOrClpParam found = parameters[iParam];
    462                     CbcOrClpParameterType type = found.type();
    463                     int valid;
    464                     numberGoodCommands++;
    465                     if (type == CBC_PARAM_GENERALQUERY) {
    466                          std::cout << "In argument list keywords have leading - "
    467                                    ", -stdin or just - switches to stdin" << std::endl;
    468                          std::cout << "One command per line (and no -)" << std::endl;
    469                          std::cout << "abcd? gives list of possibilities, if only one + explanation" << std::endl;
    470                          std::cout << "abcd?? adds explanation, if only one fuller help" << std::endl;
    471                          std::cout << "abcd without value (where expected) gives current value" << std::endl;
    472                          std::cout << "abcd value sets value" << std::endl;
    473                          std::cout << "Commands are:" << std::endl;
    474                          int maxAcross = 10;
    475                          bool evenHidden = false;
    476                          int printLevel =
    477                               parameters[whichParam(CLP_PARAM_STR_ALLCOMMANDS,
    478                                                     numberParameters, parameters)].currentOptionAsInteger();
    479                          int convertP[] = {2, 1, 0};
    480                          printLevel = convertP[printLevel];
    481                          if ((verbose & 8) != 0) {
    482                               // even hidden
    483                               evenHidden = true;
    484                               verbose &= ~8;
    485                          }
    486 #ifdef COIN_HAS_ASL
    487                          if (verbose < 4 && usingAmpl)
    488                            verbose += 4;
    489 #endif
    490                          if (verbose)
    491                               maxAcross = 1;
    492                          int limits[] = {1, 101, 201, 301, 401};
    493                          std::vector<std::string> types;
    494                          types.push_back("Double parameters:");
    495                          types.push_back("Int parameters:");
    496                          types.push_back("Keyword parameters:");
    497                          types.push_back("Actions or string parameters:");
    498                          int iType;
    499                          for (iType = 0; iType < 4; iType++) {
    500                               int across = 0;
    501                               int lengthLine = 0;
    502                               if ((verbose % 4) != 0)
    503                                    std::cout << std::endl;
    504                               std::cout << types[iType] << std::endl;
    505                               if ((verbose & 2) != 0)
    506                                    std::cout << std::endl;
    507                               for ( iParam = 0; iParam < numberParameters; iParam++ ) {
    508                                    int type = parameters[iParam].type();
    509                                    //printf("%d type %d limits %d %d display %d\n",iParam,
    510                                    //   type,limits[iType],limits[iType+1],parameters[iParam].displayThis());
    511                                    if ((parameters[iParam].displayThis() >= printLevel || evenHidden) &&
    512                                              type >= limits[iType]
    513                                              && type < limits[iType+1]) {
    514                                         if (!across) {
    515                                              if ((verbose & 2) != 0)
    516                                                   std::cout << "Command ";
    517                                         }
    518                                         int length = parameters[iParam].lengthMatchName() + 1;
    519                                         if (lengthLine + length > 80) {
    520                                              std::cout << std::endl;
    521                                              across = 0;
    522                                              lengthLine = 0;
    523                                         }
    524                                         std::cout << " " << parameters[iParam].matchName();
    525                                         lengthLine += length ;
    526                                         across++;
    527                                         if (across == maxAcross) {
    528                                              across = 0;
    529                                              if (verbose) {
    530                                                   // put out description as well
    531                                                   if ((verbose & 1) != 0)
    532                                                        std::cout << parameters[iParam].shortHelp();
    533                                                   std::cout << std::endl;
    534                                                   if ((verbose & 2) != 0) {
    535                                                        std::cout << "---- description" << std::endl;
    536                                                        parameters[iParam].printLongHelp();
    537                                                        std::cout << "----" << std::endl << std::endl;
    538                                                   }
    539                                              } else {
    540                                                   std::cout << std::endl;
    541                                              }
    542                                         }
    543                                    }
    544                               }
    545                               if (across)
    546                                    std::cout << std::endl;
    547                          }
    548                     } else if (type == CBC_PARAM_FULLGENERALQUERY) {
    549                          std::cout << "Full list of commands is:" << std::endl;
    550                          int maxAcross = 5;
    551                          int limits[] = {1, 101, 201, 301, 401};
    552                          std::vector<std::string> types;
    553                          types.push_back("Double parameters:");
    554                          types.push_back("Int parameters:");
    555                          types.push_back("Keyword parameters and others:");
    556                          types.push_back("Actions:");
    557                          int iType;
    558                          for (iType = 0; iType < 4; iType++) {
    559                               int across = 0;
    560                               std::cout << types[iType] << std::endl;
    561                               for ( iParam = 0; iParam < numberParameters; iParam++ ) {
    562                                    int type = parameters[iParam].type();
    563                                    if (type >= limits[iType]
    564                                              && type < limits[iType+1]) {
    565                                         if (!across)
    566                                              std::cout << "  ";
    567                                         std::cout << parameters[iParam].matchName() << "  ";
    568                                         across++;
    569                                         if (across == maxAcross) {
    570                                              std::cout << std::endl;
    571                                              across = 0;
    572                                         }
    573                                    }
    574                               }
    575                               if (across)
    576                                    std::cout << std::endl;
    577                          }
    578                     } else if (type < 101) {
    579                          // get next field as double
    580                          double value = CoinReadGetDoubleField(argc, argv, &valid);
    581                          if (!valid) {
    582                               parameters[iParam].setDoubleParameter(thisModel, value);
    583                          } else if (valid == 1) {
    584                               std::cout << " is illegal for double parameter " << parameters[iParam].name() << " value remains " <<
    585                                         parameters[iParam].doubleValue() << std::endl;
    586                          } else {
    587                               std::cout << parameters[iParam].name() << " has value " <<
    588                                         parameters[iParam].doubleValue() << std::endl;
    589                          }
    590                     } else if (type < 201) {
    591                          // get next field as int
    592                          int value = CoinReadGetIntField(argc, argv, &valid);
    593                          if (!valid) {
    594                               if (parameters[iParam].type() == CLP_PARAM_INT_PRESOLVEPASS)
    595                                    preSolve = value;
    596                               else if (parameters[iParam].type() == CLP_PARAM_INT_IDIOT)
    597                                    doIdiot = value;
    598                               else if (parameters[iParam].type() == CLP_PARAM_INT_SPRINT)
    599                                    doSprint = value;
    600                               else if (parameters[iParam].type() == CLP_PARAM_INT_OUTPUTFORMAT)
    601                                    outputFormat = value;
    602                               else if (parameters[iParam].type() == CLP_PARAM_INT_SLPVALUE)
    603                                    slpValue = value;
    604                               else if (parameters[iParam].type() == CLP_PARAM_INT_CPP)
    605                                    cppValue = value;
    606                               else if (parameters[iParam].type() == CLP_PARAM_INT_PRESOLVEOPTIONS)
    607                                    presolveOptions = value;
    608                               else if (parameters[iParam].type() == CLP_PARAM_INT_PRINTOPTIONS)
    609                                    printOptions = value;
    610                               else if (parameters[iParam].type() == CLP_PARAM_INT_SUBSTITUTION)
    611                                    substitution = value;
    612                               else if (parameters[iParam].type() == CLP_PARAM_INT_DUALIZE)
    613                                    dualize = value;
    614                               else if (parameters[iParam].type() == CLP_PARAM_INT_VERBOSE)
    615                                    verbose = value;
    616                               parameters[iParam].setIntParameter(thisModel, value);
    617                          } else if (valid == 1) {
    618                               std::cout << " is illegal for integer parameter " << parameters[iParam].name() << " value remains " <<
    619                                         parameters[iParam].intValue() << std::endl;
    620                          } else {
    621                               std::cout << parameters[iParam].name() << " has value " <<
    622                                         parameters[iParam].intValue() << std::endl;
    623                          }
    624                     } else if (type < 301) {
    625                          // one of several strings
    626                          std::string value = CoinReadGetString(argc, argv);
    627                          int action = parameters[iParam].parameterOption(value);
    628                          if (action < 0) {
    629                               if (value != "EOL") {
    630                                    // no match
    631                                    parameters[iParam].printOptions();
    632                               } else {
    633                                    // print current value
    634                                    std::cout << parameters[iParam].name() << " has value " <<
    635                                              parameters[iParam].currentOption() << std::endl;
    636                               }
    637                          } else {
    638                               parameters[iParam].setCurrentOption(action);
    639                               // for now hard wired
    640                               switch (type) {
    641                               case CLP_PARAM_STR_DIRECTION:
    642                                 if (action == 0) {
    643                                         models[iModel].setOptimizationDirection(1);
    644 #ifdef ABC_INHERIT
    645                                         thisModel->setOptimizationDirection(1);
    646 #endif
    647                                 }  else if (action == 1) {
    648                                         models[iModel].setOptimizationDirection(-1);
    649 #ifdef ABC_INHERIT
    650                                         thisModel->setOptimizationDirection(-1);
    651 #endif
    652                                 }  else {
    653                                         models[iModel].setOptimizationDirection(0);
    654 #ifdef ABC_INHERIT
    655                                         thisModel->setOptimizationDirection(0);
    656 #endif
    657                                 }
    658                                    break;
    659                               case CLP_PARAM_STR_DUALPIVOT:
    660                                    if (action == 0) {
    661                                         ClpDualRowSteepest steep(3);
    662                                         thisModel->setDualRowPivotAlgorithm(steep);
    663 #ifdef ABC_INHERIT
    664                                         AbcDualRowSteepest steep2(3);
    665                                         models[iModel].setDualRowPivotAlgorithm(steep2);
    666 #endif
    667                                    } else if (action == 1) {
    668                                         //ClpDualRowDantzig dantzig;
    669                                         ClpDualRowDantzig dantzig;
    670                                         thisModel->setDualRowPivotAlgorithm(dantzig);
    671 #ifdef ABC_INHERIT
    672                                         AbcDualRowDantzig dantzig2;
    673                                         models[iModel].setDualRowPivotAlgorithm(dantzig2);
    674 #endif
    675                                    } else if (action == 2) {
    676                                         // partial steep
    677                                         ClpDualRowSteepest steep(2);
    678                                         thisModel->setDualRowPivotAlgorithm(steep);
    679 #ifdef ABC_INHERIT
    680                                         AbcDualRowSteepest steep2(2);
    681                                         models[iModel].setDualRowPivotAlgorithm(steep2);
    682 #endif
    683                                    } else {
    684                                         ClpDualRowSteepest steep;
    685                                         thisModel->setDualRowPivotAlgorithm(steep);
    686 #ifdef ABC_INHERIT
    687                                         AbcDualRowSteepest steep2;
    688                                         models[iModel].setDualRowPivotAlgorithm(steep2);
    689 #endif
    690                                    }
    691                                    break;
    692                               case CLP_PARAM_STR_PRIMALPIVOT:
    693                                    if (action == 0) {
    694                                         ClpPrimalColumnSteepest steep(3);
    695                                         thisModel->setPrimalColumnPivotAlgorithm(steep);
    696                                    } else if (action == 1) {
    697                                         ClpPrimalColumnSteepest steep(0);
    698                                         thisModel->setPrimalColumnPivotAlgorithm(steep);
    699                                    } else if (action == 2) {
    700                                         ClpPrimalColumnDantzig dantzig;
    701                                         thisModel->setPrimalColumnPivotAlgorithm(dantzig);
    702                                    } else if (action == 3) {
    703                                         ClpPrimalColumnSteepest steep(4);
    704                                         thisModel->setPrimalColumnPivotAlgorithm(steep);
    705                                    } else if (action == 4) {
    706                                         ClpPrimalColumnSteepest steep(1);
    707                                         thisModel->setPrimalColumnPivotAlgorithm(steep);
    708                                    } else if (action == 5) {
    709                                         ClpPrimalColumnSteepest steep(2);
    710                                         thisModel->setPrimalColumnPivotAlgorithm(steep);
    711                                    } else if (action == 6) {
    712                                         ClpPrimalColumnSteepest steep(10);
    713                                         thisModel->setPrimalColumnPivotAlgorithm(steep);
    714                                    }
    715                                    break;
    716                               case CLP_PARAM_STR_SCALING:
    717                                    thisModel->scaling(action);
    718                                    break;
    719                               case CLP_PARAM_STR_AUTOSCALE:
    720                                    thisModel->setAutomaticScaling(action != 0);
    721                                    break;
    722                               case CLP_PARAM_STR_SPARSEFACTOR:
    723                                    thisModel->setSparseFactorization((1 - action) != 0);
    724                                    break;
    725                               case CLP_PARAM_STR_BIASLU:
    726                                    thisModel->factorization()->setBiasLU(action);
    727                                    break;
    728                               case CLP_PARAM_STR_PERTURBATION:
    729                                    if (action == 0)
    730                                         thisModel->setPerturbation(50);
    731                                    else
    732                                         thisModel->setPerturbation(100);
    733                                    break;
    734                               case CLP_PARAM_STR_ERRORSALLOWED:
    735                                    allowImportErrors = action;
    736                                    break;
    737                               case CLP_PARAM_STR_ABCWANTED:
    738                                    models[iModel].setAbcState(action);
    739                                    break;
    740                               case CLP_PARAM_STR_INTPRINT:
    741                                    printMode = action;
    742                                    break;
    743                               case CLP_PARAM_STR_KEEPNAMES:
    744                                    keepImportNames = 1 - action;
    745                                    break;
    746                               case CLP_PARAM_STR_PRESOLVE:
    747                                    if (action == 0)
    748                                         preSolve = DEFAULT_PRESOLVE_PASSES;
    749                                    else if (action == 1)
    750                                         preSolve = 0;
    751                                    else if (action == 2)
    752                                         preSolve = 10;
    753                                    else
    754                                         preSolveFile = true;
    755                                    break;
    756                               case CLP_PARAM_STR_PFI:
    757                                    thisModel->factorization()->setForrestTomlin(action == 0);
    758                                    break;
    759                               case CLP_PARAM_STR_FACTORIZATION:
    760                                    models[iModel].factorization()->forceOtherFactorization(action);
    761 #ifdef ABC_INHERIT
    762                                    thisModel->factorization()->forceOtherFactorization(action);
    763 #endif
    764                                    break;
    765                               case CLP_PARAM_STR_CRASH:
    766                                    doCrash = action;
    767                                    break;
    768                               case CLP_PARAM_STR_VECTOR:
    769                                    doVector = action;
    770                                    break;
    771                               case CLP_PARAM_STR_MESSAGES:
    772                                    models[iModel].messageHandler()->setPrefix(action != 0);
    773 #ifdef ABC_INHERIT
    774                                    thisModel->messageHandler()->setPrefix(action != 0);
    775 #endif
    776                                    break;
    777                               case CLP_PARAM_STR_CHOLESKY:
    778                                    choleskyType = action;
    779                                    break;
    780                               case CLP_PARAM_STR_GAMMA:
    781                                    gamma = action;
    782                                    break;
    783                               case CLP_PARAM_STR_BARRIERSCALE:
    784                                    scaleBarrier = action;
    785                                    break;
    786                               case CLP_PARAM_STR_KKT:
    787                                    doKKT = action;
    788                                    break;
    789                               case CLP_PARAM_STR_CROSSOVER:
    790                                    crossover = action;
    791                                    break;
    792                               default:
    793                                    //abort();
    794                                    break;
    795                               }
    796                          }
    797                     } else {
    798                          // action
    799                          if (type == CLP_PARAM_ACTION_EXIT) {
    800 #ifdef COIN_HAS_ASL
    801                            if (usingAmpl) {
    802                              writeAmpl(&info);
    803                              freeArrays2(&info);
    804                              freeArgs(&info);
    805                            }
    806 #endif
    807                            break; // stop all
    808                          }
    809                          switch (type) {
    810                          case CLP_PARAM_ACTION_DUALSIMPLEX:
    811                          case CLP_PARAM_ACTION_PRIMALSIMPLEX:
    812                          case CLP_PARAM_ACTION_EITHERSIMPLEX:
    813                          case CLP_PARAM_ACTION_BARRIER:
    814                               // synonym for dual
    815                          case CBC_PARAM_ACTION_BAB:
    816                               if (goodModels[iModel]) {
    817                                 if (type==CLP_PARAM_ACTION_EITHERSIMPLEX||
    818                                     type==CBC_PARAM_ACTION_BAB)
    819                                   models[iModel].setMoreSpecialOptions(16384|
    820                                                                        models[iModel].moreSpecialOptions());
    821                                    double objScale =
    822                                         parameters[whichParam(CLP_PARAM_DBL_OBJSCALE2, numberParameters, parameters)].doubleValue();
    823                                    if (objScale != 1.0) {
    824                                         int iColumn;
    825                                         int numberColumns = models[iModel].numberColumns();
    826                                         double * dualColumnSolution =
    827                                              models[iModel].dualColumnSolution();
    828                                         ClpObjective * obj = models[iModel].objectiveAsObject();
    829                                         assert(dynamic_cast<ClpLinearObjective *> (obj));
    830                                         double offset;
    831                                         double * objective = obj->gradient(NULL, NULL, offset, true);
    832                                         for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    833                                              dualColumnSolution[iColumn] *= objScale;
    834                                              objective[iColumn] *= objScale;;
    835                                         }
    836                                         int iRow;
    837                                         int numberRows = models[iModel].numberRows();
    838                                         double * dualRowSolution =
    839                                              models[iModel].dualRowSolution();
    840                                         for (iRow = 0; iRow < numberRows; iRow++)
    841                                              dualRowSolution[iRow] *= objScale;
    842                                         models[iModel].setObjectiveOffset(objScale * models[iModel].objectiveOffset());
    843                                    }
    844                                    ClpSolve::SolveType method;
    845                                    ClpSolve::PresolveType presolveType;
    846                                    ClpSolve solveOptions;
    847 #ifndef ABC_INHERIT
    848                                    ClpSimplex * model2 = models + iModel;
    849 #else
    850                                    AbcSimplex * model2 = models + iModel;
    851                                    if (type==CLP_PARAM_ACTION_EITHERSIMPLEX||
    852                                        type==CBC_PARAM_ACTION_BAB)
    853                                      solveOptions.setSpecialOption(3,0); // allow +-1
    854 #endif
    855                                    if (dualize==4) {
    856                                      solveOptions.setSpecialOption(4, 77);
    857                                      dualize=0;
    858                                    }
    859                                    if (dualize) {
    860                                         bool tryIt = true;
    861                                         double fractionColumn = 1.0;
    862                                         double fractionRow = 1.0;
    863                                         if (dualize == 3) {
    864                                              dualize = 1;
    865                                              int numberColumns = model2->numberColumns();
    866                                              int numberRows = model2->numberRows();
    867 #ifndef ABC_INHERIT
    868                                              if (numberRows < 50000 || 5 * numberColumns > numberRows) {
    869 #else
    870                                              if (numberRows < 500 || 4 * numberColumns > numberRows) {
    871 #endif
    872                                                   tryIt = false;
    873                                              } else {
    874                                                   fractionColumn = 0.1;
    875                                                   fractionRow = 0.3;
    876                                              }
    877                                         }
    878                                         if (tryIt) {
    879                                           ClpSimplex * thisModel=model2;
    880                                              thisModel = static_cast<ClpSimplexOther *> (thisModel)->dualOfModel(fractionRow, fractionColumn);
    881                                              if (thisModel) {
    882                                                   printf("Dual of model has %d rows and %d columns\n",
    883                                                          thisModel->numberRows(), thisModel->numberColumns());
    884                                                   thisModel->setOptimizationDirection(1.0);
    885 #ifndef ABC_INHERIT
    886                                                   model2=thisModel;
    887 #else
    888                                                   int abcState=model2->abcState();
    889                                                   model2=new AbcSimplex(*thisModel);
    890                                                   model2->setAbcState(abcState);
    891                                                   delete thisModel;
    892 #endif
    893                                              } else {
    894                                                   thisModel = models + iModel;
    895                                                   dualize = 0;
    896                                              }
    897                                         } else {
    898                                              dualize = 0;
    899                                         }
    900                                    }
    901                                    if (preSolveFile)
    902                                         presolveOptions |= 0x40000000;
    903                                    solveOptions.setPresolveActions(presolveOptions);
    904                                    solveOptions.setSubstitution(substitution);
    905                                    if (preSolve != DEFAULT_PRESOLVE_PASSES && preSolve) {
    906                                         presolveType = ClpSolve::presolveNumber;
    907                                         if (preSolve < 0) {
    908                                              preSolve = - preSolve;
    909                                              if (preSolve <= 100) {
    910                                                   presolveType = ClpSolve::presolveNumber;
    911                                                   printf("Doing %d presolve passes - picking up non-costed slacks\n",
    912                                                          preSolve);
    913                                                   solveOptions.setDoSingletonColumn(true);
    914                                              } else {
    915                                                   preSolve -= 100;
    916                                                   presolveType = ClpSolve::presolveNumberCost;
    917                                                   printf("Doing %d presolve passes - picking up costed slacks\n",
    918                                                          preSolve);
    919                                              }
    920                                         }
    921                                    } else if (preSolve) {
    922                                         presolveType = ClpSolve::presolveOn;
    923                                    } else {
    924                                         presolveType = ClpSolve::presolveOff;
    925                                    }
    926                                    solveOptions.setPresolveType(presolveType, preSolve);
    927                                    if (type == CLP_PARAM_ACTION_DUALSIMPLEX ||
    928                                              type == CBC_PARAM_ACTION_BAB) {
    929                                         method = ClpSolve::useDual;
    930                                    } else if (type == CLP_PARAM_ACTION_PRIMALSIMPLEX) {
    931                                         method = ClpSolve::usePrimalorSprint;
    932                                    } else if (type == CLP_PARAM_ACTION_EITHERSIMPLEX) {
    933                                         method = ClpSolve::automatic;
    934                                    } else {
    935                                         method = ClpSolve::useBarrier;
    936 #ifdef ABC_INHERIT
    937                                         if (doIdiot > 0)
    938                                              solveOptions.setSpecialOption(1, 2, doIdiot); // dense threshold
    939 #endif
    940                                         if (crossover == 1) {
    941                                              method = ClpSolve::useBarrierNoCross;
    942                                         } else if (crossover == 2) {
    943                                              ClpObjective * obj = models[iModel].objectiveAsObject();
    944                                              if (obj->type() > 1) {
    945                                                   method = ClpSolve::useBarrierNoCross;
    946                                                   presolveType = ClpSolve::presolveOff;
    947                                                   solveOptions.setPresolveType(presolveType, preSolve);
    948                                              }
    949                                         }
    950                                    }
    951                                    solveOptions.setSolveType(method);
    952                                    solveOptions.setSpecialOption(5, printOptions & 1);
    953                                    if (doVector) {
    954                                         ClpMatrixBase * matrix = models[iModel].clpMatrix();
    955                                         if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
    956                                              ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
    957                                              clpMatrix->makeSpecialColumnCopy();
    958                                         }
    959                                    }
    960                                    if (method == ClpSolve::useDual) {
    961                                         // dual
    962                                         if (doCrash)
    963                                              solveOptions.setSpecialOption(0, 1, doCrash); // crash
    964                                         else if (doIdiot)
    965                                              solveOptions.setSpecialOption(0, 2, doIdiot); // possible idiot
    966                                    } else if (method == ClpSolve::usePrimalorSprint) {
    967                                         // primal
    968                                         // if slp turn everything off
    969                                         if (slpValue > 0) {
    970                                              doCrash = false;
    971                                              doSprint = 0;
    972                                              doIdiot = -1;
    973                                              solveOptions.setSpecialOption(1, 10, slpValue); // slp
    974                                              method = ClpSolve::usePrimal;
    975                                         }
    976                                         if (doCrash) {
    977                                              solveOptions.setSpecialOption(1, 1, doCrash); // crash
    978                                         } else if (doSprint > 0) {
    979                                              // sprint overrides idiot
    980                                              solveOptions.setSpecialOption(1, 3, doSprint); // sprint
    981                                         } else if (doIdiot > 0) {
    982                                              solveOptions.setSpecialOption(1, 2, doIdiot); // idiot
    983                                         } else if (slpValue <= 0) {
    984                                              if (doIdiot == 0) {
    985                                                   if (doSprint == 0)
    986                                                        solveOptions.setSpecialOption(1, 4); // all slack
    987                                                   else
    988                                                        solveOptions.setSpecialOption(1, 9); // all slack or sprint
    989                                              } else {
    990                                                   if (doSprint == 0)
    991                                                        solveOptions.setSpecialOption(1, 8); // all slack or idiot
    992                                                   else
    993                                                        solveOptions.setSpecialOption(1, 7); // initiative
    994                                              }
    995                                         }
    996                                         if (basisHasValues == -1)
    997                                              solveOptions.setSpecialOption(1, 11); // switch off values
    998                                    } else if (method == ClpSolve::useBarrier || method == ClpSolve::useBarrierNoCross) {
    999                                         int barrierOptions = choleskyType;
    1000                                         if (scaleBarrier) {
    1001                                              if ((scaleBarrier & 1) != 0)
    1002                                                   barrierOptions |= 8;
    1003                                              barrierOptions |= 2048 * (scaleBarrier >> 1);
    1004                                         }
    1005                                         if (doKKT)
    1006                                              barrierOptions |= 16;
    1007                                         if (gamma)
    1008                                              barrierOptions |= 32 * gamma;
    1009                                         if (crossover == 3)
    1010                                              barrierOptions |= 256; // try presolve in crossover
    1011                                         solveOptions.setSpecialOption(4, barrierOptions);
    1012                                    }
    1013                                    int status;
    1014                                    if (cppValue >= 0) {
    1015                                         // generate code
    1016                                         FILE * fp = fopen("user_driver.cpp", "w");
    1017                                         if (fp) {
    1018                                              // generate enough to do solveOptions
    1019                                              model2->generateCpp(fp);
    1020                                              solveOptions.generateCpp(fp);
    1021                                              fclose(fp);
    1022                                              // now call generate code
    1023                                              generateCode("user_driver.cpp", cppValue);
    1024                                         } else {
    1025                                              std::cout << "Unable to open file user_driver.cpp" << std::endl;
    1026                                         }
    1027                                    }
    1028 #ifdef CLP_MULTIPLE_FACTORIZATIONS
    1029                                    int denseCode = parameters[whichParam(CBC_PARAM_INT_DENSE, numberParameters, parameters)].intValue();
    1030                                    if (denseCode!=-1)
    1031                                      model2->factorization()->setGoDenseThreshold(denseCode);
    1032                                    int smallCode = parameters[whichParam(CBC_PARAM_INT_SMALLFACT, numberParameters, parameters)].intValue();
    1033                                    if (smallCode!=-1)
    1034                                      model2->factorization()->setGoSmallThreshold(smallCode);
    1035                                    model2->factorization()->goDenseOrSmall(model2->numberRows());
    1036 #endif
    1037                                    try {
    1038                                         status = model2->initialSolve(solveOptions);
    1039 #ifdef COIN_HAS_ASL
    1040                             if (usingAmpl) {
    1041                                 double value = model2->getObjValue() * model2->getObjSense();
    1042                                 char buf[300];
    1043                                 int pos = 0;
    1044                                 int iStat = model2->status();
    1045                                 if (iStat == 0) {
    1046                                     pos += sprintf(buf + pos, "optimal," );
    1047                                 } else if (iStat == 1) {
    1048                                     // infeasible
    1049                                     pos += sprintf(buf + pos, "infeasible,");
    1050                                 } else if (iStat == 2) {
    1051                                     // unbounded
    1052                                     pos += sprintf(buf + pos, "unbounded,");
    1053                                 } else if (iStat == 3) {
    1054                                     pos += sprintf(buf + pos, "stopped on iterations or time,");
    1055                                 } else if (iStat == 4) {
    1056                                     iStat = 7;
    1057                                     pos += sprintf(buf + pos, "stopped on difficulties,");
    1058                                 } else if (iStat == 5) {
    1059                                     iStat = 3;
    1060                                     pos += sprintf(buf + pos, "stopped on ctrl-c,");
    1061                                 } else if (iStat == 6) {
    1062                                     // bab infeasible
    1063                                     pos += sprintf(buf + pos, "integer infeasible,");
    1064                                     iStat = 1;
    1065                                 } else {
    1066                                     pos += sprintf(buf + pos, "status unknown,");
    1067                                     iStat = 6;
    1068                                 }
    1069                                 info.problemStatus = iStat;
    1070                                 info.objValue = value;
    1071                                 pos += sprintf(buf + pos, " objective %.*g", ampl_obj_prec(),
    1072                                                value);
    1073                                 sprintf(buf + pos, "\n%d iterations",
    1074                                         model2->getIterationCount());
    1075                                 free(info.primalSolution);
    1076                                 int numberColumns = model2->numberColumns();
    1077                                 info.primalSolution = reinterpret_cast<double *> (malloc(numberColumns * sizeof(double)));
    1078                                 CoinCopyN(model2->primalColumnSolution(), numberColumns, info.primalSolution);
    1079                                 int numberRows = model2->numberRows();
    1080                                 free(info.dualSolution);
    1081                                 info.dualSolution = reinterpret_cast<double *> (malloc(numberRows * sizeof(double)));
    1082                                 CoinCopyN(model2->dualRowSolution(), numberRows, info.dualSolution);
    1083                                 CoinWarmStartBasis * basis = model2->getBasis();
    1084                                 free(info.rowStatus);
    1085                                 info.rowStatus = reinterpret_cast<int *> (malloc(numberRows * sizeof(int)));
    1086                                 free(info.columnStatus);
    1087                                 info.columnStatus = reinterpret_cast<int *> (malloc(numberColumns * sizeof(int)));
    1088                                 // Put basis in
    1089                                 int i;
    1090                                 // free,basic,ub,lb are 0,1,2,3
    1091                                 for (i = 0; i < numberRows; i++) {
    1092                                     CoinWarmStartBasis::Status status = basis->getArtifStatus(i);
    1093                                     info.rowStatus[i] = status;
    1094                                 }
    1095                                 for (i = 0; i < numberColumns; i++) {
    1096                                     CoinWarmStartBasis::Status status = basis->getStructStatus(i);
    1097                                     info.columnStatus[i] = status;
    1098                                 }
    1099                                 // put buffer into info
    1100                                 strcpy(info.buffer, buf);
    1101                                 delete basis;
    1102                             }
    1103 #endif
    1104 #ifndef NDEBUG
    1105                                         // if infeasible check ray
    1106                                         if (model2->status()==1) {
    1107                                           ClpSimplex * simplex = model2;
    1108                                           if(simplex->ray()) {
    1109                                             // make sure we use non-scaled versions
    1110                                             ClpPackedMatrix * saveMatrix = simplex->swapScaledMatrix(NULL);
    1111                                             double * saveScale = simplex->swapRowScale(NULL);
    1112                                             // could use existing arrays
    1113                                             int numberRows=simplex->numberRows();
    1114                                             int numberColumns=simplex->numberColumns();
    1115                                             double * farkas = new double [2*numberColumns+numberRows];
    1116                                             double * bound = farkas + numberColumns;
    1117                                             double * effectiveRhs = bound + numberColumns;
    1118                                             // get ray as user would
    1119                                             double * ray = simplex->infeasibilityRay();
    1120                                             // get farkas row
    1121                                             memset(farkas,0,(2*numberColumns+numberRows)*sizeof(double));
    1122                                             simplex->transposeTimes(-1.0,ray,farkas);
    1123                                             // Put nonzero bounds in bound
    1124                                             const double * columnLower = simplex->columnLower();
    1125                                             const double * columnUpper = simplex->columnUpper();
    1126                                             int numberBad=0;
    1127                                             for (int i=0;i<numberColumns;i++) {
    1128                                               double value = farkas[i];
    1129                                               double boundValue=0.0;
    1130                                               if (simplex->getStatus(i)==ClpSimplex::basic) {
    1131                                                 // treat as zero if small
    1132                                                 if (fabs(value)<1.0e-8) {
    1133                                                   value=0.0;
    1134                                                   farkas[i]=0.0;
    1135                                                 }
    1136                                                 if (value) {
    1137                                                   //printf("basic %d direction %d farkas %g\n",
    1138                                                   //       i,simplex->directionOut(),value);
    1139                                                   if (value<0.0)
    1140                                                     boundValue=columnLower[i];
    1141                                                   else
    1142                                                     boundValue=columnUpper[i];
    1143                                                 }
    1144                                               } else if (fabs(value)>1.0e-10) {
    1145                                                 if (value<0.0)
    1146                                                   boundValue=columnLower[i];
    1147                                                 else
    1148                                                   boundValue=columnUpper[i];
    1149                                               }
    1150                                               bound[i]=boundValue;
    1151                                               if (fabs(boundValue)>1.0e10)
    1152                                                 numberBad++;
    1153                                             }
    1154                                             const double * rowLower = simplex->rowLower();
    1155                                             const double * rowUpper = simplex->rowUpper();
    1156                                             //int pivotRow = simplex->spareIntArray_[3];
    1157                                             //bool badPivot=pivotRow<0;
    1158                                             for (int i=0;i<numberRows;i++) {
    1159                                               double value = ray[i];
    1160                                               double rhsValue=0.0;
    1161                                               if (simplex->getRowStatus(i)==ClpSimplex::basic) {
    1162                                                 // treat as zero if small
    1163                                                 if (fabs(value)<1.0e-8) {
    1164                                                   value=0.0;
    1165                                                   ray[i]=0.0;
    1166                                                 }
    1167                                                 if (value) {
    1168                                                   //printf("row basic %d direction %d ray %g\n",
    1169                                                   //       i,simplex->directionOut(),value);
    1170                                                   if (value<0.0)
    1171                                                     rhsValue=rowLower[i];
    1172                                                   else
    1173                                                     rhsValue=rowUpper[i];
    1174                                                 }
    1175                                               } else if (fabs(value)>1.0e-10) {
    1176                                                 if (value<0.0)
    1177                                                   rhsValue=rowLower[i];
    1178                                                 else
    1179                                                   rhsValue=rowUpper[i];
    1180                                               }
    1181                                               effectiveRhs[i]=rhsValue;
    1182                                               if (fabs(effectiveRhs[i])>1.0e10)
    1183                                                 printf("Large rhs row %d %g\n",
    1184                                                        i,effectiveRhs[i]);
    1185                                             }
    1186                                             simplex->times(-1.0,bound,effectiveRhs);
    1187                                             double bSum=0.0;
    1188                                             for (int i=0;i<numberRows;i++) {
    1189                                               bSum += effectiveRhs[i]*ray[i];
    1190                                               if (fabs(effectiveRhs[i])>1.0e10)
    1191                                                 printf("Large rhs row %d %g after\n",
    1192                                                        i,effectiveRhs[i]);
    1193                                             }
    1194                                             if (numberBad||bSum>1.0e-6) {
    1195                                               printf("Bad infeasibility ray %g  - %d bad\n",
    1196                                                      bSum,numberBad);
    1197                                             } else {
    1198                                               //printf("Good ray - infeasibility %g\n",
    1199                                               //     -bSum);
    1200                                             }
    1201                                             delete [] ray;
    1202                                             delete [] farkas;
    1203                                             simplex->swapRowScale(saveScale);
    1204                                             simplex->swapScaledMatrix(saveMatrix);
    1205                                           } else {
    1206                                             //printf("No dual ray\n");
    1207                                           }
    1208                                         }
    1209 #endif
    1210                                    } catch (CoinError e) {
    1211                                         e.print();
    1212                                         status = -1;
    1213                                    }
    1214                                    if (dualize) {
    1215                                      ClpSimplex * thisModel=models+iModel;
    1216                                         int returnCode = static_cast<ClpSimplexOther *> (thisModel)->restoreFromDual(model2);
    1217                                         if (model2->status() == 3)
    1218                                              returnCode = 0;
    1219                                         delete model2;
    1220                                         if (returnCode && dualize != 2) {
    1221                                              currentModel = models + iModel;
    1222                                              // register signal handler
    1223                                              signal(SIGINT, signal_handler);
    1224                                              thisModel->primal(1);
    1225                                              currentModel = NULL;
    1226                                         }
    1227                                         // switch off (user can switch back on)
    1228                                         parameters[whichParam(CLP_PARAM_INT_DUALIZE,
    1229                                                               numberParameters, parameters)].setIntValue(dualize);
    1230                                    }
    1231                                    if (status >= 0)
    1232                                         basisHasValues = 1;
    1233                               } else {
    1234                                    std::cout << "** Current model not valid" << std::endl;
    1235                               }
    1236                               break;
    1237                          case CLP_PARAM_ACTION_STATISTICS:
    1238                               if (goodModels[iModel]) {
    1239                                    // If presolve on look at presolved
    1240                                    bool deleteModel2 = false;
    1241                                    ClpSimplex * model2 = models + iModel;
    1242                                    if (preSolve) {
    1243                                         ClpPresolve pinfo;
    1244                                         int presolveOptions2 = presolveOptions&~0x40000000;
    1245                                         if ((presolveOptions2 & 0xffff) != 0)
    1246                                              pinfo.setPresolveActions(presolveOptions2);
    1247                                         pinfo.setSubstitution(substitution);
    1248                                         if ((printOptions & 1) != 0)
    1249                                              pinfo.statistics();
    1250                                         double presolveTolerance =
    1251                                              parameters[whichParam(CLP_PARAM_DBL_PRESOLVETOLERANCE, numberParameters, parameters)].doubleValue();
    1252                                         model2 =
    1253                                              pinfo.presolvedModel(models[iModel], presolveTolerance,
    1254                                                                   true, preSolve);
    1255                                         if (model2) {
    1256                                              printf("Statistics for presolved model\n");
    1257                                              deleteModel2 = true;
    1258                                         } else {
    1259                                              printf("Presolved model looks infeasible - will use unpresolved\n");
    1260                                              model2 = models + iModel;
    1261                                         }
    1262                                    } else {
    1263                                         printf("Statistics for unpresolved model\n");
    1264                                         model2 =  models + iModel;
    1265                                    }
    1266                                    statistics(models + iModel, model2);
    1267                                    if (deleteModel2)
    1268                                         delete model2;
    1269                               } else {
    1270                                    std::cout << "** Current model not valid" << std::endl;
    1271                               }
    1272                               break;
    1273                          case CLP_PARAM_ACTION_TIGHTEN:
    1274                               if (goodModels[iModel]) {
    1275                                    int numberInfeasibilities = models[iModel].tightenPrimalBounds();
    1276                                    if (numberInfeasibilities)
    1277                                         std::cout << "** Analysis indicates model infeasible" << std::endl;
    1278                               } else {
    1279                                    std::cout << "** Current model not valid" << std::endl;
    1280                               }
    1281                               break;
    1282                          case CLP_PARAM_ACTION_PLUSMINUS:
    1283                               if (goodModels[iModel]) {
    1284                                    ClpMatrixBase * saveMatrix = models[iModel].clpMatrix();
    1285                                    ClpPackedMatrix* clpMatrix =
    1286                                         dynamic_cast< ClpPackedMatrix*>(saveMatrix);
    1287                                    if (clpMatrix) {
    1288                                         ClpPlusMinusOneMatrix * newMatrix = new ClpPlusMinusOneMatrix(*(clpMatrix->matrix()));
    1289                                         if (newMatrix->getIndices()) {
    1290                                              models[iModel].replaceMatrix(newMatrix);
    1291                                              delete saveMatrix;
    1292                                              std::cout << "Matrix converted to +- one matrix" << std::endl;
    1293                                         } else {
    1294                                              std::cout << "Matrix can not be converted to +- 1 matrix" << std::endl;
    1295                                         }
    1296                                    } else {
    1297                                         std::cout << "Matrix not a ClpPackedMatrix" << std::endl;
    1298                                    }
    1299                               } else {
    1300                                    std::cout << "** Current model not valid" << std::endl;
    1301                               }
    1302                               break;
    1303                          case CLP_PARAM_ACTION_NETWORK:
    1304                               if (goodModels[iModel]) {
    1305                                    ClpMatrixBase * saveMatrix = models[iModel].clpMatrix();
    1306                                    ClpPackedMatrix* clpMatrix =
    1307                                         dynamic_cast< ClpPackedMatrix*>(saveMatrix);
    1308                                    if (clpMatrix) {
    1309                                         ClpNetworkMatrix * newMatrix = new ClpNetworkMatrix(*(clpMatrix->matrix()));
    1310                                         if (newMatrix->getIndices()) {
    1311                                              models[iModel].replaceMatrix(newMatrix);
    1312                                              delete saveMatrix;
    1313                                              std::cout << "Matrix converted to network matrix" << std::endl;
    1314                                         } else {
    1315                                              std::cout << "Matrix can not be converted to network matrix" << std::endl;
    1316                                         }
    1317                                    } else {
    1318                                         std::cout << "Matrix not a ClpPackedMatrix" << std::endl;
    1319                                    }
    1320                               } else {
    1321                                    std::cout << "** Current model not valid" << std::endl;
    1322                               }
    1323                               break;
    1324                          case CLP_PARAM_ACTION_IMPORT: {
    1325                               // get next field
    1326                               field = CoinReadGetString(argc, argv);
    1327                               if (field == "$") {
    1328                                    field = parameters[iParam].stringValue();
    1329                               } else if (field == "EOL") {
    1330                                    parameters[iParam].printString();
    1331                                    break;
    1332                               } else {
    1333                                    parameters[iParam].setStringValue(field);
    1334                               }
    1335                               std::string fileName;
    1336                               bool canOpen = false;
    1337                               // See if gmpl file
    1338                               int gmpl = 0;
    1339                               std::string gmplData;
    1340                               if (field == "-") {
    1341                                    // stdin
    1342                                    canOpen = true;
    1343                                    fileName = "-";
    1344                               } else {
    1345                                    // See if .lp
    1346                                    {
    1347                                         const char * c_name = field.c_str();
    1348                                         size_t length = strlen(c_name);
    1349                                         if (length > 3 && !strncmp(c_name + length - 3, ".lp", 3))
    1350                                              gmpl = -1; // .lp
    1351                                    }
    1352                                    bool absolutePath;
    1353                                    if (dirsep == '/') {
    1354                                         // non Windows (or cygwin)
    1355                                         absolutePath = (field[0] == '/');
    1356                                    } else {
    1357                                         //Windows (non cycgwin)
    1358                                         absolutePath = (field[0] == '\\');
    1359                                         // but allow for :
    1360                                         if (strchr(field.c_str(), ':'))
    1361                                              absolutePath = true;
    1362                                    }
    1363                                    if (absolutePath) {
    1364                                         fileName = field;
    1365                                         size_t length = field.size();
    1366                                         size_t percent = field.find('%');
    1367                                         if (percent < length && percent > 0) {
    1368                                              gmpl = 1;
    1369                                              fileName = field.substr(0, percent);
    1370                                              gmplData = field.substr(percent + 1);
    1371                                              if (percent < length - 1)
    1372                                                   gmpl = 2; // two files
    1373                                              printf("GMPL model file %s and data file %s\n",
    1374                                                     fileName.c_str(), gmplData.c_str());
    1375                                         }
    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                                         // See if gmpl (model & data) - or even lp file
    1388                                         size_t length = field.size();
    1389                                         size_t percent = field.find('%');
    1390                                         if (percent < length && percent > 0) {
    1391                                              gmpl = 1;
    1392                                              fileName = directory + field.substr(0, percent);
    1393                                              gmplData = directory + field.substr(percent + 1);
    1394                                              if (percent < length - 1)
    1395                                                   gmpl = 2; // two files
    1396                                              printf("GMPL model file %s and data file %s\n",
    1397                                                     fileName.c_str(), gmplData.c_str());
    1398                                         }
    1399                                    }
    1400                                    std::string name = fileName;
    1401                                    if (fileCoinReadable(name)) {
    1402                                         // can open - lets go for it
    1403                                         canOpen = true;
    1404                                         if (gmpl == 2) {
    1405                                              FILE *fp;
    1406                                              fp = fopen(gmplData.c_str(), "r");
    1407                                              if (fp) {
    1408                                                   fclose(fp);
    1409                                              } else {
    1410                                                   canOpen = false;
    1411                                                   std::cout << "Unable to open file " << gmplData << std::endl;
    1412                                              }
    1413                                         }
    1414                                    } else {
    1415                                         std::cout << "Unable to open file " << fileName << std::endl;
    1416                                    }
    1417                               }
    1418                               if (canOpen) {
    1419                                    int status;
    1420                                    if (!gmpl)
    1421                                         status = models[iModel].readMps(fileName.c_str(),
    1422                                                                         keepImportNames != 0,
    1423                                                                         allowImportErrors != 0);
    1424                                    else if (gmpl > 0)
    1425                                         status = models[iModel].readGMPL(fileName.c_str(),
    1426                                                                          (gmpl == 2) ? gmplData.c_str() : NULL,
    1427                                                                          keepImportNames != 0);
    1428                                    else
    1429 #ifdef KILL_ZERO_READLP
    1430                                      status = models[iModel].readLp(fileName.c_str(), models[iModel].getSmallElementValue());
    1431 #else
    1432                                      status = models[iModel].readLp(fileName.c_str(), 1.0e-12);
    1433 #endif
    1434                                    if (!status || (status > 0 && allowImportErrors)) {
    1435                                         goodModels[iModel] = true;
    1436                                         // sets to all slack (not necessary?)
    1437                                         thisModel->createStatus();
    1438                                         time2 = CoinCpuTime();
    1439                                         totalTime += time2 - time1;
    1440                                         time1 = time2;
    1441                                         // Go to canned file if just input file
    1442                                         if (CbcOrClpRead_mode == 2 && argc == 2) {
    1443                                              // only if ends .mps
    1444                                              char * find = const_cast<char *>(strstr(fileName.c_str(), ".mps"));
    1445                                              if (find && find[4] == '\0') {
    1446                                                   find[1] = 'p';
    1447                                                   find[2] = 'a';
    1448                                                   find[3] = 'r';
    1449                                                   FILE *fp = fopen(fileName.c_str(), "r");
    1450                                                   if (fp) {
    1451                                                        CbcOrClpReadCommand = fp; // Read from that file
    1452                                                        CbcOrClpRead_mode = -1;
    1453                                                   }
    1454                                              }
    1455                                         }
    1456                                    } else {
    1457                                         // errors
    1458                                         std::cout << "There were " << status <<
    1459                                                   " errors on input" << std::endl;
    1460                                    }
    1461                               }
    1462                          }
    1463                          break;
    1464                          case CLP_PARAM_ACTION_EXPORT:
    1465                               if (goodModels[iModel]) {
    1466                                    double objScale =
    1467                                         parameters[whichParam(CLP_PARAM_DBL_OBJSCALE2, numberParameters, parameters)].doubleValue();
    1468                                    if (objScale != 1.0) {
    1469                                         int iColumn;
    1470                                         int numberColumns = models[iModel].numberColumns();
    1471                                         double * dualColumnSolution =
    1472                                              models[iModel].dualColumnSolution();
    1473                                         ClpObjective * obj = models[iModel].objectiveAsObject();
    1474                                         assert(dynamic_cast<ClpLinearObjective *> (obj));
    1475                                         double offset;
    1476                                         double * objective = obj->gradient(NULL, NULL, offset, true);
    1477                                         for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1478                                              dualColumnSolution[iColumn] *= objScale;
    1479                                              objective[iColumn] *= objScale;;
    1480                                         }
    1481                                         int iRow;
    1482                                         int numberRows = models[iModel].numberRows();
    1483                                         double * dualRowSolution =
    1484                                              models[iModel].dualRowSolution();
    1485                                         for (iRow = 0; iRow < numberRows; iRow++)
    1486                                              dualRowSolution[iRow] *= objScale;
    1487                                         models[iModel].setObjectiveOffset(objScale * models[iModel].objectiveOffset());
    1488                                    }
    1489                                    // get next field
    1490                                    field = CoinReadGetString(argc, argv);
    1491                                    if (field == "$") {
    1492                                         field = parameters[iParam].stringValue();
    1493                                    } else if (field == "EOL") {
    1494                                         parameters[iParam].printString();
    1495                                         break;
    1496                                    } else {
    1497                                         parameters[iParam].setStringValue(field);
    1498                                    }
    1499                                    std::string fileName;
    1500                                    bool canOpen = false;
    1501                                    if (field[0] == '/' || field[0] == '\\') {
    1502                                         fileName = field;
    1503                                    } else if (field[0] == '~') {
    1504                                         char * environVar = getenv("HOME");
    1505                                         if (environVar) {
    1506                                              std::string home(environVar);
    1507                                              field = field.erase(0, 1);
    1508                                              fileName = home + field;
    1509                                         } else {
    1510                                              fileName = field;
    1511                                         }
    1512                                    } else {
    1513                                         fileName = directory + field;
    1514                                    }
    1515                                    FILE *fp = fopen(fileName.c_str(), "w");
    1516                                    if (fp) {
    1517                                         // can open - lets go for it
    1518                                         fclose(fp);
    1519                                         canOpen = true;
    1520                                    } else {
    1521                                         std::cout << "Unable to open file " << fileName << std::endl;
    1522                                    }
    1523                                    if (canOpen) {
    1524                                         // If presolve on then save presolved
    1525                                         bool deleteModel2 = false;
    1526                                         ClpSimplex * model2 = models + iModel;
    1527                                         if (dualize && dualize < 3) {
    1528                                              model2 = static_cast<ClpSimplexOther *> (model2)->dualOfModel();
    1529                                              printf("Dual of model has %d rows and %d columns\n",
    1530                                                     model2->numberRows(), model2->numberColumns());
    1531                                              model2->setOptimizationDirection(1.0);
    1532                                              preSolve = 0; // as picks up from model
    1533                                         }
    1534                                         if (preSolve) {
    1535                                              ClpPresolve pinfo;
    1536                                              int presolveOptions2 = presolveOptions&~0x40000000;
    1537                                              if ((presolveOptions2 & 0xffff) != 0)
    1538                                                   pinfo.setPresolveActions(presolveOptions2);
    1539                                              pinfo.setSubstitution(substitution);
    1540                                              if ((printOptions & 1) != 0)
    1541                                                   pinfo.statistics();
    1542                                              double presolveTolerance =
    1543                                                   parameters[whichParam(CLP_PARAM_DBL_PRESOLVETOLERANCE, numberParameters, parameters)].doubleValue();
    1544                                              model2 =
    1545                                                   pinfo.presolvedModel(models[iModel], presolveTolerance,
    1546                                                                        true, preSolve, false, false);
    1547                                              if (model2) {
    1548                                                   printf("Saving presolved model on %s\n",
    1549                                                          fileName.c_str());
    1550                                                   deleteModel2 = true;
    1551                                              } else {
    1552                                                   printf("Presolved model looks infeasible - saving original on %s\n",
    1553                                                          fileName.c_str());
    1554                                                   deleteModel2 = false;
    1555                                                   model2 = models + iModel;
    1556 
    1557                                              }
    1558                                         } else {
    1559                                              printf("Saving model on %s\n",
    1560                                                     fileName.c_str());
    1561                                         }
    1562 #if 0
    1563                                         // Convert names
    1564                                         int iRow;
    1565                                         int numberRows = model2->numberRows();
    1566                                         int iColumn;
    1567                                         int numberColumns = model2->numberColumns();
    1568 
    1569                                         char ** rowNames = NULL;
    1570                                         char ** columnNames = NULL;
    1571                                         if (model2->lengthNames()) {
    1572                                              rowNames = new char * [numberRows];
    1573                                              for (iRow = 0; iRow < numberRows; iRow++) {
    1574                                                   rowNames[iRow] =
    1575                                                        CoinStrdup(model2->rowName(iRow).c_str());
    1576 #ifdef STRIPBLANKS
    1577                                                   char * xx = rowNames[iRow];
    1578                                                   int i;
    1579                                                   int length = strlen(xx);
    1580                                                   int n = 0;
    1581                                                   for (i = 0; i < length; i++) {
    1582                                                        if (xx[i] != ' ')
    1583                                                             xx[n++] = xx[i];
    1584                                                   }
    1585                                                   xx[n] = '\0';
    1586 #endif
    1587                                              }
    1588 
    1589                                              columnNames = new char * [numberColumns];
    1590                                              for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1591                                                   columnNames[iColumn] =
    1592                                                        CoinStrdup(model2->columnName(iColumn).c_str());
    1593 #ifdef STRIPBLANKS
    1594                                                   char * xx = columnNames[iColumn];
    1595                                                   int i;
    1596                                                   int length = strlen(xx);
    1597                                                   int n = 0;
    1598                                                   for (i = 0; i < length; i++) {
    1599                                                        if (xx[i] != ' ')
    1600                                                             xx[n++] = xx[i];
    1601                                                   }
    1602                                                   xx[n] = '\0';
    1603 #endif
    1604                                              }
    1605                                         }
    1606                                         CoinMpsIO writer;
    1607                                         writer.setMpsData(*model2->matrix(), COIN_DBL_MAX,
    1608                                                           model2->getColLower(), model2->getColUpper(),
    1609                                                           model2->getObjCoefficients(),
    1610                                                           (const char*) 0 /*integrality*/,
    1611                                                           model2->getRowLower(), model2->getRowUpper(),
    1612                                                           columnNames, rowNames);
    1613                                         // Pass in array saying if each variable integer
    1614                                         writer.copyInIntegerInformation(model2->integerInformation());
    1615                                         writer.setObjectiveOffset(model2->objectiveOffset());
    1616                                         writer.writeMps(fileName.c_str(), 0, 1, 1);
    1617                                         if (rowNames) {
    1618                                              for (iRow = 0; iRow < numberRows; iRow++) {
    1619                                                   free(rowNames[iRow]);
    1620                                              }
    1621                                              delete [] rowNames;
    1622                                              for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1623                                                   free(columnNames[iColumn]);
    1624                                              }
    1625                                              delete [] columnNames;
    1626                                         }
    1627 #else
    1628                                         model2->writeMps(fileName.c_str(), (outputFormat - 1) / 2, 1 + ((outputFormat - 1) & 1));
    1629 #endif
    1630                                         if (deleteModel2)
    1631                                              delete model2;
    1632                                         time2 = CoinCpuTime();
    1633                                         totalTime += time2 - time1;
    1634                                         time1 = time2;
    1635                                    }
    1636                               } else {
    1637                                    std::cout << "** Current model not valid" << std::endl;
    1638                               }
    1639                               break;
    1640                          case CLP_PARAM_ACTION_BASISIN:
    1641                               if (goodModels[iModel]) {
    1642                                    // get next field
    1643                                    field = CoinReadGetString(argc, argv);
    1644                                    if (field == "$") {
    1645                                         field = parameters[iParam].stringValue();
    1646                                    } else if (field == "EOL") {
    1647                                         parameters[iParam].printString();
    1648                                         break;
    1649                                    } else {
    1650                                         parameters[iParam].setStringValue(field);
    1651                                    }
    1652                                    std::string fileName;
    1653                                    bool canOpen = false;
    1654                                    if (field == "-") {
    1655                                         // stdin
    1656                                         canOpen = true;
    1657                                         fileName = "-";
    1658                                    } else {
    1659                                         if (field[0] == '/' || field[0] == '\\') {
    1660                                              fileName = field;
    1661                                         } else if (field[0] == '~') {
    1662                                              char * environVar = getenv("HOME");
    1663                                              if (environVar) {
    1664                                                   std::string home(environVar);
    1665                                                   field = field.erase(0, 1);
    1666                                                   fileName = home + field;
    1667                                              } else {
    1668                                                   fileName = field;
    1669                                              }
    1670                                         } else {
    1671                                              fileName = directory + field;
    1672                                         }
    1673                                         FILE *fp = fopen(fileName.c_str(), "r");
    1674                                         if (fp) {
    1675                                              // can open - lets go for it
    1676                                              fclose(fp);
    1677                                              canOpen = true;
    1678                                         } else {
    1679                                              std::cout << "Unable to open file " << fileName << std::endl;
    1680                                         }
    1681                                    }
    1682                                    if (canOpen) {
    1683                                         int values = thisModel->readBasis(fileName.c_str());
    1684                                         if (values == 0)
    1685                                              basisHasValues = -1;
    1686                                         else
    1687                                              basisHasValues = 1;
    1688                                    }
    1689                               } else {
    1690                                    std::cout << "** Current model not valid" << std::endl;
    1691                               }
    1692                               break;
    1693                          case CLP_PARAM_ACTION_PRINTMASK:
    1694                               // get next field
    1695                          {
    1696                               std::string name = CoinReadGetString(argc, argv);
    1697                               if (name != "EOL") {
    1698                                    parameters[iParam].setStringValue(name);
    1699                                    printMask = name;
    1700                               } else {
    1701                                    parameters[iParam].printString();
    1702                               }
    1703                          }
    1704                          break;
    1705                          case CLP_PARAM_ACTION_BASISOUT:
    1706                               if (goodModels[iModel]) {
    1707                                    // get next field
    1708                                    field = CoinReadGetString(argc, argv);
    1709                                    if (field == "$") {
    1710                                         field = parameters[iParam].stringValue();
    1711                                    } else if (field == "EOL") {
    1712                                         parameters[iParam].printString();
    1713                                         break;
    1714                                    } else {
    1715                                         parameters[iParam].setStringValue(field);
    1716                                    }
    1717                                    std::string fileName;
    1718                                    bool canOpen = false;
    1719                                    if (field[0] == '/' || field[0] == '\\') {
    1720                                         fileName = field;
    1721                                    } else if (field[0] == '~') {
    1722                                         char * environVar = getenv("HOME");
    1723                                         if (environVar) {
    1724                                              std::string home(environVar);
    1725                                              field = field.erase(0, 1);
    1726                                              fileName = home + field;
    1727                                         } else {
    1728                                              fileName = field;
    1729                                         }
    1730                                    } else {
    1731                                         fileName = directory + field;
    1732                                    }
    1733                                    FILE *fp = fopen(fileName.c_str(), "w");
    1734                                    if (fp) {
    1735                                         // can open - lets go for it
    1736                                         fclose(fp);
    1737                                         canOpen = true;
    1738                                    } else {
    1739                                         std::cout << "Unable to open file " << fileName << std::endl;
    1740                                    }
    1741                                    if (canOpen) {
    1742                                         ClpSimplex * model2 = models + iModel;
    1743                                         model2->writeBasis(fileName.c_str(), outputFormat > 1, outputFormat - 2);
    1744                                         time2 = CoinCpuTime();
    1745                                         totalTime += time2 - time1;
    1746                                         time1 = time2;
    1747                                    }
    1748                               } else {
    1749                                    std::cout << "** Current model not valid" << std::endl;
    1750                               }
    1751                               break;
    1752                          case CLP_PARAM_ACTION_PARAMETRICS:
    1753                               if (goodModels[iModel]) {
    1754                                    // get next field
    1755                                    field = CoinReadGetString(argc, argv);
    1756                                    if (field == "$") {
    1757                                         field = parameters[iParam].stringValue();
    1758                                    } else if (field == "EOL") {
    1759                                         parameters[iParam].printString();
    1760                                         break;
    1761                                    } else {
    1762                                         parameters[iParam].setStringValue(field);
    1763                                    }
    1764                                    std::string fileName;
    1765                                    //bool canOpen = false;
    1766                                    if (field[0] == '/' || field[0] == '\\') {
    1767                                         fileName = field;
    1768                                    } else if (field[0] == '~') {
    1769                                         char * environVar = getenv("HOME");
    1770                                         if (environVar) {
    1771                                              std::string home(environVar);
    1772                                              field = field.erase(0, 1);
    1773                                              fileName = home + field;
    1774                                         } else {
    1775                                              fileName = field;
    1776                                         }
    1777                                    } else {
    1778                                         fileName = directory + field;
    1779                                    }
    1780                                    ClpSimplex * model2 = models + iModel;
    1781                                    static_cast<ClpSimplexOther *> (model2)->parametrics(fileName.c_str());
    1782                                    time2 = CoinCpuTime();
    1783                                    totalTime += time2 - time1;
    1784                                    time1 = time2;
    1785                               } else {
    1786                                    std::cout << "** Current model not valid" << std::endl;
    1787                               }
    1788                               break;
    1789                          case CLP_PARAM_ACTION_SAVE: {
    1790                               // get next field
    1791                               field = CoinReadGetString(argc, argv);
    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                               bool canOpen = false;
    1802                               if (field[0] == '/' || field[0] == '\\') {
    1803                                    fileName = field;
    1804                               } else if (field[0] == '~') {
    1805                                    char * environVar = getenv("HOME");
    1806                                    if (environVar) {
    1807                                         std::string home(environVar);
    1808                                         field = field.erase(0, 1);
    1809                                         fileName = home + field;
    1810                                    } else {
    1811                                         fileName = field;
    1812                                    }
    1813                               } else {
    1814                                    fileName = directory + field;
    1815                               }
    1816                               FILE *fp = fopen(fileName.c_str(), "wb");
    1817                               if (fp) {
    1818                                    // can open - lets go for it
    1819                                    fclose(fp);
    1820                                    canOpen = true;
    1821                               } else {
    1822                                    std::cout << "Unable to open file " << fileName << std::endl;
    1823                               }
    1824                               if (canOpen) {
    1825                                    int status;
    1826                                    // If presolve on then save presolved
    1827                                    bool deleteModel2 = false;
    1828                                    ClpSimplex * model2 = models + iModel;
    1829                                    if (preSolve) {
    1830                                         ClpPresolve pinfo;
    1831                                         double presolveTolerance =
    1832                                              parameters[whichParam(CLP_PARAM_DBL_PRESOLVETOLERANCE, numberParameters, parameters)].doubleValue();
    1833                                         model2 =
    1834                                              pinfo.presolvedModel(models[iModel], presolveTolerance,
    1835                                                                   false, preSolve);
    1836                                         if (model2) {
    1837                                              printf("Saving presolved model on %s\n",
    1838                                                     fileName.c_str());
    1839                                              deleteModel2 = true;
    1840                                         } else {
    1841                                              printf("Presolved model looks infeasible - saving original on %s\n",
    1842                                                     fileName.c_str());
    1843                                              deleteModel2 = false;
    1844                                              model2 = models + iModel;
    1845 
    1846                                         }
    1847                                    } else {
    1848                                         printf("Saving model on %s\n",
    1849                                                fileName.c_str());
    1850                                    }
    1851                                    status = model2->saveModel(fileName.c_str());
    1852                                    if (deleteModel2)
    1853                                         delete model2;
    1854                                    if (!status) {
    1855                                         goodModels[iModel] = true;
    1856                                         time2 = CoinCpuTime();
    1857                                         totalTime += time2 - time1;
    1858                                         time1 = time2;
    1859                                    } else {
    1860                                         // errors
    1861                                         std::cout << "There were errors on output" << std::endl;
    1862                                    }
    1863                               }
    1864                          }
    1865                          break;
    1866                          case CLP_PARAM_ACTION_RESTORE: {
    1867                               // get next field
    1868                               field = CoinReadGetString(argc, argv);
    1869                               if (field == "$") {
    1870                                    field = parameters[iParam].stringValue();
    1871                               } else if (field == "EOL") {
    1872                                    parameters[iParam].printString();
    1873                                    break;
    1874                               } else {
    1875                                    parameters[iParam].setStringValue(field);
    1876                               }
    1877                               std::string fileName;
    1878                               bool canOpen = false;
    1879                               if (field[0] == '/' || field[0] == '\\') {
    1880                                    fileName = field;
    1881                               } else if (field[0] == '~') {
    1882                                    char * environVar = getenv("HOME");
    1883                                    if (environVar) {
    1884                                         std::string home(environVar);
    1885                                         field = field.erase(0, 1);
    1886                                         fileName = home + field;
    1887                                    } else {
    1888                                         fileName = field;
    1889                                    }
    1890                               } else {
    1891                                    fileName = directory + field;
    1892                               }
    1893                               FILE *fp = fopen(fileName.c_str(), "rb");
    1894                               if (fp) {
    1895                                    // can open - lets go for it
    1896                                    fclose(fp);
    1897                                    canOpen = true;
    1898                               } else {
    1899                                    std::cout << "Unable to open file " << fileName << std::endl;
    1900                               }
    1901                               if (canOpen) {
    1902                                    int status = models[iModel].restoreModel(fileName.c_str());
    1903                                    if (!status) {
    1904                                         goodModels[iModel] = true;
    1905                                         time2 = CoinCpuTime();
    1906                                         totalTime += time2 - time1;
    1907                                         time1 = time2;
    1908                                    } else {
    1909                                         // errors
    1910                                         std::cout << "There were errors on input" << std::endl;
    1911                                    }
    1912                               }
    1913                          }
    1914                          break;
    1915                          case CLP_PARAM_ACTION_MAXIMIZE:
    1916                               models[iModel].setOptimizationDirection(-1);
    1917 #ifdef ABC_INHERIT
    1918                               thisModel->setOptimizationDirection(-1);
    1919 #endif
    1920                               break;
    1921                          case CLP_PARAM_ACTION_MINIMIZE:
    1922                               models[iModel].setOptimizationDirection(1);
    1923 #ifdef ABC_INHERIT
    1924                               thisModel->setOptimizationDirection(1);
    1925 #endif
    1926                               break;
    1927                          case CLP_PARAM_ACTION_ALLSLACK:
    1928                               thisModel->allSlackBasis(true);
    1929 #ifdef ABC_INHERIT
    1930                               models[iModel].allSlackBasis();
    1931 #endif
    1932                               break;
    1933                          case CLP_PARAM_ACTION_REVERSE:
    1934                               if (goodModels[iModel]) {
    1935                                    int iColumn;
    1936                                    int numberColumns = models[iModel].numberColumns();
    1937                                    double * dualColumnSolution =
    1938                                         models[iModel].dualColumnSolution();
    1939                                    ClpObjective * obj = models[iModel].objectiveAsObject();
    1940                                    assert(dynamic_cast<ClpLinearObjective *> (obj));
    1941                                    double offset;
    1942                                    double * objective = obj->gradient(NULL, NULL, offset, true);
    1943                                    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    1944                                         dualColumnSolution[iColumn] = -dualColumnSolution[iColumn];
    1945                                         objective[iColumn] = -objective[iColumn];
    1946                                    }
    1947                                    int iRow;
    1948                                    int numberRows = models[iModel].numberRows();
    1949                                    double * dualRowSolution =
    1950                                         models[iModel].dualRowSolution();
    1951                                    for (iRow = 0; iRow < numberRows; iRow++) {
    1952                                         dualRowSolution[iRow] = -dualRowSolution[iRow];
    1953                                    }
    1954                                    models[iModel].setObjectiveOffset(-models[iModel].objectiveOffset());
    1955                               }
    1956                               break;
    1957                          case CLP_PARAM_ACTION_DIRECTORY: {
    1958                               std::string name = CoinReadGetString(argc, argv);
    1959                               if (name != "EOL") {
    1960                                    size_t length = name.length();
    1961                                    if (length > 0 && name[length-1] == dirsep) {
    1962                                         directory = name;
    1963                                    } else {
    1964                                         directory = name + dirsep;
    1965                                    }
    1966                                    parameters[iParam].setStringValue(directory);
    1967                               } else {
    1968                                    parameters[iParam].printString();
    1969                               }
    1970                          }
    1971                          break;
    1972                          case CLP_PARAM_ACTION_DIRSAMPLE: {
    1973                               std::string name = CoinReadGetString(argc, argv);
    1974                               if (name != "EOL") {
    1975                                    size_t length = name.length();
    1976                                    if (length > 0 && name[length-1] == dirsep) {
    1977                                         dirSample = name;
    1978                                    } else {
    1979                                         dirSample = name + dirsep;
    1980                                    }
    1981                                    parameters[iParam].setStringValue(dirSample);
    1982                               } else {
    1983                                    parameters[iParam].printString();
    1984                               }
    1985                          }
    1986                          break;
    1987                          case CLP_PARAM_ACTION_DIRNETLIB: {
    1988                               std::string name = CoinReadGetString(argc, argv);
    1989                               if (name != "EOL") {
    1990                                    size_t length = name.length();
    1991                                    if (length > 0 && name[length-1] == dirsep) {
    1992                                         dirNetlib = name;
    1993                                    } else {
    1994                                         dirNetlib = name + dirsep;
    1995                                    }
    1996                                    parameters[iParam].setStringValue(dirNetlib);
    1997                               } else {
    1998                                    parameters[iParam].printString();
    1999                               }
    2000                          }
    2001                          break;
    2002                          case CBC_PARAM_ACTION_DIRMIPLIB: {
    2003                               std::string name = CoinReadGetString(argc, argv);
    2004                               if (name != "EOL") {
    2005                                    size_t length = name.length();
    2006                                    if (length > 0 && name[length-1] == dirsep) {
    2007                                         dirMiplib = name;
    2008                                    } else {
    2009                                         dirMiplib = name + dirsep;
    2010                                    }
    2011                                    parameters[iParam].setStringValue(dirMiplib);
    2012                               } else {
    2013                                    parameters[iParam].printString();
    2014                               }
    2015                          }
    2016                          break;
    2017                          case CLP_PARAM_ACTION_STDIN:
    2018                               CbcOrClpRead_mode = -1;
    2019                               break;
    2020                          case CLP_PARAM_ACTION_NETLIB_DUAL:
    2021                          case CLP_PARAM_ACTION_NETLIB_EITHER:
    2022                          case CLP_PARAM_ACTION_NETLIB_BARRIER:
    2023                          case CLP_PARAM_ACTION_NETLIB_PRIMAL:
    2024                          case CLP_PARAM_ACTION_NETLIB_TUNE: {
    2025                               // create fields for unitTest
    2026                               const char * fields[4];
    2027                               int nFields = 4;
    2028                               fields[0] = "fake main from unitTest";
    2029                               std::string mpsfield = "-dirSample=";
    2030                               mpsfield += dirSample.c_str();
    2031                               fields[1] = mpsfield.c_str();
    2032                               std::string netfield = "-dirNetlib=";
    2033                               netfield += dirNetlib.c_str();
    2034                               fields[2] = netfield.c_str();
    2035                               fields[3] = "-netlib";
    2036                               int algorithm;
    2037                               if (type == CLP_PARAM_ACTION_NETLIB_DUAL) {
    2038                                    std::cerr << "Doing netlib with dual algorithm" << std::endl;
    2039                                    algorithm = 0;
    2040                                    models[iModel].setMoreSpecialOptions(models[iModel].moreSpecialOptions()|32768);
    2041                               } else if (type == CLP_PARAM_ACTION_NETLIB_BARRIER) {
    2042                                    std::cerr << "Doing netlib with barrier algorithm" << std::endl;
    2043                                    algorithm = 2;
    2044                               } else if (type == CLP_PARAM_ACTION_NETLIB_EITHER) {
    2045                                    std::cerr << "Doing netlib with dual or primal algorithm" << std::endl;
    2046                                    algorithm = 3;
    2047                               } else if (type == CLP_PARAM_ACTION_NETLIB_TUNE) {
    2048                                    std::cerr << "Doing netlib with best algorithm!" << std::endl;
    2049                                    algorithm = 5;
    2050                                    // uncomment next to get active tuning
    2051                                    // algorithm=6;
    2052                               } else {
    2053                                    std::cerr << "Doing netlib with primal algorithm" << std::endl;
    2054                                    algorithm = 1;
    2055                               }
    2056                               //int specialOptions = models[iModel].specialOptions();
    2057                               models[iModel].setSpecialOptions(0);
    2058                               ClpSolve solveOptions;
    2059                               ClpSolve::PresolveType presolveType;
    2060                               if (preSolve)
    2061                                    presolveType = ClpSolve::presolveOn;
    2062                               else
    2063                                    presolveType = ClpSolve::presolveOff;
    2064                               solveOptions.setPresolveType(presolveType, 5);
    2065                               if (doSprint >= 0 || doIdiot >= 0) {
    2066                                    if (doSprint > 0) {
    2067                                         // sprint overrides idiot
    2068                                         solveOptions.setSpecialOption(1, 3, doSprint); // sprint
    2069                                    } else if (doIdiot > 0) {
    2070                                         solveOptions.setSpecialOption(1, 2, doIdiot); // idiot
    2071                                    } else {
    2072                                         if (doIdiot == 0) {
    2073                                              if (doSprint == 0)
    2074                                                   solveOptions.setSpecialOption(1, 4); // all slack
    2075                                              else
    2076                                                   solveOptions.setSpecialOption(1, 9); // all slack or sprint
    2077                                         } else {
    2078                                              if (doSprint == 0)
    2079                                                   solveOptions.setSpecialOption(1, 8); // all slack or idiot
    2080                                              else
    2081                                                   solveOptions.setSpecialOption(1, 7); // initiative
    2082                                         }
    2083                                    }
    2084                               }
    2085 #if FACTORIZATION_STATISTICS
    2086                               {
    2087                                 extern int loSizeX;
    2088                                 extern int hiSizeX;
    2089                                 for (int i=0;i<argc;i++) {
    2090                                   if (!strcmp(argv[i],"-losize")) {
    2091                                     int size=atoi(argv[i+1]);
    2092                                     if (size>0)
    2093                                       loSizeX=size;
    2094                                   }
    2095                                   if (!strcmp(argv[i],"-hisize")) {
    2096                                     int size=atoi(argv[i+1]);
    2097                                     if (size>loSizeX)
    2098                                       hiSizeX=size;
    2099                                   }
    2100                                 }
    2101                                 if (loSizeX!=-1||hiSizeX!=1000000)
    2102                                   printf("Solving problems %d<= and <%d\n",loSizeX,hiSizeX);
    2103                               }
    2104 #endif
    2105                               // for moment then back to models[iModel]
    2106 #ifndef ABC_INHERIT
    2107                               int specialOptions = models[iModel].specialOptions();
    2108                               mainTest(nFields, fields, algorithm, *thisModel,
    2109                                        solveOptions, specialOptions, doVector != 0);
    2110 #else
    2111                               //if (!processTune) {
    2112                               //specialOptions=0;
    2113                               //models->setSpecialOptions(models->specialOptions()&~65536);
    2114                               // }
    2115                               mainTest(nFields, fields, algorithm, *models,
    2116                                        solveOptions, 0, doVector != 0);
    2117 #endif
    2118                          }
    2119                          break;
    2120                          case CLP_PARAM_ACTION_UNITTEST: {
    2121                               // create fields for unitTest
    2122                               const char * fields[2];
    2123                               int nFields = 2;
    2124                               fields[0] = "fake main from unitTest";
    2125                               std::string dirfield = "-dirSample=";
    2126                               dirfield += dirSample.c_str();
    2127                               fields[1] = dirfield.c_str();
    2128                               int specialOptions = models[iModel].specialOptions();
    2129                               models[iModel].setSpecialOptions(0);
    2130                               int algorithm = -1;
    2131                               if (models[iModel].numberRows())
    2132                                    algorithm = 7;
    2133                               ClpSolve solveOptions;
    2134                               ClpSolve::PresolveType presolveType;
    2135                               if (preSolve)
    2136                                    presolveType = ClpSolve::presolveOn;
    2137                               else
    2138                                    presolveType = ClpSolve::presolveOff;
    2139                               solveOptions.setPresolveType(presolveType, 5);
    2140 #ifndef ABC_INHERIT
    2141                               mainTest(nFields, fields, algorithm, *thisModel,
    2142                                        solveOptions, specialOptions, doVector != 0);
    2143 #else
    2144                               mainTest(nFields, fields, algorithm, *models,
    2145                                        solveOptions, specialOptions, doVector != 0);
    2146 #endif
    2147                          }
    2148                          break;
    2149                          case CLP_PARAM_ACTION_FAKEBOUND:
    2150                               if (goodModels[iModel]) {
    2151                                    // get bound
    2152                                    double value = CoinReadGetDoubleField(argc, argv, &valid);
    2153                                    if (!valid) {
    2154                                         std::cout << "Setting " << parameters[iParam].name() <<
    2155                                                   " to DEBUG " << value << std::endl;
    2156                                         int iRow;
    2157                                         int numberRows = models[iModel].numberRows();
    2158                                         double * rowLower = models[iModel].rowLower();
    2159                                         double * rowUpper = models[iModel].rowUpper();
    2160                                         for (iRow = 0; iRow < numberRows; iRow++) {
    2161                                              // leave free ones for now
    2162                                              if (rowLower[iRow] > -1.0e20 || rowUpper[iRow] < 1.0e20) {
    2163                                                   rowLower[iRow] = CoinMax(rowLower[iRow], -value);
    2164                                                   rowUpper[iRow] = CoinMin(rowUpper[iRow], value);
    2165                                              }
    2166                                         }
    2167                                         int iColumn;
    2168                                         int numberColumns = models[iModel].numberColumns();
    2169                                         double * columnLower = models[iModel].columnLower();
    2170                                         double * columnUpper = models[iModel].columnUpper();
    2171                                         for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2172                                              // leave free ones for now
    2173                                              if (columnLower[iColumn] > -1.0e20 ||
    2174                                                        columnUpper[iColumn] < 1.0e20) {
    2175                                                   columnLower[iColumn] = CoinMax(columnLower[iColumn], -value);
    2176                                                   columnUpper[iColumn] = CoinMin(columnUpper[iColumn], value);
    2177                                              }
    2178                                         }
    2179                                    } else if (valid == 1) {
    2180                                         abort();
    2181                                    } else {
    2182                                         std::cout << "enter value for " << parameters[iParam].name() <<
    2183                                                   std::endl;
    2184                                    }
    2185                               }
    2186                               break;
    2187                          case CLP_PARAM_ACTION_REALLY_SCALE:
    2188                               if (goodModels[iModel]) {
    2189                                    ClpSimplex newModel(models[iModel],
    2190                                                        models[iModel].scalingFlag());
    2191                                    printf("model really really scaled\n");
    2192                                    models[iModel] = newModel;
    2193                               }
    2194                               break;
    2195                          case CLP_PARAM_ACTION_USERCLP:
    2196                               // Replace the sample code by whatever you want
    2197                               if (goodModels[iModel]) {
    2198                                    ClpSimplex * thisModel = &models[iModel];
    2199                                    printf("Dummy user code - model has %d rows and %d columns\n",
    2200                                           thisModel->numberRows(), thisModel->numberColumns());
    2201                               }
    2202                               break;
    2203                          case CLP_PARAM_ACTION_HELP:
    2204                               std::cout << "Coin LP version " << CLP_VERSION
    2205                                         << ", build " << __DATE__ << std::endl;
    2206                               std::cout << "Non default values:-" << std::endl;
    2207                               std::cout << "Perturbation " << models[0].perturbation() << " (default 100)"
    2208                                         << std::endl;
    2209                               CoinReadPrintit(
    2210                                    "Presolve being done with 5 passes\n\
    2211 Dual steepest edge steep/partial on matrix shape and factorization density\n\
    2212 Clpnnnn taken out of messages\n\
    2213 If Factorization frequency default then done on size of matrix\n\n\
    2214 (-)unitTest, (-)netlib or (-)netlibp will do standard tests\n\n\
    2215 You can switch to interactive mode at any time so\n\
    2216 clp watson.mps -scaling off -primalsimplex\nis the same as\n\
    2217 clp watson.mps -\nscaling off\nprimalsimplex"
    2218                               );
    2219                               break;
    2220                          case CLP_PARAM_ACTION_SOLUTION:
    2221                          case CLP_PARAM_ACTION_GMPL_SOLUTION:
    2222                               if (goodModels[iModel]) {
    2223                                    // get next field
    2224                                    field = CoinReadGetString(argc, argv);
    2225                                    bool append = false;
    2226                                    if (field == "append$") {
    2227                                      field = "$";
    2228                                      append = true;
    2229                                    }
    2230                                    if (field == "$") {
    2231                                         field = parameters[iParam].stringValue();
    2232                                    } else if (field == "EOL") {
    2233                                         parameters[iParam].printString();
    2234                                         break;
    2235                                    } else {
    2236                                         parameters[iParam].setStringValue(field);
    2237                                    }
    2238                                    std::string fileName;
    2239                                    FILE *fp = NULL;
    2240                                    if (field == "-" || field == "EOL" || field == "stdout") {
    2241                                         // stdout
    2242                                         fp = stdout;
    2243                                         fprintf(fp, "\n");
    2244                                    } else if (field == "stderr") {
    2245                                         // stderr
    2246                                         fp = stderr;
    2247                                         fprintf(fp, "\n");
    2248                                    } else {
    2249                                         if (field[0] == '/' || field[0] == '\\') {
    2250                                              fileName = field;
    2251                                         } else if (field[0] == '~') {
    2252                                              char * environVar = getenv("HOME");
    2253                                              if (environVar) {
    2254                                                   std::string home(environVar);
    2255                                                   field = field.erase(0, 1);
    2256                                                   fileName = home + field;
    2257                                              } else {
    2258                                                   fileName = field;
    2259                                              }
    2260                                         } else {
    2261                                              fileName = directory + field;
    2262                                         }
    2263                                         if (!append)
    2264                                           fp = fopen(fileName.c_str(), "w");
    2265                                         else
    2266                                           fp = fopen(fileName.c_str(), "a");
    2267                                    }
    2268                                    if (fp) {
    2269                                      // See if Glpk
    2270                                      if (type == CLP_PARAM_ACTION_GMPL_SOLUTION) {
    2271                                        int numberRows = models[iModel].getNumRows();
    2272                                        int numberColumns = models[iModel].getNumCols();
    2273                                        int numberGlpkRows=numberRows+1;
    2274 #ifdef COIN_HAS_GLPK
    2275                                        if (cbc_glp_prob) {
    2276                                          // from gmpl
    2277                                          numberGlpkRows=glp_get_num_rows(cbc_glp_prob);
    2278                                          if (numberGlpkRows!=numberRows)
    2279                                            printf("Mismatch - cbc %d rows, glpk %d\n",
    2280                                                   numberRows,numberGlpkRows);
    2281                                        }
    2282 #endif
    2283                                        fprintf(fp,"%d %d\n",numberGlpkRows,
    2284                                                numberColumns);
    2285                                        int iStat = models[iModel].status();
    2286                                        int iStat2 = GLP_UNDEF;
    2287                                        if (iStat == 0) {
    2288                                          // optimal
    2289                                          iStat2 = GLP_FEAS;
    2290                                        } else if (iStat == 1) {
    2291                                          // infeasible
    2292                                          iStat2 = GLP_NOFEAS;
    2293                                        } else if (iStat == 2) {
    2294                                          // unbounded
    2295                                          // leave as 1
    2296                                        } else if (iStat >= 3 && iStat <= 5) {
    2297                                          iStat2 = GLP_FEAS;
    2298                                        }
    2299                                        double objValue = models[iModel].getObjValue()
    2300                                          * models[iModel].getObjSense();
    2301                                        fprintf(fp,"%d 2 %g\n",iStat2,objValue);
    2302                                        if (numberGlpkRows > numberRows) {
    2303                                          // objective as row
    2304                                          fprintf(fp,"4 %g 1.0\n",objValue);
    2305                                        }
    2306                                        int lookup[6]=
    2307                                          {4,1,3,2,4,5};
    2308                                        const double * primalRowSolution =
    2309                                          models[iModel].primalRowSolution();
    2310                                        const double * dualRowSolution =
    2311                                          models[iModel].dualRowSolution();
    2312                                        for (int i=0;i<numberRows;i++) {
    2313                                          fprintf(fp,"%d %g %g\n",lookup[models[iModel].getRowStatus(i)],
    2314                                                  primalRowSolution[i],dualRowSolution[i]);
    2315                                        }
    2316                                        const double * primalColumnSolution =
    2317                                          models[iModel].primalColumnSolution();
    2318                                        const double * dualColumnSolution =
    2319                                          models[iModel].dualColumnSolution();
    2320                                        for (int i=0;i<numberColumns;i++) {
    2321                                          fprintf(fp,"%d %g %g\n",lookup[models[iModel].getColumnStatus(i)],
    2322                                                  primalColumnSolution[i],dualColumnSolution[i]);
    2323                                        }
    2324                                        fclose(fp);
    2325 #ifdef COIN_HAS_GLPK
    2326                                        if (cbc_glp_prob) {
    2327                                          glp_read_sol(cbc_glp_prob,fileName.c_str());
    2328                                          glp_mpl_postsolve(cbc_glp_tran,
    2329                                                            cbc_glp_prob,
    2330                                                            GLP_SOL);
    2331                                          // free up as much as possible
    2332                                          glp_free(cbc_glp_prob);
    2333                                          glp_mpl_free_wksp(cbc_glp_tran);
    2334                                          cbc_glp_prob = NULL;
    2335                                          cbc_glp_tran = NULL;
    2336                                          //gmp_free_mem();
    2337                                          /* check that no memory blocks are still allocated */
    2338                                          glp_free_env();
    2339                                        }
    2340 #endif
    2341                                        break;
    2342                                      }
    2343                                         // Write solution header (suggested by Luigi Poderico)
    2344                                         double objValue = models[iModel].getObjValue() * models[iModel].getObjSense();
    2345                                         int iStat = models[iModel].status();
    2346                                         if (iStat == 0) {
    2347                                              fprintf(fp, "optimal\n" );
    2348                                         } else if (iStat == 1) {
    2349                                              // infeasible
    2350                                              fprintf(fp, "infeasible\n" );
    2351                                         } else if (iStat == 2) {
    2352                                              // unbounded
    2353                                              fprintf(fp, "unbounded\n" );
    2354                                         } else if (iStat == 3) {
    2355                                              fprintf(fp, "stopped on iterations or time\n" );
    2356                                         } else if (iStat == 4) {
    2357                                              fprintf(fp, "stopped on difficulties\n" );
    2358                                         } else if (iStat == 5) {
    2359                                              fprintf(fp, "stopped on ctrl-c\n" );
    2360                                         } else {
    2361                                              fprintf(fp, "status unknown\n" );
    2362                                         }
    2363                                         char printFormat[50];
    2364                                         sprintf(printFormat,"Objective value %s\n",
    2365                                                 CLP_QUOTE(CLP_OUTPUT_FORMAT));
    2366                                         fprintf(fp, printFormat, objValue);
    2367                                         if (printMode==9) {
    2368                                           // just statistics
    2369                                           int numberRows = models[iModel].numberRows();
    2370                                           double * dualRowSolution = models[iModel].dualRowSolution();
    2371                                           double * primalRowSolution =
    2372                                             models[iModel].primalRowSolution();
    2373                                           double * rowLower = models[iModel].rowLower();
    2374                                           double * rowUpper = models[iModel].rowUpper();
    2375                                           double highestPrimal;
    2376                                           double lowestPrimal;
    2377                                           double highestDual;
    2378                                           double lowestDual;
    2379                                           double largestAway;
    2380                                           int numberAtLower;
    2381                                           int numberAtUpper;
    2382                                           int numberBetween;
    2383                                           highestPrimal=-COIN_DBL_MAX;
    2384                                           lowestPrimal=COIN_DBL_MAX;
    2385                                           highestDual=-COIN_DBL_MAX;
    2386                                           lowestDual=COIN_DBL_MAX;
    2387                                           largestAway=0.0;;
    2388                                           numberAtLower=0;
    2389                                           numberAtUpper=0;
    2390                                           numberBetween=0;
    2391                                           for (int iRow=0;iRow<numberRows;iRow++) {
    2392                                             double primal=primalRowSolution[iRow];
    2393                                             double lower=rowLower[iRow];
    2394                                             double upper=rowUpper[iRow];
    2395                                             double dual=dualRowSolution[iRow];
    2396                                             highestPrimal=CoinMax(highestPrimal,primal);
    2397                                             lowestPrimal=CoinMin(lowestPrimal,primal);
    2398                                             highestDual=CoinMax(highestDual,dual);
    2399                                             lowestDual=CoinMin(lowestDual,dual);
    2400                                             if (primal<lower+1.0e-6) {
    2401                                               numberAtLower++;
    2402                                             } else if (primal>upper-1.0e-6) {
    2403                                               numberAtUpper++;
    2404                                             } else {
    2405                                               numberBetween++;
    2406                                               largestAway=CoinMax(largestAway,
    2407                                                                   CoinMin(primal-lower,upper-primal));
    2408                                             }
    2409                                           }
    2410                                           printf("For rows %d at lower, %d between, %d at upper - lowest %g, highest %g most away %g - highest dual %g lowest %g\n",
    2411                                                  numberAtLower,numberBetween,numberAtUpper,
    2412                                                  lowestPrimal,highestPrimal,largestAway,
    2413                                                  lowestDual,highestDual);
    2414                                           int numberColumns = models[iModel].numberColumns();
    2415                                           double * dualColumnSolution = models[iModel].dualColumnSolution();
    2416                                           double * primalColumnSolution =
    2417                                             models[iModel].primalColumnSolution();
    2418                                           double * columnLower = models[iModel].columnLower();
    2419                                           double * columnUpper = models[iModel].columnUpper();
    2420                                           highestPrimal=-COIN_DBL_MAX;
    2421                                           lowestPrimal=COIN_DBL_MAX;
    2422                                           highestDual=-COIN_DBL_MAX;
    2423                                           lowestDual=COIN_DBL_MAX;
    2424                                           largestAway=0.0;;
    2425                                           numberAtLower=0;
    2426                                           numberAtUpper=0;
    2427                                           numberBetween=0;
    2428                                           for (int iColumn=0;iColumn<numberColumns;iColumn++) {
    2429                                             double primal=primalColumnSolution[iColumn];
    2430                                             double lower=columnLower[iColumn];
    2431                                             double upper=columnUpper[iColumn];
    2432                                             double dual=dualColumnSolution[iColumn];
    2433                                             highestPrimal=CoinMax(highestPrimal,primal);
    2434                                             lowestPrimal=CoinMin(lowestPrimal,primal);
    2435                                             highestDual=CoinMax(highestDual,dual);
    2436                                             lowestDual=CoinMin(lowestDual,dual);
    2437                                             if (primal<lower+1.0e-6) {
    2438                                               numberAtLower++;
    2439                                             } else if (primal>upper-1.0e-6) {
    2440                                               numberAtUpper++;
    2441                                             } else {
    2442                                               numberBetween++;
    2443                                               largestAway=CoinMax(largestAway,
    2444                                                                   CoinMin(primal-lower,upper-primal));
    2445                                             }
    2446                                           }
    2447                                           printf("For columns %d at lower, %d between, %d at upper - lowest %g, highest %g most away %g - highest dual %g lowest %g\n",
    2448                                                  numberAtLower,numberBetween,numberAtUpper,
    2449                                                  lowestPrimal,highestPrimal,largestAway,
    2450                                                  lowestDual,highestDual);
    2451                                           break;
    2452                                         }
    2453                                         // make fancy later on
    2454                                         int iRow;
    2455                                         int numberRows = models[iModel].numberRows();
    2456                                         int lengthName = models[iModel].lengthNames(); // 0 if no names
    2457                                         int lengthPrint = CoinMax(lengthName, 8);
    2458                                         // in general I don't want to pass around massive
    2459                                         // amounts of data but seems simpler here
    2460                                         std::vector<std::string> rowNames =
    2461                                              *(models[iModel].rowNames());
    2462                                         std::vector<std::string> columnNames =
    2463                                              *(models[iModel].columnNames());
    2464 
    2465                                         double * dualRowSolution = models[iModel].dualRowSolution();
    2466                                         double * primalRowSolution =
    2467                                              models[iModel].primalRowSolution();
    2468                                         double * rowLower = models[iModel].rowLower();
    2469                                         double * rowUpper = models[iModel].rowUpper();
    2470                                         double primalTolerance = models[iModel].primalTolerance();
    2471                                         bool doMask = (printMask != "" && lengthName);
    2472                                         int * maskStarts = NULL;
    2473                                         int maxMasks = 0;
    2474                                         char ** masks = NULL;
    2475                                         if (doMask) {
    2476                                              int nAst = 0;
    2477                                              const char * pMask2 = printMask.c_str();
    2478                                              char pMask[100];
    2479                                              size_t lengthMask = strlen(pMask2);
    2480                                              assert (lengthMask < 100);
    2481                                              if (*pMask2 == '"') {
    2482                                                   if (pMask2[lengthMask-1] != '"') {
    2483                                                        printf("mismatched \" in mask %s\n", pMask2);
    2484                                                        break;
    2485                                                   } else {
    2486                                                        strcpy(pMask, pMask2 + 1);
    2487                                                        *strchr(pMask, '"') = '\0';
    2488                                                   }
    2489                                              } else if (*pMask2 == '\'') {
    2490                                                   if (pMask2[lengthMask-1] != '\'') {
    2491                                                        printf("mismatched ' in mask %s\n", pMask2);
    2492                                                        break;
    2493                                                   } else {
    2494                                                        strcpy(pMask, pMask2 + 1);
    2495                                                        *strchr(pMask, '\'') = '\0';
    2496                                                   }
    2497                                              } else {
    2498                                                   strcpy(pMask, pMask2);
    2499                                              }
    2500                                              if (lengthMask > static_cast<size_t>(lengthName)) {
    2501                                                   printf("mask %s too long - skipping\n", pMask);
    2502                                                   break;
    2503                                              }
    2504                                              maxMasks = 1;
    2505                                              for (size_t iChar = 0; iChar < lengthMask; iChar++) {
    2506                                                   if (pMask[iChar] == '*') {
    2507                                                        nAst++;
    2508                                                        maxMasks *= (lengthName + 1);
    2509                                                   }
    2510                                              }
    2511                                              int nEntries = 1;
    2512                                              maskStarts = new int[lengthName+2];
    2513                                              masks = new char * [maxMasks];
    2514                                              char ** newMasks = new char * [maxMasks];
    2515                                              int i;
    2516                                              for (i = 0; i < maxMasks; i++) {
    2517                                                   masks[i] = new char[lengthName+1];
    2518                                                   newMasks[i] = new char[lengthName+1];
    2519                                              }
    2520                                              strcpy(masks[0], pMask);
    2521                                              for (int iAst = 0; iAst < nAst; iAst++) {
    2522                                                   int nOldEntries = nEntries;
    2523                                                   nEntries = 0;
    2524                                                   for (int iEntry = 0; iEntry < nOldEntries; iEntry++) {
    2525                                                        char * oldMask = masks[iEntry];
    2526                                                        char * ast = strchr(oldMask, '*');
    2527                                                        assert (ast);
    2528                                                        size_t length = strlen(oldMask) - 1;
    2529                                                        size_t nBefore = ast - oldMask;
    2530                                                        size_t nAfter = length - nBefore;
    2531                                                        // and add null
    2532                                                        nAfter++;
    2533                                                        for (size_t i = 0; i <= lengthName - length; i++) {
    2534                                                             char * maskOut = newMasks[nEntries];
    2535                                                             CoinMemcpyN(oldMask, static_cast<int>(nBefore), maskOut);
    2536                                                             for (size_t k = 0; k < i; k++)
    2537                                                                  maskOut[k+nBefore] = '?';
    2538                                                             CoinMemcpyN(ast + 1, static_cast<int>(nAfter), maskOut + nBefore + i);
    2539                                                             nEntries++;
    2540                                                             assert (nEntries <= maxMasks);
    2541                                                        }
    2542                                                   }
    2543                                                   char ** temp = masks;
    2544                                                   masks = newMasks;
    2545                                                   newMasks = temp;
    2546                                              }
    2547                                              // Now extend and sort
    2548                                              int * sort = new int[nEntries];
    2549                                              for (i = 0; i < nEntries; i++) {
    2550                                                   char * maskThis = masks[i];
    2551                                                   size_t length = strlen(maskThis);
    2552                                                   while (length > 0 && maskThis[length-1] == ' ')
    2553                                                        length--;
    2554                                                   maskThis[length] = '\0';
    2555                                                   sort[i] = static_cast<int>(length);
    2556                                              }
    2557                                              CoinSort_2(sort, sort + nEntries, masks);
    2558                                              int lastLength = -1;
    2559                                              for (i = 0; i < nEntries; i++) {
    2560                                                   int length = sort[i];
    2561                                                   while (length > lastLength)
    2562                                                        maskStarts[++lastLength] = i;
    2563                                              }
    2564                                              maskStarts[++lastLength] = nEntries;
    2565                                              delete [] sort;
    2566                                              for (i = 0; i < maxMasks; i++)
    2567                                                   delete [] newMasks[i];
    2568                                              delete [] newMasks;
    2569                                         }
    2570                                         if (printMode > 5) {
    2571                                           int numberColumns = models[iModel].numberColumns();
    2572                                           // column length unless rhs ranging
    2573                                           int number = numberColumns;
    2574                                           switch (printMode) {
    2575                                             // bound ranging
    2576                                             case 6:
    2577                                               fprintf(fp,"Bound ranging");
    2578                                               break;
    2579                                             // rhs ranging
    2580                                             case 7:
    2581                                               fprintf(fp,"Rhs ranging");
    2582                                               number = numberRows;
    2583                                               break;
    2584                                             // objective ranging
    2585                                             case 8:
    2586                                               fprintf(fp,"Objective ranging");
    2587                                               break;
    2588                                           }
    2589                                           if (lengthName)
    2590                                             fprintf(fp,",name");
    2591                                           fprintf(fp,",increase,variable,decrease,variable\n");
    2592                                           int * which = new int [ number];
    2593                                           if (printMode != 7) {
    2594                                             if (!doMask) {
    2595                                               for (int i = 0; i < number;i ++)
    2596                                                 which[i]=i;
    2597                                             } else {
    2598                                               int n = 0;
    2599                                               for (int i = 0; i < number;i ++) {
    2600                                                 if (maskMatches(maskStarts,masks,columnNames[i]))
    2601                                                   which[n++]=i;
    2602                                               }
    2603                                               if (n) {
    2604                                                 number=n;
    2605                                               } else {
    2606                                                 printf("No names match - doing all\n");
    2607                                                 for (int i = 0; i < number;i ++)
    2608                                                   which[i]=i;
    2609                                               }
    2610                                             }
    2611                                           } else {
    2612                                             if (!doMask) {
    2613                                               for (int i = 0; i < number;i ++)
    2614                                                 which[i]=i+numberColumns;
    2615                                             } else {
    2616                                               int n = 0;
    2617                                               for (int i = 0; i < number;i ++) {
    2618                                                 if (maskMatches(maskStarts,masks,rowNames[i]))
    2619                                                   which[n++]=i+numberColumns;
    2620                                               }
    2621                                               if (n) {
    2622                                                 number=n;
    2623                                               } else {
    2624                                                 printf("No names match - doing all\n");
    2625                                                 for (int i = 0; i < number;i ++)
    2626                                                   which[i]=i+numberColumns;
    2627                                               }
    2628                                             }
    2629                                           }
    2630                                           double * valueIncrease = new double [ number];
    2631                                           int * sequenceIncrease = new int [ number];
    2632                                           double * valueDecrease = new double [ number];
    2633                                           int * sequenceDecrease = new int [ number];
    2634                                           switch (printMode) {
    2635                                             // bound or rhs ranging
    2636                                             case 6:
    2637                                             case 7:
    2638                                               models[iModel].primalRanging(numberRows,
    2639                                                                            which, valueIncrease, sequenceIncrease,
    2640                                                                            valueDecrease, sequenceDecrease);
    2641                                               break;
    2642                                             // objective ranging
    2643                                             case 8:
    2644                                               models[iModel].dualRanging(number,
    2645                                                                            which, valueIncrease, sequenceIncrease,
    2646                                                                            valueDecrease, sequenceDecrease);
    2647                                               break;
    2648                                           }
    2649                                           for (int i = 0; i < number; i++) {
    2650                                             int iWhich = which[i];
    2651                                             fprintf(fp, "%d,", (iWhich<numberColumns) ? iWhich : iWhich-numberColumns);
    2652                                             if (lengthName) {
    2653                                               const char * name = (printMode==7) ? rowNames[iWhich-numberColumns].c_str() : columnNames[iWhich].c_str();
    2654                                               fprintf(fp,"%s,",name);
    2655                                             }
    2656                                             if (valueIncrease[i]<1.0e30) {
    2657                                               fprintf(fp, "%.10g,", valueIncrease[i]);
    2658                                               int outSequence = sequenceIncrease[i];
    2659                                               if (outSequence<numberColumns) {
    2660                                                 if (lengthName)
    2661                                                   fprintf(fp,"%s,",columnNames[outSequence].c_str());
    2662                                                 else
    2663                                                   fprintf(fp,"C%7.7d,",outSequence);
    2664                                               } else {
    2665                                                 outSequence -= numberColumns;
    2666                                                 if (lengthName)
    2667                                                   fprintf(fp,"%s,",rowNames[outSequence].c_str());
    2668                                                 else
    2669                                                   fprintf(fp,"R%7.7d,",outSequence);
    2670                                               }
    2671                                             } else {
    2672                                               fprintf(fp,"1.0e100,,");
    2673                                             }
    2674                                             if (valueDecrease[i]<1.0e30) {
    2675                                               fprintf(fp, "%.10g,", valueDecrease[i]);
    2676                                               int outSequence = sequenceDecrease[i];
    2677                                               if (outSequence<numberColumns) {
    2678                                                 if (lengthName)
    2679                                                   fprintf(fp,"%s",columnNames[outSequence].c_str());
    2680                                                 else
    2681                                                   fprintf(fp,"C%7.7d",outSequence);
    2682                                               } else {
    2683                                                 outSequence -= numberColumns;
    2684                                                 if (lengthName)
    2685                                                   fprintf(fp,"%s",rowNames[outSequence].c_str());
    2686                                                 else
    2687                                                   fprintf(fp,"R%7.7d",outSequence);
    2688                                               }
    2689                                             } else {
    2690                                               fprintf(fp,"1.0e100,");
    2691                                             }
    2692                                             fprintf(fp,"\n");
    2693                                           }
    2694                                           if (fp != stdout)
    2695                                             fclose(fp);
    2696                                           delete [] which;
    2697                                           delete [] valueIncrease;
    2698                                           delete [] sequenceIncrease;
    2699                                           delete [] valueDecrease;
    2700                                           delete [] sequenceDecrease;
    2701                                           if (masks) {
    2702                                             delete [] maskStarts;
    2703                                             for (int i = 0; i < maxMasks; i++)
    2704                                               delete [] masks[i];
    2705                                             delete [] masks;
    2706                                           }
    2707                                           break;
    2708                                         }
    2709                                         sprintf(printFormat," %s         %s\n",
    2710                                                 CLP_QUOTE(CLP_OUTPUT_FORMAT),
    2711                                                 CLP_QUOTE(CLP_OUTPUT_FORMAT));
    2712                                         if (printMode > 2) {
    2713                                              for (iRow = 0; iRow < numberRows; iRow++) {
    2714                                                   int type = printMode - 3;
    2715                                                   if (primalRowSolution[iRow] > rowUpper[iRow] + primalTolerance ||
    2716                                                             primalRowSolution[iRow] < rowLower[iRow] - primalTolerance) {
    2717                                                        fprintf(fp, "** ");
    2718                                                        type = 2;
    2719                                                   } else if (fabs(primalRowSolution[iRow]) > 1.0e-8) {
    2720                                                        type = 1;
    2721                                                   } else if (numberRows < 50) {
    2722                                                        type = 3;
    2723                                                   }
    2724                                                   if (doMask && !maskMatches(maskStarts, masks, rowNames[iRow]))
    2725                                                        type = 0;
    2726                                                   if (type) {
    2727                                                        fprintf(fp, "%7d ", iRow);
    2728                                                        if (lengthName) {
    2729                                                             const char * name = rowNames[iRow].c_str();
    2730                                                             size_t n = strlen(name);
    2731                                                             size_t i;
    2732                                                             for (i = 0; i < n; i++)
    2733                                                                  fprintf(fp, "%c", name[i]);
    2734                                                             for (; i < static_cast<size_t>(lengthPrint); i++)
    2735                                                                  fprintf(fp, " ");
    2736                                                        }
    2737                                                        fprintf(fp, printFormat, primalRowSolution[iRow],
    2738                                                                dualRowSolution[iRow]);
    2739                                                   }
    2740                                              }
    2741                                         }
    2742                                         int iColumn;
    2743                                         int numberColumns = models[iModel].numberColumns();
    2744                                         double * dualColumnSolution =
    2745                                              models[iModel].dualColumnSolution();
    2746                                         double * primalColumnSolution =
    2747                                              models[iModel].primalColumnSolution();
    2748                                         double * columnLower = models[iModel].columnLower();
    2749                                         double * columnUpper = models[iModel].columnUpper();
    2750                                         for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2751                                              int type = (printMode > 3) ? 1 : 0;
    2752                                              if (primalColumnSolution[iColumn] > columnUpper[iColumn] + primalTolerance ||
    2753                                                        primalColumnSolution[iColumn] < columnLower[iColumn] - primalTolerance) {
    2754                                                   fprintf(fp, "** ");
    2755                                                   type = 2;
    2756                                              } else if (fabs(primalColumnSolution[iColumn]) > 1.0e-8) {
    2757                                                   type = 1;
    2758                                              } else if (numberColumns < 50) {
    2759                                                   type = 3;
    2760                                              }
    2761                                              if (doMask && !maskMatches(maskStarts, masks,
    2762                                                                         columnNames[iColumn]))
    2763                                                   type = 0;
    2764                                              if (type) {
    2765                                                   fprintf(fp, "%7d ", iColumn);
    2766                                                   if (lengthName) {
    2767                                                        const char * name = columnNames[iColumn].c_str();
    2768                                                        size_t n = strlen(name);
    2769                                                        size_t i;
    2770                                                        for (i = 0; i < n; i++)
    2771                                                             fprintf(fp, "%c", name[i]);
    2772                                                        for (; i < static_cast<size_t>(lengthPrint); i++)
    2773                                                             fprintf(fp, " ");
    2774                                                   }
    2775                                                   fprintf(fp, printFormat,
    2776                                                           primalColumnSolution[iColumn],
    2777                                                           dualColumnSolution[iColumn]);
    2778                                              }
    2779                                         }
    2780                                         if (fp != stdout)
    2781                                              fclose(fp);
    2782                                         if (masks) {
    2783                                              delete [] maskStarts;
    2784                                              for (int i = 0; i < maxMasks; i++)
    2785                                                   delete [] masks[i];
    2786                                              delete [] masks;
    2787                                         }
    2788                                    } else {
    2789                                         std::cout << "Unable to open file " << fileName << std::endl;
    2790                                    }
    2791                               } else {
    2792                                    std::cout << "** Current model not valid" << std::endl;
    2793 
    2794                               }
    2795 
    2796                               break;
    2797                          case CLP_PARAM_ACTION_SAVESOL:
    2798                               if (goodModels[iModel]) {
    2799                                    // get next field
    2800                                    field = CoinReadGetString(argc, argv);
    2801                                    if (field == "$") {
    2802                                         field = parameters[iParam].stringValue();
    2803                                    } else if (field == "EOL") {
    2804                                         parameters[iParam].printString();
    2805                                         break;
    2806                                    } else {
    2807                                         parameters[iParam].setStringValue(field);
    2808                                    }
    2809                                    std::string fileName;
    2810                                    if (field[0] == '/' || field[0] == '\\') {
    2811                                         fileName = field;
    2812                                    } else if (field[0] == '~') {
    2813                                         char * environVar = getenv("HOME");
    2814                                         if (environVar) {
    2815                                              std::string home(environVar);
    2816                                              field = field.erase(0, 1);
    2817                                              fileName = home + field;
    2818                                         } else {
    2819                                              fileName = field;
    2820                                         }
    2821                                    } else {
    2822                                         fileName = directory + field;
    2823                                    }
    2824                                    saveSolution(models + iModel, fileName);
    2825                               } else {
    2826                                    std::cout << "** Current model not valid" << std::endl;
    2827 
    2828                               }
    2829                               break;
    2830                          case CLP_PARAM_ACTION_ENVIRONMENT:
    2831                               CbcOrClpEnvironmentIndex = 0;
    2832                               break;
    2833                          default:
    2834                               abort();
    2835                          }
    2836                     }
    2837                } else if (!numberMatches) {
    2838                     std::cout << "No match for " << field << " - ? for list of commands"
    2839                               << std::endl;
    2840                } else if (numberMatches == 1) {
    2841                     if (!numberQuery) {
    2842                          std::cout << "Short match for " << field << " - completion: ";
    2843                          std::cout << parameters[firstMatch].matchName() << std::endl;
    2844                     } else if (numberQuery) {
    2845                          std::cout << parameters[firstMatch].matchName() << " : ";
    2846                          std::cout << parameters[firstMatch].shortHelp() << std::endl;
    2847                          if (numberQuery >= 2)
    2848                               parameters[firstMatch].printLongHelp();
    2849                     }
    2850                } else {
    2851                     if (!numberQuery)
    2852                          std::cout << "Multiple matches for " << field << " - possible completions:"
    2853                                    << std::endl;
    2854                     else
    2855                          std::cout << "Completions of " << field << ":" << std::endl;
    2856                     for ( iParam = 0; iParam < numberParameters; iParam++ ) {
    2857                          int match = parameters[iParam].matches(field);
    2858                          if (match && parameters[iParam].displayThis()) {
    2859                               std::cout << parameters[iParam].matchName();
    2860                               if (numberQuery >= 2)
    2861                                    std::cout << " : " << parameters[iParam].shortHelp();
    2862                               std::cout << std::endl;
    2863                          }
    2864                     }
    2865                }
    2866           }
    2867           delete [] models;
    2868           delete [] goodModels;
    2869      }
    2870 #ifdef COIN_HAS_GLPK
    2871      if (cbc_glp_prob) {
    2872        // free up as much as possible
    2873        glp_free(cbc_glp_prob);
    2874        glp_mpl_free_wksp(cbc_glp_tran);
    2875        glp_free_env();
    2876        cbc_glp_prob = NULL;
    2877        cbc_glp_tran = NULL;
    2878      }
    2879 #endif
    2880      // By now all memory should be freed
    2881 #ifdef DMALLOC
    2882      dmalloc_log_unfreed();
    2883      dmalloc_shutdown();
    2884 #endif
    2885      return 0;
    2886 }
    2887 static void breakdown(const char * name, int numberLook, const double * region)
    2888 {
    2889      double range[] = {
    2890           -COIN_DBL_MAX,
    2891           -1.0e15, -1.0e11, -1.0e8, -1.0e5, -1.0e4, -1.0e3, -1.0e2, -1.0e1,
    2892           -1.0,
    2893           -1.0e-1, -1.0e-2, -1.0e-3, -1.0e-4, -1.0e-5, -1.0e-8, -1.0e-11, -1.0e-15,
    2894           0.0,
    2895           1.0e-15, 1.0e-11, 1.0e-8, 1.0e-5, 1.0e-4, 1.0e-3, 1.0e-2, 1.0e-1,
    2896           1.0,
    2897           1.0e1, 1.0e2, 1.0e3, 1.0e4, 1.0e5, 1.0e8, 1.0e11, 1.0e15,
    2898           COIN_DBL_MAX
    2899      };
    2900      int nRanges = static_cast<int> (sizeof(range) / sizeof(double));
    2901      int * number = new int[nRanges];
    2902      memset(number, 0, nRanges * sizeof(int));
    2903      int * numberExact = new int[nRanges];
    2904      memset(numberExact, 0, nRanges * sizeof(int));
    2905      int i;
    2906      for ( i = 0; i < numberLook; i++) {
    2907           double value = region[i];
    2908           for (int j = 0; j < nRanges; j++) {
    2909                if (value == range[j]) {
    2910                     numberExact[j]++;
    2911                     break;
    2912                } else if (value < range[j]) {
    2913                     number[j]++;
    2914                     break;
    2915                }
    2916           }
    2917      }
    2918      printf("\n%s has %d entries\n", name, numberLook);
    2919      for (i = 0; i < nRanges; i++) {
    2920           if (number[i])
    2921                printf("%d between %g and %g", number[i], range[i-1], range[i]);
    2922           if (numberExact[i]) {
    2923                if (number[i])
    2924                     printf(", ");
    2925                printf("%d exactly at %g", numberExact[i], range[i]);
    2926           }
    2927           if (number[i] + numberExact[i])
    2928                printf("\n");
    2929      }
    2930      delete [] number;
    2931      delete [] numberExact;
    2932 }
    2933 void sortOnOther(int * column,
    2934                  const CoinBigIndex * rowStart,
    2935                  int * order,
    2936                  int * other,
    2937                  int nRow,
    2938                  int nInRow,
    2939                  int where)
    2940 {
    2941      if (nRow < 2 || where >= nInRow)
    2942           return;
    2943      // do initial sort
    2944      int kRow;
    2945      int iRow;
    2946      for ( kRow = 0; kRow < nRow; kRow++) {
    2947           iRow = order[kRow];
    2948           other[kRow] = column[rowStart[iRow] + where];
    2949      }
    2950      CoinSort_2(other, other + nRow, order);
    2951      int first = 0;
    2952      iRow = order[0];
    2953      int firstC = column[rowStart[iRow] + where];
    2954      kRow = 1;
    2955      while (kRow < nRow) {
    2956           int lastC = 9999999;;
    2957           for (; kRow < nRow + 1; kRow++) {
    2958                if (kRow < nRow) {
    2959                     iRow = order[kRow];
    2960                     lastC = column[rowStart[iRow] + where];
    2961                } else {
    2962                     lastC = 9999999;
    2963                }
    2964                if (lastC > firstC)
    2965                     break;
    2966           }
    2967           // sort
    2968           sortOnOther(column, rowStart, order + first, other, kRow - first,
    2969                       nInRow, where + 1);
    2970           firstC = lastC;
    2971           first = kRow;
    2972      }
    2973 }
    2974 static void statistics(ClpSimplex * originalModel, ClpSimplex * model)
    2975 {
    2976      int numberColumns = originalModel->numberColumns();
    2977      const char * integerInformation  = originalModel->integerInformation();
    2978      const double * columnLower = originalModel->columnLower();
    2979      const double * columnUpper = originalModel->columnUpper();
    2980      int numberIntegers = 0;
    2981      int numberBinary = 0;
    2982      int iRow, iColumn;
    2983      if (integerInformation) {
    2984           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    2985                if (integerInformation[iColumn]) {
    2986                     if (columnUpper[iColumn] > columnLower[iColumn]) {
    2987                          numberIntegers++;
    2988                          if (columnLower[iColumn] == 0.0 && columnUpper[iColumn] == 1)
    2989                               numberBinary++;
    2990                     }
    2991                }
    2992           }
    2993           printf("Original problem has %d integers (%d of which binary)\n",
    2994                  numberIntegers,numberBinary);
    2995      }
    2996      numberColumns = model->numberColumns();
    2997      int numberRows = model->numberRows();
    2998      columnLower = model->columnLower();
    2999      columnUpper = model->columnUpper();
    3000      const double * rowLower = model->rowLower();
    3001      const double * rowUpper = model->rowUpper();
    3002      const double * objective = model->objective();
    3003      if (model->integerInformation()) {
    3004        const char * integerInformation  = model->integerInformation();
    3005        int numberIntegers = 0;
    3006        int numberBinary = 0;
    3007        double * obj = new double [numberColumns];
    3008        int * which = new int [numberColumns];
    3009        for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
    3010          if (columnUpper[iColumn] > columnLower[iColumn]) {
    3011            if (integerInformation[iColumn]) {
    3012              numberIntegers++;
    3013              if (columnLower[iColumn] == 0.0 && columnUpper[iColumn] == 1)
    3014                numberBinary++;
    3015            }
    3016          }
    3017        }
    3018        if(numberColumns != originalModel->numberColumns())
    3019          printf("Presolved problem has %d integers (%d of which binary)\n",
    3020                 numberIntegers,numberBinary);
    3021        for (int ifInt=0;ifInt<2;ifInt++) {
    3022          for (int ifAbs=0;ifAbs<2;ifAbs++) {
    3023            int numberSort=0;
    3024            int numberZero=0;
    3025            int numberDifferentObj=0;
    3026            for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
    3027              if (columnUpper[iColumn] > columnLower[iColumn]) {
    3028                if (!ifInt||integerInformation[iColumn]) {
    3029                  obj[numberSort]=(ifAbs) ? fabs(objective[iColumn]) :
    3030                    objective[iColumn];
    3031                  which[numberSort++]=iColumn;
    3032                  if (!objective[iColumn])
    3033                    numberZero++;
    3034                }
    3035              }
    3036            }
    3037            CoinSort_2(obj,obj+numberSort,which);
    3038            double last=obj[0];
    3039            for (int jColumn = 1; jColumn < numberSort; jColumn++) {
    3040              if (fabs(obj[jColumn]-last)>1.0e-12) {
    3041                numberDifferentObj++;
    3042                last=obj[jColumn];
    3043              }
    3044            }
    3045            numberDifferentObj++;
    3046            printf("==== ");
    3047            if (ifInt)
    3048              printf("for integers ");
    3049            if (!ifAbs)
    3050              printf("%d zero objective ",numberZero);
    3051            else
    3052              printf("absolute objective values ");
    3053            printf("%d different\n",numberDifferentObj);
    3054            bool saveModel=false;
    3055            int target=model->logLevel();
    3056            if (target>10000) {
    3057              if (ifInt&&!ifAbs)
    3058                saveModel=true;
    3059              target-=10000;
    3060            }
    3061 
    3062            if (target<=100)
    3063              target=12;
    3064            else
    3065              target-=100;
    3066            if (numberDifferentObj<target) {
    3067              int iLast=0;
    3068              double last=obj[0];
    3069              for (int jColumn = 1; jColumn < numberSort; jColumn++) {
    3070                if (fabs(obj[jColumn]-last)>1.0e-12) {
    3071                  printf("%d variables have objective of %g\n",
    3072                         jColumn-iLast,last);
    3073                  iLast=jColumn;
    3074                  last=obj[jColumn];
    3075                }
    3076              }
    3077              printf("%d variables have objective of %g\n",
    3078                     numberSort-iLast,last);
    3079              if (saveModel) {
    3080                int spaceNeeded=numberSort+numberDifferentObj;
    3081                int * columnAdd = new int[spaceNeeded+numberDifferentObj+1];
    3082                double * elementAdd = new double[spaceNeeded];
    3083                int * rowAdd = new int[2*numberDifferentObj+1];
    3084                int * newIsInteger = rowAdd+numberDifferentObj+1;
    3085                double * objectiveNew = new double[3*numberDifferentObj];
    3086                double * lowerNew = objectiveNew+numberDifferentObj;
    3087                double * upperNew = lowerNew+numberDifferentObj;
    3088                memset(columnAdd+spaceNeeded,0,
    3089                       (numberDifferentObj+1)*sizeof(int));
    3090                ClpSimplex tempModel=*model;
    3091                int iLast=0;
    3092                double last=obj[0];
    3093                numberDifferentObj=0;
    3094                int numberElements=0;
    3095                rowAdd[0]=0;
    3096                double * objective = tempModel.objective();
    3097                for (int jColumn = 1; jColumn < numberSort+1; jColumn++) {
    3098                  if (jColumn==numberSort||fabs(obj[jColumn]-last)>1.0e-12) {
    3099                    // not if just one
    3100                    if (jColumn-iLast>1) {
    3101                      bool allInteger=integerInformation!=NULL;
    3102                      int iColumn=which[iLast];
    3103                      objectiveNew[numberDifferentObj]=objective[iColumn];
    3104                      double lower=0.0;
    3105                      double upper=0.0;
    3106                      for (int kColumn=iLast;kColumn<jColumn;kColumn++) {
    3107                        iColumn=which[kColumn];
    3108                        objective[iColumn]=0.0;
    3109                        double lowerValue=columnLower[iColumn];
    3110                        double upperValue=columnUpper[iColumn];
    3111                        double elementValue=-1.0;
    3112                        if (objectiveNew[numberDifferentObj]*objective[iColumn]<0.0) {
    3113                          lowerValue=-columnUpper[iColumn];
    3114                          upperValue=-columnLower[iColumn];
    3115                          elementValue=1.0;
    3116                        }
    3117                        columnAdd[numberElements]=iColumn;
    3118                        elementAdd[numberElements++]=elementValue;
    3119                        if (integerInformation&&!integerInformation[iColumn])
    3120                          allInteger=false;
    3121                        if (lower!=-COIN_DBL_MAX) {
    3122                          if (lowerValue!=-COIN_DBL_MAX)
    3123                            lower += lowerValue;
    3124                          else
    3125                            lower=-COIN_DBL_MAX;
    3126                        }
    3127                        if (upper!=COIN_DBL_MAX) {
    3128                          if (upperValue!=COIN_DBL_MAX)
    3129                            upper += upperValue;
    3130                          else
    3131                            upper=COIN_DBL_MAX;
    3132                        }
    3133                      }
    3134                      columnAdd[numberElements]=numberColumns+numberDifferentObj;
    3135                      elementAdd[numberElements++]=1.0;
    3136                      newIsInteger[numberDifferentObj]= (allInteger) ? 1 : 0;
    3137                      lowerNew[numberDifferentObj]=lower;
    3138                      upperNew[numberDifferentObj]=upper;
    3139                      numberDifferentObj++;
    3140                      rowAdd[numberDifferentObj]=numberElements;
    3141                    }
    3142                    iLast=jColumn;
    3143                    last=obj[jColumn];
    3144                  }
    3145                }
    3146                // add columns
    3147                tempModel.addColumns(numberDifferentObj, lowerNew, upperNew,
    3148                                     objectiveNew,
    3149                                     columnAdd+spaceNeeded, NULL, NULL);
    3150                // add constraints and make integer if all integer in group
    3151                for (int iObj=0; iObj < numberDifferentObj; iObj++) {
    3152                  lowerNew[iObj]=0.0;
    3153                  upperNew[iObj]=0.0;
    3154                  if (newIsInteger[iObj])
    3155                    tempModel.setInteger(numberColumns+iObj);
    3156                }
    3157                tempModel.addRows(numberDifferentObj, lowerNew, upperNew,
    3158                                  rowAdd,columnAdd,elementAdd);
    3159                delete [] columnAdd;
    3160                delete [] elementAdd;
    3161                delete [] rowAdd;
    3162                delete [] objectiveNew;
    3163                // save
    3164                std::string tempName = model->problemName();
    3165                if (ifInt)
    3166                  tempName += "_int";
    3167                if (ifAbs)
    3168                  tempName += "_abs";
    3169                tempName += ".mps";
    3170                tempModel.writeMps(tempName.c_str());
    3171              }
    3172            }
    3173          }
    3174        }
    3175        delete [] which;
    3176        delete [] obj;
    3177        printf("===== end objective counts\n");
    3178      }
    3179      CoinPackedMatrix * matrix = model->matrix();
    3180      CoinBigIndex numberElements = matrix->getNumElements();
    3181      const int * columnLength = matrix->getVectorLengths();
    3182      //const CoinBigIndex * columnStart = matrix->getVectorStarts();
    3183      const double * elementByColumn = matrix->getElements();
    3184      int * number = new int[numberRows+1];
    3185      memset(number, 0, (numberRows + 1)*sizeof(int));
    3186      int numberObjSingletons = 0;
    3187      /* cType
    3188         0 0/inf, 1 0/up, 2 lo/inf, 3 lo/up, 4 free, 5 fix, 6 -inf/0, 7 -inf/up,
    3189         8 0/1
    3190      */
    3191      int cType[9];
    3192      std::string cName[] = {"0.0->inf,", "0.0->up,", "lo->inf,", "lo->up,", "free,", "fixed,", "-inf->0.0,",
    3193                             "-inf->up,", "0.0->1.0"
    3194                            };
    3195      int nObjective = 0;
    3196      memset(cType, 0, sizeof(cType));
    3197      for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    3198           int length = columnLength[iColumn];
    3199           if (length == 1 && objective[iColumn])
    3200                numberObjSingletons++;
    3201           number[length]++;
    3202           if (objective[iColumn])
    3203                nObjective++;
    3204           if (columnLower[iColumn] > -1.0e20) {
    3205                if (columnLower[iColumn] == 0.0) {
    3206                     if (columnUpper[iColumn] > 1.0e20)
    3207                          cType[0]++;
    3208                     else if (columnUpper[iColumn] == 1.0)
    3209                          cType[8]++;
    3210                     else if (columnUpper[iColumn] == 0.0)
    3211                          cType[5]++;
    3212                     else
    3213                          cType[1]++;
    3214                } else {
    3215                     if (columnUpper[iColumn] > 1.0e20)
    3216                          cType[2]++;
    3217                     else if (columnUpper[iColumn] == columnLower[iColumn])
    3218                          cType[5]++;
    3219                     else
    3220                          cType[3]++;
    3221                }
    3222           } else {
    3223                if (columnUpper[iColumn] > 1.0e20)
    3224                     cType[4]++;
    3225                else if (columnUpper[iColumn] == 0.0)
    3226                     cType[6]++;
    3227                else
    3228                     cType[7]++;
    3229           }
    3230      }
    3231      /* rType
    3232         0 E 0, 1 E 1, 2 E -1, 3 E other, 4 G 0, 5 G 1, 6 G other,
    3233         7 L 0,  8 L 1, 9 L other, 10 Range 0/1, 11 Range other, 12 free
    3234      */
    3235      int rType[13];
    3236      std::string rName[] = {"E 0.0,", "E 1.0,", "E -1.0,", "E other,", "G 0.0,", "G 1.0,", "G other,",
    3237                             "L 0.0,", "L 1.0,", "L other,", "Range 0.0->1.0,", "Range other,", "Free"
    3238                            };
    3239      memset(rType, 0, sizeof(rType));
    3240      for (iRow = 0; iRow < numberRows; iRow++) {
    3241           if (rowLower[iRow] > -1.0e20) {
    3242                if (rowLower[iRow] == 0.0) {
    3243                     if (rowUpper[iRow] > 1.0e20)
    3244                          rType[4]++;
    3245                     else if (rowUpper[iRow] == 1.0)
    3246                          rType[10]++;
    3247                     else if (rowUpper[iRow] == 0.0)
    3248                          rType[0]++;
    3249                     else
    3250                          rType[11]++;
    3251                } else if (rowLower[iRow] == 1.0) {
    3252                     if (rowUpper[iRow] > 1.0e20)
    3253                          rType[5]++;
    3254                     else if (rowUpper[iRow] == rowLower[iRow])
    3255                          rType[1]++;
    3256                     else
    3257                          rType[11]++;
    3258                } else if (rowLower[iRow] == -1.0) {
    3259                     if (rowUpper[iRow] > 1.0e20)
    3260                          rType[6]++;
    3261                     else if (rowUpper[iRow] == rowLower[iRow])
    3262                          rType[2]++;
    3263                     else
    3264                          rType[11]++;
    3265                } else {
    3266                     if (rowUpper[iRow] > 1.0e20)
    3267                          rType[6]++;
    3268                     else if (rowUpper[iRow] == rowLower[iRow])
    3269                          rType[3]++;
    3270                     else
    3271                          rType[11]++;
    3272                }
    3273           } else {
    3274                if (rowUpper[iRow] > 1.0e20)
    3275                     rType[12]++;
    3276                else if (rowUpper[iRow] == 0.0)
    3277                     rType[7]++;
    3278                else if (rowUpper[iRow] == 1.0)
    3279                     rType[8]++;
    3280                else
    3281                     rType[9]++;
    3282           }
    3283      }
    3284      // Basic statistics
    3285      printf("\n\nProblem has %d rows, %d columns (%d with objective) and %d elements\n",
    3286             numberRows, numberColumns, nObjective, numberElements);
    3287      if (number[0] + number[1]) {
    3288           printf("There are ");
    3289           if (numberObjSingletons)
    3290                printf("%d singletons with objective ", numberObjSingletons);
    3291           int numberNoObj = number[1] - numberObjSingletons;
    3292           if (numberNoObj)
    3293                printf("%d singletons with no objective ", numberNoObj);
    3294           if (number[0])
    3295                printf("** %d columns have no entries", number[0]);
    3296           printf("\n");
    3297      }
    3298      printf("Column breakdown:\n");
    3299      int k;
    3300      for (k = 0; k < static_cast<int> (sizeof(cType) / sizeof(int)); k++) {
    3301           printf("%d of type %s ", cType[k], cName[k].c_str());
    3302           if (((k + 1) % 3) == 0)
    3303                printf("\n");
    3304      }
    3305      if ((k % 3) != 0)
    3306           printf("\n");
    3307      printf("Row breakdown:\n");
    3308      for (k = 0; k < static_cast<int> (sizeof(rType) / sizeof(int)); k++) {
    3309           printf("%d of type %s ", rType[k], rName[k].c_str());
    3310           if (((k + 1) % 3) == 0)
    3311                printf("\n");
    3312      }
    3313      if ((k % 3) != 0)
    3314           printf("\n");
    3315      //#define SYM
    3316 #ifndef SYM
    3317      if (model->logLevel() < 2)
    3318           return ;
    3319 #endif
    3320      int kMax = model->logLevel() > 3 ? 1000000 : 10;
    3321      k = 0;
    3322      for (iRow = 1; iRow <= numberRows; iRow++) {
    3323           if (number[iRow]) {
    3324                k++;
    3325                printf("%d columns have %d entries\n", number[iRow], iRow);
    3326                if (k == kMax)
    3327                     break;
    3328           }
    3329      }
    3330      if (k < numberRows) {
    3331           int kk = k;
    3332           k = 0;
    3333           for (iRow = numberRows; iRow >= 1; iRow--) {
    3334                if (number[iRow]) {
    3335                     k++;
    3336                     if (k == kMax)
    3337                          break;
    3338                }
    3339           }
    3340           if (k > kk) {
    3341                printf("\n    .........\n\n");
    3342                iRow = k;
    3343                k = 0;
    3344                for (; iRow < numberRows; iRow++) {
    3345                     if (number[iRow]) {
    3346                          k++;
    3347                          printf("%d columns have %d entries\n", number[iRow], iRow);
    3348                          if (k == kMax)
    3349                               break;
    3350                     }
    3351                }
    3352           }
    3353      }
    3354      delete [] number;
    3355      printf("\n\n");
    3356      if (model->logLevel() == 63
    3357 #ifdef SYM
    3358                || true
    3359 #endif
    3360         ) {
    3361           // get column copy
    3362           CoinPackedMatrix columnCopy = *matrix;
    3363           const int * columnLength = columnCopy.getVectorLengths();
    3364           number = new int[numberRows+1];
    3365           memset(number, 0, (numberRows + 1)*sizeof(int));
    3366           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    3367                int length = columnLength[iColumn];
    3368                number[length]++;
    3369           }
    3370           k = 0;
    3371           for (iRow = 1; iRow <= numberRows; iRow++) {
    3372                if (number[iRow]) {
    3373                     k++;
    3374                }
    3375           }
    3376           int * row = columnCopy.getMutableIndices();
    3377           const CoinBigIndex * columnStart = columnCopy.getVectorStarts();
    3378           double * element = columnCopy.getMutableElements();
    3379           int * order = new int[numberColumns];
    3380           int * other = new int[numberColumns];
    3381           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    3382                int length = columnLength[iColumn];
    3383                order[iColumn] = iColumn;
    3384                other[iColumn] = length;
    3385                CoinBigIndex start = columnStart[iColumn];
    3386                CoinSort_2(row + start, row + start + length, element + start);
    3387           }
    3388           CoinSort_2(other, other + numberColumns, order);
    3389           int jColumn = number[0] + number[1];
    3390           for (iRow = 2; iRow <= numberRows; iRow++) {
    3391                if (number[iRow]) {
    3392                     printf("XX %d columns have %d entries\n", number[iRow], iRow);
    3393                     int kColumn = jColumn + number[iRow];
    3394                     sortOnOther(row, columnStart,
    3395                                 order + jColumn, other, number[iRow], iRow, 0);
    3396                     // Now print etc
    3397                     if (iRow < 500000) {
    3398                          for (int lColumn = jColumn; lColumn < kColumn; lColumn++) {
    3399                               iColumn = order[lColumn];
    3400                               CoinBigIndex start = columnStart[iColumn];
    3401                               if (model->logLevel() == 63) {
    3402                                    printf("column %d %g <= ", iColumn, columnLower[iColumn]);
    3403                                    for (CoinBigIndex i = start; i < start + iRow; i++)
    3404                                         printf("( %d, %g) ", row[i], element[i]);
    3405                                    printf("<= %g\n", columnUpper[iColumn]);
    3406                               }
    3407                          }
    3408                     }
    3409                     jColumn = kColumn;
    3410                }
    3411           }
    3412           delete [] order;
    3413           delete [] other;
    3414           delete [] number;
    3415      }
    3416      // get row copy
    3417      CoinPackedMatrix rowCopy = *matrix;
    3418      rowCopy.reverseOrdering();
    3419      const int * rowLength = rowCopy.getVectorLengths();
    3420      number = new int[numberColumns+1];
    3421      memset(number, 0, (numberColumns + 1)*sizeof(int));
    3422      if (model->logLevel() > 3) {
    3423        // get column copy
    3424        CoinPackedMatrix columnCopy = *matrix;
    3425        const int * columnLength = columnCopy.getVectorLengths();
    3426        const int * row = columnCopy.getIndices();
    3427        const CoinBigIndex * columnStart = columnCopy.getVectorStarts();
    3428        const double * element = columnCopy.getElements();
    3429        const double * elementByRow = rowCopy.getElements();
    3430        const int * rowStart = rowCopy.getVectorStarts();
    3431        const int * column = rowCopy.getIndices();
    3432        int nPossibleZeroCost=0;
    3433        int nPossibleNonzeroCost=0;
    3434        for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
    3435          int length = columnLength[iColumn];
    3436          if (columnLower[iColumn]<-1.0e30&&columnUpper[iColumn]>1.0e30) {
    3437            if (length==1) {
    3438              printf("Singleton free %d - cost %g\n",iColumn,objective[iColumn]);
    3439            } else if (length==2) {
    3440              int iRow0=row[columnStart[iColumn]];
    3441              int iRow1=row[columnStart[iColumn]+1];
    3442              double element0=element[columnStart[iColumn]];
    3443              double element1=element[columnStart[iColumn]+1];
    3444              int n0=rowLength[iRow0];
    3445              int n1=rowLength[iRow1];
    3446              printf("Doubleton free %d - cost %g - %g in %srow with %d entries and %g in %srow with %d entries\n",
    3447                     iColumn,objective[iColumn],element0,(rowLower[iRow0]==rowUpper[iRow0]) ? "==" : "",n0,
    3448                     element1,(rowLower[iRow1]==rowUpper[iRow1]) ? "==" : "",n1);
    3449 
    3450            }
    3451          }
    3452          if (length==1) {
    3453            int iRow=row[columnStart[iColumn]];
    3454            double value=COIN_DBL_MAX;
    3455            for (int i=rowStart[iRow];i<rowStart[iRow]+rowLength[iRow];i++) {
    3456              int jColumn=column[i];
    3457              if (jColumn!=iColumn) {
    3458                if (value!=elementByRow[i]) {
    3459                  if (value==COIN_DBL_MAX) {
    3460                    value=elementByRow[i];
    3461                  } else {
    3462                    value = -COIN_DBL_MAX;
    3463                    break;
    3464                  }
    3465                }
    3466              }
    3467            }
    3468            if (!objective[iColumn]) {
    3469              if (model->logLevel() > 4)
    3470              printf("Singleton %d with no objective in row with %d elements - rhs %g,%g\n",iColumn,rowLength[iRow],rowLower[iRow],rowUpper[iRow]);
    3471              nPossibleZeroCost++;
    3472            } else if (value!=-COIN_DBL_MAX) {
    3473              if (model->logLevel() > 4)
    3474                printf("Singleton %d (%s) with objective in row %d (%s) with %d equal elements - rhs %g,%g\n",iColumn,model->getColumnName(iColumn).c_str(),
    3475                       iRow,model->getRowName(iRow).c_str(),
    3476                       rowLength[iRow],rowLower[iRow],rowUpper[iRow]);
    3477              nPossibleNonzeroCost++;
    3478            }
    3479          }
    3480        }
    3481        if (nPossibleZeroCost||nPossibleNonzeroCost)
    3482          printf("%d singletons with zero cost, %d with valid cost\n",
    3483                 nPossibleZeroCost,nPossibleNonzeroCost);
    3484        // look for DW
    3485        int * blockStart = new int [3*(numberRows+numberColumns)+1];
    3486        int * columnBlock = blockStart+numberRows;
    3487        int * nextColumn = columnBlock+numberColumns;
    3488        int * blockCount = nextColumn+numberColumns;
    3489        int * blockEls = blockCount+numberRows+1;
    3490        int * countIntegers = blockEls+numberRows;
    3491        memset(countIntegers,0,numberColumns*sizeof(int));
    3492        int direction[2]={-1,1};
    3493        int bestBreak=-1;
    3494        double bestValue=0.0;
    3495        int iPass=0;
    3496        int halfway=(numberRows+1)/2;
    3497        int firstMaster=-1;
    3498        int lastMaster=-2;
    3499        while (iPass<2) {
    3500          int increment=direction[iPass];
    3501          int start= increment>0 ? 0 : numberRows-1;
    3502          int stop=increment>0 ? numberRows : -1;
    3503          int numberBlocks=0;
    3504          int thisBestBreak=-1;
    3505          double thisBestValue=COIN_DBL_MAX;
    3506          int numberRowsDone=0;
    3507          int numberMarkedColumns=0;
    3508          int maximumBlockSize=0;
    3509          for (int i=0;i<numberRows+2*numberColumns;i++)
    3510            blockStart[i]=-1;
    3511          for (int i=0;i<numberRows+1;i++)
    3512            blockCount[i]=0;
    3513          for (int iRow=start;iRow!=stop;iRow+=increment) {
    3514            int iBlock = -1;
    3515            for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    3516              int iColumn=column[j];
    3517              int whichColumnBlock=columnBlock[iColumn];
    3518              if (whichColumnBlock>=0) {
    3519                // column marked
    3520                if (iBlock<0) {
    3521                  // put row in that block
    3522                  iBlock=whichColumnBlock;
    3523                } else if (iBlock!=whichColumnBlock) {
    3524                  // merge
    3525                  blockCount[iBlock]+=blockCount[whichColumnBlock];
    3526                  blockCount[whichColumnBlock]=0;
    3527                  int jColumn=blockStart[whichColumnBlock];
    3528                  while (jColumn>=0) {
    3529                    columnBlock[jColumn]=iBlock;
    3530                    iColumn=jColumn;
    3531                    jColumn=nextColumn[jColumn];
    3532                  }
    3533                  nextColumn[iColumn]=blockStart[iBlock];
    3534                  blockStart[iBlock]=blockStart[whichColumnBlock];
    3535                  blockStart[whichColumnBlock]=-1;
    3536                }
    3537              }
    3538            }
    3539            int n=numberMarkedColumns;
    3540            if (iBlock<0) {
    3541              //new block
    3542              if (rowLength[iRow]) {
    3543                numberBlocks++;
    3544                iBlock=numberBlocks;
    3545                int jColumn=column[rowStart[iRow]];
    3546                columnBlock[jColumn]=iBlock;
    3547                blockStart[iBlock]=jColumn;
    3548                numberMarkedColumns++;
    3549                for (CoinBigIndex j=rowStart[iRow]+1;j<rowStart[iRow]+rowLength[iRow];j++) {
    3550                  int iColumn=column[j];
    3551                  columnBlock[iColumn]=iBlock;
    3552                  numberMarkedColumns++;
    3553                  nextColumn[jColumn]=iColumn;
    3554                  jColumn=iColumn;
    3555                }
    3556                blockCount[iBlock]=numberMarkedColumns-n;
    3557              } else {
    3558                // empty
    3559                iBlock=numberRows;
    3560              }
    3561            } else {
    3562              // put in existing block
    3563              int jColumn=blockStart[iBlock];
    3564              for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    3565                int iColumn=column[j];
    3566                assert (columnBlock[iColumn]<0||columnBlock[iColumn]==iBlock);
    3567                if (columnBlock[iColumn]<0) {
    3568                  columnBlock[iColumn]=iBlock;
    3569                  numberMarkedColumns++;
    3570                  nextColumn[iColumn]=jColumn;
    3571                  jColumn=iColumn;
    3572                }
    3573              }
    3574              blockStart[iBlock]=jColumn;
    3575              blockCount[iBlock]+=numberMarkedColumns-n;
    3576            }
    3577            maximumBlockSize=CoinMax(maximumBlockSize,blockCount[iBlock]);
    3578            numberRowsDone++;
    3579            if (thisBestValue*numberRowsDone > maximumBlockSize&&numberRowsDone>halfway) {
    3580              thisBestBreak=iRow;
    3581              thisBestValue=static_cast<double>(maximumBlockSize)/
    3582                static_cast<double>(numberRowsDone);
    3583            }
    3584          }
    3585          if (thisBestBreak==stop)
    3586            thisBestValue=COIN_DBL_MAX;
    3587          iPass++;
    3588          if (iPass==1) {
    3589            bestBreak=thisBestBreak;
    3590            bestValue=thisBestValue;
    3591          } else {
    3592            if (bestValue<thisBestValue) {
    3593              firstMaster=0;
    3594              lastMaster=bestBreak;
    3595            } else {
    3596              firstMaster=thisBestBreak; // ? +1
    3597              lastMaster=numberRows;
    3598            }
    3599          }
    3600        }
    3601        if (firstMaster<lastMaster) {
    3602          printf("%d master rows %d <= < %d\n",lastMaster-firstMaster,
    3603                 firstMaster,lastMaster);
    3604          for (int i=0;i<numberRows+2*numberColumns;i++)
    3605            blockStart[i]=-1;
    3606          for (int i=firstMaster;i<lastMaster;i++)
    3607            blockStart[i]=-2;
    3608          int firstRow=0;
    3609          int numberBlocks=-1;
    3610          while (true) {
    3611            for (;firstRow<numberRows;firstRow++) {
    3612              if (blockStart[firstRow]==-1)
    3613                break;
    3614            }
    3615            if (firstRow==numberRows)
    3616              break;
    3617            int nRows=0;
    3618            numberBlocks++;
    3619            int numberStack=1;
    3620            blockCount[0] = firstRow;
    3621            while (numberStack) {
    3622              int iRow=blockCount[--numberStack];
    3623              for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
    3624                int iColumn=column[j];
    3625                int iBlock=columnBlock[iColumn];
    3626                if (iBlock<0) {
    3627                  columnBlock[iColumn]=numberBlocks;
    3628                  for (CoinBigIndex k=columnStart[iColumn];
    3629                       k<columnStart[iColumn]+columnLength[iColumn];k++) {
    3630                    int jRow=row[k];
    3631                    int rowBlock=blockStart[jRow];
    3632                    if (rowBlock==-1) {
    3633                      nRows++;
    3634                      blockStart[jRow]=numberBlocks;
    3635                      blockCount[numberStack++]=jRow;
    3636                    }
    3637                  }
    3638                }
    3639              }
    3640            }
    3641            if (!nRows) {
    3642              // empty!!
    3643              numberBlocks--;
    3644            }
    3645            firstRow++;
    3646          }
    3647          // adjust
    3648          numberBlocks++;
    3649          for (int i=0;i<numberBlocks;i++) {
    3650            blockCount[i]=0;
    3651            nextColumn[i]=0;
    3652          }
    3653          int numberEmpty=0;
    3654          int numberMaster=0;
    3655          memset(blockEls,0,numberBlocks*sizeof(int));
    3656          for (int iRow = 0; iRow < numberRows; iRow++) {
    3657            int iBlock=blockStart[iRow];
    3658            if (iBlock>=0) {
    3659              blockCount[iBlock]++;
    3660              blockEls[iBlock]+=rowLength[iRow];
    3661            } else {
    3662              if (iBlock==-2)
    3663                numberMaster++;
    3664              else
    3665                numberEmpty++;
    3666            }
    3667          }
    3668          int numberEmptyColumns=0;
    3669          int numberMasterColumns=0;
    3670          int numberMasterIntegers=0;
    3671          for (int iColumn = 0; iColumn < numberColumns; iColumn++) {
    3672            int iBlock=columnBlock[iColumn];
    3673            bool isInteger = (model->isInteger(iColumn));
    3674            if (iBlock>=0) {
    3675              nextColumn[iBlock]++;
    3676              if (isInteger)
    3677                countIntegers[iBlock]++;
    3678            } else {
    3679              if (isInteger)
    3680                numberMasterIntegers++;
    3681              if (columnLength[iColumn])
    3682                numberMasterColumns++;
    3683              else
    3684                numberEmptyColumns++;
    3685            }
    3686          }
    3687          int largestRows=0;
    3688          int largestColumns=0;
    3689          for (int i=0;i<numberBlocks;i++) {
    3690            if (blockCount[i]+nextColumn[i]>largestRows+largestColumns) {
    3691              largestRows=blockCount[i];
    3692              largestColumns=nextColumn[i];
    3693            }
    3694          }
    3695          bool useful=true;
    3696          if (numberMaster>halfway||largestRows*3>numberRows)
    3697            useful=false;
    3698          printf("%s %d blocks (largest %d,%d), %d master rows (%d empty) out of %d, %d master columns (%d empty, %d integer) out of %d\n",
    3699                 useful ? "**Useful" : "NoGood",
    3700                 numberBlocks,largestRows,largestColumns,numberMaster,numberEmpty,numberRows,
    3701                 numberMasterColumns,numberEmptyColumns,numberMasterIntegers,
    3702                 numberColumns);
    3703          for (int i=0;i<numberBlocks;i++)
    3704            printf("Block %d has %d rows and %d columns (%d elements, %d integers)\n",
    3705                   i,blockCount[i],nextColumn[i],blockEls[i],countIntegers[i]);
    3706          FILE * fpBlocks = fopen("blocks.data","wb");
    3707          printf("Blocks data on file blocks.data\n");
    3708          int stats[3];
    3709          stats[0]=numberRows;
    3710          stats[1]=numberColumns;
    3711          stats[2]=numberBlocks;
    3712          size_t numberWritten;
    3713          numberWritten=fwrite(stats,sizeof(int),3,fpBlocks);
    3714          assert (numberWritten==3);
    3715          numberWritten=fwrite(blockStart,sizeof(int),numberRows,fpBlocks);
    3716          assert (numberWritten==numberRows);
    3717          numberWritten=fwrite(columnBlock,sizeof(int),numberColumns,fpBlocks);
    3718          assert (numberWritten==numberColumns);
    3719          fclose(fpBlocks);
    3720          if (model->logLevel() == 17) {
    3721            int * whichRows=new int[numberRows+numberColumns];
    3722            int * whichColumns=whichRows+numberRows;
    3723            char name[20];
    3724            for (int iBlock=0;iBlock<numberBlocks;iBlock++) {
    3725              sprintf(name,"block%d.mps",iBlock);
    3726              int nRows=0;
    3727              for (int iRow=0;iRow<numberRows;iRow++) {
    3728                if (blockStart[iRow]==iBlock)
    3729                  whichRows[nRows++]=iRow;
    3730              }
    3731              int nColumns=0;
    3732              for (int iColumn=0;iColumn<numberColumns;iColumn++) {
    3733                if (columnBlock[iColumn]==iBlock)
    3734                  whichColumns[nColumns++]=iColumn;
    3735              }
    3736              ClpSimplex subset(model,nRows,whichRows,nColumns,whichColumns);
    3737              for (int jRow=0;jRow<nRows;jRow++) {
    3738                int iRow = whichRows[jRow];
    3739                std::string name = model->getRowName(iRow);
    3740                subset.setRowName(jRow,name);
    3741              }
    3742              for (int jColumn=0;jColumn<nColumns;jColumn++) {
    3743                int iColumn = whichColumns[jColumn];
    3744                if (model->isInteger(iColumn))
    3745                  subset.setInteger(jColumn);
    3746                std::string name = model->getColumnName(iColumn);
    3747                subset.setColumnName(jColumn,name);
    3748              }
    3749              subset.writeMps(name,0,1);
    3750            }
    3751            delete [] whichRows;
    3752          }
    3753        }
    3754        delete [] blockStart;
    3755      }
    3756      for (iRow = 0; iRow < numberRows; iRow++) {
    3757           int length = rowLength[iRow];
    3758           number[length]++;
    3759      }
    3760      if (number[0])
    3761           printf("** %d rows have no entries\n", number[0]);
    3762      k = 0;
    3763      for (iColumn = 1; iColumn <= numberColumns; iColumn++) {
    3764           if (number[iColumn]) {
    3765                k++;
    3766                printf("%d rows have %d entries\n", number[iColumn], iColumn);
    3767                if (k == kMax)
    3768                     break;
    3769           }
    3770      }
    3771      if (k < numberColumns) {
    3772           int kk = k;
    3773           k = 0;
    3774           for (iColumn = numberColumns; iColumn >= 1; iColumn--) {
    3775                if (number[iColumn]) {
    3776                     k++;
    3777                     if (k == kMax)
    3778                          break;
    3779                }
    3780           }
    3781           if (k > kk) {
    3782                printf("\n    .........\n\n");
    3783                iColumn = k;
    3784                k = 0;
    3785                for (; iColumn < numberColumns; iColumn++) {
    3786                     if (number[iColumn]) {
    3787                          k++;
    3788                          printf("%d rows have %d entries\n", number[iColumn], iColumn);
    3789                          if (k == kMax)
    3790                               break;
    3791                     }
    3792                }
    3793           }
    3794      }
    3795      if (model->logLevel() == 63
    3796 #ifdef SYM
    3797                || true
    3798 #endif
    3799         ) {
    3800           int * column = rowCopy.getMutableIndices();
    3801           const CoinBigIndex * rowStart = rowCopy.getVectorStarts();
    3802           double * element = rowCopy.getMutableElements();
    3803           int * order = new int[numberRows];
    3804           int * other = new int[numberRows];
    3805           for (iRow = 0; iRow < numberRows; iRow++) {
    3806                int length = rowLength[iRow];
    3807                order[iRow] = iRow;
    3808                other[iRow] = length;
    3809                CoinBigIndex start = rowStart[iRow];
    3810                CoinSort_2(column + start, column + start + length, element + start);
    3811           }
    3812           CoinSort_2(other, other + numberRows, order);
    3813           int jRow = number[0] + number[1];
    3814           double * weight = new double[numberRows];
    3815           double * randomColumn = new double[numberColumns+1];
    3816           double * randomRow = new double [numberRows+1];
    3817           int * sortRow = new int [numberRows];
    3818           int * possibleRow = new int [numberRows];
    3819           int * backRow = new int [numberRows];
    3820           int * stackRow = new int [numberRows];
    3821           int * sortColumn = new int [numberColumns];
    3822           int * possibleColumn = new int [numberColumns];
    3823           int * backColumn = new int [numberColumns];
    3824           int * backColumn2 = new int [numberColumns];
    3825           int * mapRow = new int [numberRows];
    3826           int * mapColumn = new int [numberColumns];
    3827           int * stackColumn = new int [numberColumns];
    3828           double randomLower = CoinDrand48();
    3829           double randomUpper = CoinDrand48();
    3830           double randomInteger = CoinDrand48();
    3831           int * startAdd = new int[numberRows+1];
    3832           int * columnAdd = new int [2*numberElements];
    3833           double * elementAdd = new double[2*numberElements];
    3834           int nAddRows = 0;
    3835           startAdd[0] = 0;
    3836           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    3837                randomColumn[iColumn] = CoinDrand48();
    3838                backColumn2[iColumn] = -1;
    3839           }
    3840           for (iColumn = 2; iColumn <= numberColumns; iColumn++) {
    3841                if (number[iColumn]) {
    3842                     printf("XX %d rows have %d entries\n", number[iColumn], iColumn);
    3843                     int kRow = jRow + number[iColumn];
    3844                     sortOnOther(column, rowStart,
    3845                                 order + jRow, other, number[iColumn], iColumn, 0);
    3846                     // Now print etc
    3847                     if (iColumn < 500000) {
    3848                          int nLook = 0;
    3849                          for (int lRow = jRow; lRow < kRow; lRow++) {
    3850                               iRow = order[lRow];
    3851                               CoinBigIndex start = rowStart[iRow];
    3852                               if (model->logLevel() == 63) {
    3853                                    printf("row %d %g <= ", iRow, rowLower[iRow]);
    3854                                    for (CoinBigIndex i = start; i < start + iColumn; i++)
    3855                                         printf("( %d, %g) ", column[i], element[i]);
    3856                                    printf("<= %g\n", rowUpper[iRow]);
    3857                               }
    3858                               int first = column[start];
    3859                               double sum = 0.0;
    3860                               for (CoinBigIndex i = start; i < start + iColumn; i++) {
    3861                                    int jColumn = column[i];
    3862                                    double value = element[i];
    3863                                    jColumn -= first;
    3864                                    assert (jColumn >= 0);
    3865                                    sum += value * randomColumn[jColumn];
    3866                               }
    3867                               if (rowLower[iRow] > -1.0e30 && rowLower[iRow])
    3868                                    sum += rowLower[iRow] * randomLower;
    3869                               else if (!rowLower[iRow])
    3870                                    sum += 1.234567e-7 * randomLower;
    3871                               if (rowUpper[iRow] < 1.0e30 && rowUpper[iRow])
    3872                                    sum += rowUpper[iRow] * randomUpper;
    3873                               else if (!rowUpper[iRow])
    3874                                    sum += 1.234567e-7 * randomUpper;
    3875                               sortRow[nLook] = iRow;
    3876                               randomRow[nLook++] = sum;
    3877                               // best way is to number unique elements and bounds and use
    3878                               if (fabs(sum) > 1.0e4)
    3879                                    sum *= 1.0e-6;
    3880                               weight[iRow] = sum;
    3881                          }
    3882                          assert (nLook <= numberRows);
    3883                          CoinSort_2(randomRow, randomRow + nLook, sortRow);
    3884                          randomRow[nLook] = COIN_DBL_MAX;
    3885                          double last = -COIN_DBL_MAX;
    3886                          int iLast = -1;
    3887                          for (int iLook = 0; iLook < nLook + 1; iLook++) {
    3888                               if (randomRow[iLook] > last) {
    3889                                    if (iLast >= 0) {
    3890                                         int n = iLook - iLast;
    3891                                         if (n > 1) {
    3892                                              //printf("%d rows possible?\n",n);
    3893                                         }
    3894                                    }
    3895                                    iLast = iLook;
    3896                                    last = randomRow[iLook];
    3897                               }
    3898                          }
    3899                     }
    3900                     jRow = kRow;
    3901                }
    3902           }
    3903           CoinPackedMatrix columnCopy = *matrix;
    3904           const int * columnLength = columnCopy.getVectorLengths();
    3905           const int * row = columnCopy.getIndices();
    3906           const CoinBigIndex * columnStart = columnCopy.getVectorStarts();
    3907           const double * elementByColumn = columnCopy.getElements();
    3908           for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    3909                int length = columnLength[iColumn];
    3910                CoinBigIndex start = columnStart[iColumn];
    3911                double sum = objective[iColumn];
    3912                if (columnLower[iColumn] > -1.0e30 && columnLower[iColumn])
    3913                     sum += columnLower[iColumn] * randomLower;
    3914                else if (!columnLower[iColumn])
    3915                     sum += 1.234567e-7 * randomLower;
    3916                if (columnUpper[iColumn] < 1.0e30 && columnUpper[iColumn])
    3917                     sum += columnUpper[iColumn] * randomUpper;
    3918                else if (!columnUpper[iColumn])
    3919                     sum += 1.234567e-7 * randomUpper;
    3920                if (model->isInteger(iColumn))
    3921                     sum += 9.87654321e-6 * randomInteger;
    3922                for (CoinBigIndex i = start; i < start + length; i++) {
    3923                     int iRow = row[i];
    3924                     sum += elementByColumn[i] * weight[iRow];
    3925                }
    3926                sortColumn[iColumn] = iColumn;
    3927                randomColumn[iColumn] = sum;
    3928           }
    3929           {
    3930                CoinSort_2(randomColumn, randomColumn + numberColumns, sortColumn);
    3931                for (iColumn = 0; iColumn < numberColumns; iColumn++) {
    3932                     int i = sortColumn[iColumn];
    3933                     backColumn[i] = iColumn;
    3934                }
    3935                randomColumn[numberColumns] = COIN_DBL_MAX;
    3936                double last = -COIN_DBL_MAX;
    3937                int iLast = -1;
    3938                for (int iLook = 0; iLook < numberColumns + 1; iLook++) {
    3939                     if (randomColumn[iLook] > last) {
    3940                          if (iLast >= 0) {
    3941                               int n = iLook - iLast;
    3942                               if (n > 1) {
    3943                                    //printf("%d columns possible?\n",n);
    3944                               }
    3945                               for (int i = iLast; i < iLook; i++) {
    3946                                    possibleColumn[sortColumn[i]] = n;
    3947                               }
    3948                          }
    3949                          iLast = iLook;
    3950                          last = randomColumn[iLook];
    3951                     }
    3952                }
    3953                for (iRow = 0; iRow < numberRows; iRow++) {
    3954                     CoinBigIndex start = rowStart[iRow];
    3955                     double sum = 0.0;
    3956                     int length = rowLength[iRow];
    3957                     for (CoinBigIndex i = start; i < start + length; i++) {
    3958                          int jColumn = column[i];
    3959                          double value = element[i];
    3960                          jColumn = backColumn[jColumn];
    3961                          sum += value * randomColumn[jColumn];
    3962                          //if (iColumn==23089||iRow==23729)
    3963                          //printf("row %d cola %d colb %d value %g rand %g sum %g\n",
    3964                          //   iRow,jColumn,column[i],value,randomColumn[jColumn],sum);
    3965                     }
    3966                     sortRow[iRow] = iRow;
    3967                     randomRow[iRow] = weight[iRow];
    3968                     randomRow[iRow] = sum;
    3969                }
    3970                CoinSort_2(randomRow, randomRow + numberRows, sortRow);
    3971                for (iRow = 0; iRow < numberRows; iRow++) {
    3972                     int i = sortRow[iRow];
    3973                     backRow[i] = iRow;
    3974                }
    3975                randomRow[numberRows] = COIN_DBL_MAX;
    3976                last = -COIN_DBL_MAX;
    3977                iLast = -1;
    3978                // Do backward indices from order
    3979                for (iRow = 0; iRow < numberRows; iRow++) {
    3980                     other[order[iRow]] = iRow;
    3981                }
    3982                for (int iLook = 0; iLook < numberRows + 1; iLook++) {
    3983                     if (randomRow[iLook] > last) {
    3984                          if (iLast >= 0) {
    3985                               int n = iLook - iLast;
    3986                               if (n > 1) {
    3987                                    //printf("%d rows possible?\n",n);
    3988                                    // Within group sort as for original "order"
    3989                                    for (int i = iLast; i < iLook; i++) {
    3990                                         int jRow = sortRow[i];
    3991                                         order[i] = other[jRow];
    3992                                    }
    3993                                    CoinSort_2(order + iLast, order + iLook, sortRow + iLast);
    3994                               }
    3995                               for (int i = iLast; i < iLook; i++) {
    3996                                    possibleRow[sortRow[i]] = n;
    3997                               }
    3998                          }
    3999                          iLast = iLook;
    4000                          last = randomRow[iLook];
    4001                     }
    4002                }
    4003                // Temp out
    4004                for (int iLook = 0; iLook < numberRows - 1000000; iLook++) {
    4005                     iRow = sortRow[iLook];
    4006                     CoinBigIndex start = rowStart[iRow];
    4007                     int length = rowLength[iRow];
    4008                     int numberPossible = possibleRow[iRow];
    4009                     for (CoinBigIndex i = start; i < start + length; i++) {
    4010                          int jColumn = column[i];
    4011                          if (possibleColumn[jColumn] != numberPossible)
    4012                               numberPossible = -1;
    4013                     }
    4014                     int n = numberPossible;
    4015                     if (numberPossible > 1) {
    4016                          //printf("pppppossible %d\n",numberPossible);
    4017                          for (int jLook = iLook + 1; jLook < iLook + numberPossible; jLook++) {
    4018                               int jRow = sortRow[jLook];
    4019                               CoinBigIndex start2 = rowStart[jRow];
    4020                               assert (numberPossible == possibleRow[jRow]);
    4021                               assert(length == rowLength[jRow]);
    4022                               for (CoinBigIndex i = start2; i < start2 + length; i++) {
    4023                                    int jColumn = column[i];
    4024                                    if (possibleColumn[jColumn] != numberPossible)
    4025                                         numberPossible = -1;
    4026                               }
    4027                          }
    4028                          if (numberPossible < 2) {
    4029                               // switch off
    4030                               for (int jLook = iLook; jLook < iLook + n; jLook++)
    4031                                    possibleRow[sortRow[jLook]] = -1;
    4032                          }
    4033                          // skip rest
    4034                          iLook += n - 1;
    4035                     } else {
    4036                          possibleRow[iRow] = -1;
    4037                     }
    4038                }
    4039                for (int iLook = 0; iLook < numberRows; iLook++) {
    4040                     iRow = sortRow[iLook];
    4041                     int numberPossible = possibleRow[iRow];
    4042                     // Only if any integers
    4043                     int numberIntegers = 0;
    4044                     CoinBigIndex start = rowStart[iRow];
    4045                     int length = rowLength[iRow];
    4046                     for (CoinBigIndex i = start; i < start + length; i++) {
    4047                          int jColumn = column[i];
    4048                          if (model->isInteger(jColumn))
    4049                               numberIntegers++;
    4050                     }
    4051                     if (numberPossible > 1 && !numberIntegers) {
    4052                          //printf("possible %d - but no integers\n",numberPossible);
    4053                     }
    4054                     if (numberPossible > 1 && (numberIntegers || false)) {
    4055                          //
    4056                          printf("possible %d - %d integers\n", numberPossible, numberIntegers);
    4057                          int lastLook = iLook;
    4058                          int nMapRow = -1;
    4059                          for (int jLook = iLook + 1; jLook < iLook + numberPossible; jLook++) {
    4060                               // stop if too many failures
    4061                               if (jLook > iLook + 10 && nMapRow < 0)
    4062                                    break;
    4063                               // Create identity mapping
    4064                               int i;
    4065                               for (i = 0; i < numberRows; i++)
    4066                                    mapRow[i] = i;
    4067                               for (i = 0; i < numberColumns; i++)
    4068                                    mapColumn[i] = i;
    4069                               int offset = jLook - iLook;
    4070                               int nStackC = 0;
    4071                               // build up row and column mapping
    4072                               int nStackR = 1;
    4073                               stackRow[0] = iLook;
    4074                               bool good = true;
    4075                               while (nStackR) {
    4076                                    nStackR--;
    4077                                    int look1 = stackRow[nStackR];
    4078                                    int look2 = look1 + offset;
    4079                                    assert (randomRow[look1] == randomRow[look2]);
    4080                                    int row1 = sortRow[look1];
    4081                                    int row2 = sortRow[look2];
    4082                                    assert (mapRow[row1] == row1);
    4083                                    assert (mapRow[row2] == row2);
    4084                                    mapRow[row1] = row2;
    4085                                    mapRow[row2] = row1;
    4086                                    CoinBigIndex start1 = rowStart[row1];
    4087                                    CoinBigIndex offset2 = rowStart[row2] - start1;
    4088                                    int length = rowLength[row1];
    4089                                    assert( length == rowLength[row2]);
    4090                                    for (CoinBigIndex i = start1; i < start1 + length; i++) {
    4091                                         int jColumn1 = column[i];
    4092                                         int jColumn2 = column[i+offset2];
    4093                                         if (randomColumn[backColumn[jColumn1]] !=
    4094                                                   randomColumn[backColumn[jColumn2]]) {
    4095                                              good = false;
    4096                                              break;
    4097                                         }
    4098                                         if (mapColumn[jColumn1] == jColumn1) {
    4099                                              // not touched
    4100                                              assert (mapColumn[jColumn2] == jColumn2);
    4101                                              if (jColumn1 != jColumn2) {
    4102                                                   // Put on stack
    4103                                                   mapColumn[jColumn1] = jColumn2;
    4104                                                   mapColumn[jColumn2] = jColumn1;
    4105                                                   stackColumn[nStackC++] = jColumn1;
    4106                                              }
    4107                                         } else {
    4108                                              if (mapColumn[jColumn1] != jColumn2 ||
    4109                                                        mapColumn[jColumn2] != jColumn1) {
    4110                                                   // bad
    4111                                                   good = false;
    4112                                                   printf("bad col\n");
    4113                                                   break;
    4114                                              }
    4115                                         }
    4116                                    }
    4117                                    if (!good)
    4118                                         break;
    4119                                    while (nStackC) {
    4120                                         nStackC--;
    4121                                         int iColumn = stackColumn[nStackC];
    4122                                         int iColumn2 = mapColumn[iColumn];
    4123                                         assert (iColumn != iColumn2);
    4124                                         int length = columnLength[iColumn];
    4125                                         assert (length == columnLength[iColumn2]);
    4126                                         CoinBigIndex start = columnStart[iColumn];
    4127                                         CoinBigIndex offset2 = columnStart[iColumn2] - start;
    4128                                         for (CoinBigIndex i = start; i < start + length; i++) {
    4129                                              int iRow = row[i];
    4130                                              int iRow2 = row[i+offset2];
    4131                                              if (mapRow[iRow] == iRow) {
    4132                                                   // First (but be careful)
    4133                                                   if (iRow != iRow2) {
    4134                                                        //mapRow[iRow]=iRow2;
    4135                                                        //mapRow[iRow2]=iRow;
    4136                                                        int iBack = backRow[iRow];
    4137                                                        int iBack2 = backRow[iRow2];
    4138                                                        if (randomRow[iBack] == randomRow[iBack2] &&
    4139                                                                  iBack2 - iBack == offset) {
    4140                                                             stackRow[nStackR++] = iBack;
    4141                                                        } else {
    4142                                                             //printf("randomRow diff - weights %g %g\n",
    4143                                                             //     weight[iRow],weight[iRow2]);
    4144                                                             // bad
    4145                                                             good = false;
    4146                                                             break;
    4147                                                        }
    4148                                                   }
    4149                                              } else {
    4150                                                   if (mapRow[iRow] != iRow2 ||
    4151                                                             mapRow[iRow2] != iRow) {
    4152                                                        // bad
    4153                                                        good = false;
    4154                                                        printf("bad row\n");
    4155                                                        break;
    4156                                                   }
    4157                                              }
    4158                                         }
    4159                                         if (!good)
    4160                                              break;
    4161                                    }
    4162                               }
    4163                               // then check OK
    4164                               if (good) {
    4165                                    for (iRow = 0; iRow < numberRows; iRow++) {
    4166                                         CoinBigIndex start = rowStart[iRow];
    4167                                         int length = rowLength[iRow];
    4168                                         if (mapRow[iRow] == iRow) {
    4169                                              for (CoinBigIndex i = start; i < start + length; i++) {
    4170                                                   int jColumn = column[i];
    4171                                                   backColumn2[jColumn] = i - start;
    4172                                              }
    4173                                              for (CoinBigIndex i = start; i < start + length; i++) {
    4174                                                   int jColumn = column[i];
    4175                                                   if (mapColumn[jColumn] != jColumn) {
    4176                                                        int jColumn2 = mapColumn[jColumn];
    4177                                                        CoinBigIndex i2 = backColumn2[jColumn2];
    4178                                                        if (i2 < 0) {
    4179                                                             good = false;
    4180                                                        } else if (element[i] != element[i2+start]) {
    4181                                                             good = false;
    4182                                                        }
    4183                                                   }
    4184                                              }
    4185                                              for (CoinBigIndex i = start; i < start + length; i++) {
    4186                                                   int jColumn = column[i];
    4187                                                   backColumn2[jColumn] = -1;
    4188                                              }
    4189                                         } else {
    4190                                              int row2 = mapRow[iRow];
    4191                                              assert (iRow = mapRow[row2]);
    4192                                              if (rowLower[iRow] != rowLower[row2] ||
    4193