source: branches/devel/Cbc/src/CbcSolver.cpp @ 600

Last change on this file since 600 was 600, checked in by forrest, 12 years ago

unit test

File size: 228.6 KB
<
Line 
1// Copyright (C) 2007, International Business Machines
2// Corporation and others.  All Rights Reserved.
3   
4#include "CbcConfig.h"
5#include "CoinPragma.hpp"
6
7#include <cassert>
8#include <cstdio>
9#include <cmath>
10#include <cfloat>
11#include <cstring>
12#include <iostream>
13
14
15#include "CoinPragma.hpp"
16#include "CoinHelperFunctions.hpp"
17// Same version as CBC
18#define CBCVERSION "1.03.00"
19
20#include "CoinMpsIO.hpp"
21#include "CoinModel.hpp"
22
23#include "ClpFactorization.hpp"
24#include "ClpQuadraticObjective.hpp"
25#include "CoinTime.hpp"
26#include "ClpSimplex.hpp"
27#include "ClpSimplexOther.hpp"
28#include "ClpSolve.hpp"
29#include "ClpPackedMatrix.hpp"
30#include "ClpPlusMinusOneMatrix.hpp"
31#include "ClpNetworkMatrix.hpp"
32#include "ClpDualRowSteepest.hpp"
33#include "ClpDualRowDantzig.hpp"
34#include "ClpLinearObjective.hpp"
35#include "ClpPrimalColumnSteepest.hpp"
36#include "ClpPrimalColumnDantzig.hpp"
37#include "ClpPresolve.hpp"
38#include "CbcOrClpParam.hpp"
39#include "OsiRowCutDebugger.hpp"
40#include "OsiChooseVariable.hpp"
41//#define CLP_MALLOC_STATISTICS
42#ifdef CLP_MALLOC_STATISTICS
43#include <malloc.h>
44#include <exception>
45#include <new>
46static double malloc_times=0.0;
47static double malloc_total=0.0;
48static int malloc_amount[]={0,32,128,256,1024,4096,16384,65536,262144,INT_MAX};
49static int malloc_n=10;
50double malloc_counts[10]={0,0,0,0,0,0,0,0,0,0};
51void * operator new (size_t size) throw (std::bad_alloc)
52{
53  malloc_times ++;
54  malloc_total += size;
55  int i;
56  for (i=0;i<malloc_n;i++) {
57    if ((int) size<=malloc_amount[i]) {
58      malloc_counts[i]++;
59      break;
60    }
61  }
62  void * p =malloc(size);
63  return p;
64}
65void operator delete (void *p) throw()
66{
67  free(p);
68}
69static void malloc_stats2()
70{
71  double average = malloc_total/malloc_times;
72  printf("count %g bytes %g - average %g\n",malloc_times,malloc_total,average);
73  for (int i=0;i<malloc_n;i++) 
74    printf("%g ",malloc_counts[i]);
75  printf("\n");
76  malloc_times=0.0;
77  malloc_total=0.0;
78  memset(malloc_counts,0,sizeof(malloc_counts));
79}
80#endif
81#ifdef DMALLOC
82#include "dmalloc.h"
83#endif
84#ifdef WSSMP_BARRIER
85#define FOREIGN_BARRIER
86#endif
87#ifdef UFL_BARRIER
88#define FOREIGN_BARRIER
89#endif
90#ifdef TAUCS_BARRIER
91#define FOREIGN_BARRIER
92#endif
93#include "CoinWarmStartBasis.hpp"
94
95#include "OsiSolverInterface.hpp"
96#include "OsiCuts.hpp"
97#include "OsiRowCut.hpp"
98#include "OsiColCut.hpp"
99#ifndef COIN_HAS_LINK
100#ifdef COIN_HAS_ASL
101#define COIN_HAS_LINK
102#endif
103#endif
104#ifdef COIN_HAS_LINK
105#include "CbcLinked.hpp"
106#endif
107#include "CglPreProcess.hpp"
108#include "CglCutGenerator.hpp"
109#include "CglGomory.hpp"
110#include "CglProbing.hpp"
111#include "CglKnapsackCover.hpp"
112#include "CglRedSplit.hpp"
113#include "CglClique.hpp"
114#include "CglFlowCover.hpp"
115#include "CglMixedIntegerRounding2.hpp"
116#include "CglTwomir.hpp"
117#include "CglDuplicateRow.hpp"
118#include "CglStored.hpp"
119#include "CglLandP.hpp"
120
121#include "CbcModel.hpp"
122#include "CbcHeuristic.hpp"
123#include "CbcHeuristicLocal.hpp"
124#include "CbcHeuristicGreedy.hpp"
125#include "CbcHeuristicFPump.hpp"
126#include "CbcHeuristicRINS.hpp"
127#include "CbcTreeLocal.hpp"
128#include "CbcCompareActual.hpp"
129#include "CbcBranchActual.hpp"
130#include  "CbcOrClpParam.hpp"
131#include  "CbcCutGenerator.hpp"
132#include  "CbcStrategy.hpp"
133#include "CbcBranchLotsize.hpp"
134
135#include "OsiClpSolverInterface.hpp"
136#ifdef COIN_HAS_ASL
137#include "Cbc_ampl.h"
138static bool usingAmpl=false;
139#endif
140static double totalTime=0.0;
141static void statistics(ClpSimplex * originalModel, ClpSimplex * model);
142static bool maskMatches(const int * starts, char ** masks,
143                        std::string & check);
144static void generateCode(CbcModel * model, const char * fileName,int type,int preProcess);
145#ifdef NDEBUG
146#undef NDEBUG
147#endif
148//#############################################################################
149// To use USERCBC or USERCLP uncomment the following define and add in your fake main program in CoinSolve.cpp
150//#define USER_HAS_FAKE_MAIN
151//  Start any fake main program
152#ifdef USER_HAS_FAKE_MAIN
153void fakeMain (ClpSimplex & model,OsiSolverInterface & osiSolver, CbcModel & babSolver);
154// Clp stuff to reduce problem
155void fakeMain2 (ClpSimplex & model,OsiClpSolverInterface & osiSolver,int options);
156#endif
157//  End any fake main program
158//#############################################################################
159
160// Allow for interrupts
161// But is this threadsafe ? (so switched off by option)
162
163#include "CoinSignal.hpp"
164static CbcModel * currentBranchModel = NULL;
165
166extern "C" {
167   static void signal_handler(int whichSignal)
168   {
169      if (currentBranchModel!=NULL) 
170         currentBranchModel->setMaximumNodes(0); // stop at next node
171      return;
172   }
173}
174
175int CbcOrClpRead_mode=1;
176FILE * CbcOrClpReadCommand=stdin;
177static bool noPrinting=false;
178static int * analyze(OsiClpSolverInterface * solverMod, int & numberChanged, double & increment,
179                     bool changeInt)
180{
181  OsiSolverInterface * solver = solverMod->clone();
182  if (0) {
183    // just get increment
184    CbcModel model(*solver);
185    model.analyzeObjective();
186    double increment2=model.getCutoffIncrement();
187    printf("initial cutoff increment %g\n",increment2);
188  }
189  const double *objective = solver->getObjCoefficients() ;
190  const double *lower = solver->getColLower() ;
191  const double *upper = solver->getColUpper() ;
192  int numberColumns = solver->getNumCols() ;
193  int numberRows = solver->getNumRows();
194  double direction = solver->getObjSense();
195  int iRow,iColumn;
196
197  // Row copy
198  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
199  const double * elementByRow = matrixByRow.getElements();
200  const int * column = matrixByRow.getIndices();
201  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
202  const int * rowLength = matrixByRow.getVectorLengths();
203
204  // Column copy
205  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
206  const double * element = matrixByCol.getElements();
207  const int * row = matrixByCol.getIndices();
208  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
209  const int * columnLength = matrixByCol.getVectorLengths();
210
211  const double * rowLower = solver->getRowLower();
212  const double * rowUpper = solver->getRowUpper();
213
214  char * ignore = new char [numberRows];
215  int * changed = new int[numberColumns];
216  int * which = new int[numberRows];
217  double * changeRhs = new double[numberRows];
218  memset(changeRhs,0,numberRows*sizeof(double));
219  memset(ignore,0,numberRows);
220  numberChanged=0;
221  int numberInteger=0;
222  for (iColumn=0;iColumn<numberColumns;iColumn++) {
223    if (upper[iColumn] > lower[iColumn]+1.0e-8&&solver->isInteger(iColumn)) 
224      numberInteger++;
225  }
226  bool finished=false;
227  while (!finished) {
228    int saveNumberChanged = numberChanged;
229    for (iRow=0;iRow<numberRows;iRow++) {
230      int numberContinuous=0;
231      double value1=0.0,value2=0.0;
232      bool allIntegerCoeff=true;
233      double sumFixed=0.0;
234      int jColumn1=-1,jColumn2=-1;
235      for (CoinBigIndex j=rowStart[iRow];j<rowStart[iRow]+rowLength[iRow];j++) {
236        int jColumn = column[j];
237        double value = elementByRow[j];
238        if (upper[jColumn] > lower[jColumn]+1.0e-8) {
239          if (!solver->isInteger(jColumn)) {
240            if (numberContinuous==0) {
241              jColumn1=jColumn;
242              value1=value;
243            } else {
244              jColumn2=jColumn;
245              value2=value;
246            }
247            numberContinuous++;
248          } else {
249            if (fabs(value-floor(value+0.5))>1.0e-12)
250              allIntegerCoeff=false;
251          }
252        } else {
253          sumFixed += lower[jColumn]*value;
254        }
255      }
256      double low = rowLower[iRow];
257      if (low>-1.0e20) {
258        low -= sumFixed;
259        if (fabs(low-floor(low+0.5))>1.0e-12)
260          allIntegerCoeff=false;
261      }
262      double up = rowUpper[iRow];
263      if (up<1.0e20) {
264        up -= sumFixed;
265        if (fabs(up-floor(up+0.5))>1.0e-12)
266          allIntegerCoeff=false;
267      }
268      if (!allIntegerCoeff)
269        continue; // can't do
270      if (numberContinuous==1) {
271        // see if really integer
272        // This does not allow for complicated cases
273        if (low==up) {
274          if (fabs(value1)>1.0e-3) {
275            value1 = 1.0/value1;
276            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
277              // integer
278              changed[numberChanged++]=jColumn1;
279              solver->setInteger(jColumn1);
280              if (upper[jColumn1]>1.0e20)
281                solver->setColUpper(jColumn1,1.0e20);
282              if (lower[jColumn1]<-1.0e20)
283                solver->setColLower(jColumn1,-1.0e20);
284            }
285          }
286        } else {
287          if (fabs(value1)>1.0e-3) {
288            value1 = 1.0/value1;
289            if (fabs(value1-floor(value1+0.5))<1.0e-12) {
290              // This constraint will not stop it being integer
291              ignore[iRow]=1;
292            }
293          }
294        }
295      } else if (numberContinuous==2) {
296        if (low==up) {
297          /* need general theory - for now just look at 2 cases -
298             1 - +- 1 one in column and just costs i.e. matching objective
299             2 - +- 1 two in column but feeds into G/L row which will try and minimize
300          */
301          if (fabs(value1)==1.0&&value1*value2==-1.0&&!lower[jColumn1]
302              &&!lower[jColumn2]) {
303            int n=0;
304            int i;
305            double objChange=direction*(objective[jColumn1]+objective[jColumn2]);
306            double bound = CoinMin(upper[jColumn1],upper[jColumn2]);
307            bound = CoinMin(bound,1.0e20);
308            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
309              int jRow = row[i];
310              double value = element[i];
311              if (jRow!=iRow) {
312                which[n++]=jRow;
313                changeRhs[jRow]=value;
314              }
315            }
316            for ( i=columnStart[jColumn1];i<columnStart[jColumn1]+columnLength[jColumn1];i++) {
317              int jRow = row[i];
318              double value = element[i];
319              if (jRow!=iRow) {
320                if (!changeRhs[jRow]) {
321                  which[n++]=jRow;
322                  changeRhs[jRow]=value;
323                } else {
324                  changeRhs[jRow]+=value;
325                }
326              }
327            }
328            if (objChange>=0.0) {
329              // see if all rows OK
330              bool good=true;
331              for (i=0;i<n;i++) {
332                int jRow = which[i];
333                double value = changeRhs[jRow];
334                if (value) {
335                  value *= bound;
336                  if (rowLength[jRow]==1) {
337                    if (value>0.0) {
338                      double rhs = rowLower[jRow];
339                      if (rhs>0.0) {
340                        double ratio =rhs/value;
341                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
342                          good=false;
343                      }
344                    } else {
345                      double rhs = rowUpper[jRow];
346                      if (rhs<0.0) {
347                        double ratio =rhs/value;
348                        if (fabs(ratio-floor(ratio+0.5))>1.0e-12)
349                          good=false;
350                      }
351                    }
352                  } else if (rowLength[jRow]==2) {
353                    if (value>0.0) {
354                      if (rowLower[jRow]>-1.0e20)
355                        good=false;
356                    } else {
357                      if (rowUpper[jRow]<1.0e20)
358                        good=false;
359                    }
360                  } else {
361                    good=false;
362                  }
363                }
364              }
365              if (good) {
366                // both can be integer
367                changed[numberChanged++]=jColumn1;
368                solver->setInteger(jColumn1);
369                if (upper[jColumn1]>1.0e20)
370                  solver->setColUpper(jColumn1,1.0e20);
371                if (lower[jColumn1]<-1.0e20)
372                  solver->setColLower(jColumn1,-1.0e20);
373                changed[numberChanged++]=jColumn2;
374                solver->setInteger(jColumn2);
375                if (upper[jColumn2]>1.0e20)
376                  solver->setColUpper(jColumn2,1.0e20);
377                if (lower[jColumn2]<-1.0e20)
378                  solver->setColLower(jColumn2,-1.0e20);
379              }
380            }
381            // clear
382            for (i=0;i<n;i++) {
383              changeRhs[which[i]]=0.0;
384            }
385          }
386        }
387      }
388    }
389    for (iColumn=0;iColumn<numberColumns;iColumn++) {
390      if (upper[iColumn] > lower[iColumn]+1.0e-8&&!solver->isInteger(iColumn)) {
391        double value;
392        value = upper[iColumn];
393        if (value<1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
394          continue;
395        value = lower[iColumn];
396        if (value>-1.0e20&&fabs(value-floor(value+0.5))>1.0e-12) 
397          continue;
398        bool integer=true;
399        for (CoinBigIndex j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
400          int iRow = row[j];
401          if (!ignore[iRow]) {
402            integer=false;
403            break;
404          }
405        }
406        if (integer) {
407          // integer
408          changed[numberChanged++]=iColumn;
409          solver->setInteger(iColumn);
410          if (upper[iColumn]>1.0e20)
411            solver->setColUpper(iColumn,1.0e20);
412          if (lower[iColumn]<-1.0e20)
413            solver->setColLower(iColumn,-1.0e20);
414        }
415      }
416    }
417    finished = numberChanged==saveNumberChanged;
418  }
419  delete [] which;
420  delete [] changeRhs;
421  delete [] ignore;
422  if (numberInteger&&!noPrinting)
423    printf("%d integer variables",numberInteger);
424  if (changeInt) {
425    if (!noPrinting) {
426      if (numberChanged)
427        printf(" and %d variables made integer\n",numberChanged);
428      else
429        printf("\n");
430    }
431    delete [] ignore;
432    //increment=0.0;
433    if (!numberChanged) {
434      delete [] changed;
435      delete solver;
436      return NULL;
437    } else {
438      for (iColumn=0;iColumn<numberColumns;iColumn++) {
439        if (solver->isInteger(iColumn))
440          solverMod->setInteger(iColumn);
441      }
442      delete solver;
443      return changed;
444    }
445  } else {
446    if (!noPrinting) {
447      if (numberChanged)
448        printf(" and %d variables could be made integer\n",numberChanged);
449      else
450        printf("\n");
451    }
452    // just get increment
453    CbcModel model(*solver);
454    if (noPrinting)
455      model.setLogLevel(0);
456    model.analyzeObjective();
457    double increment2=model.getCutoffIncrement();
458    if (increment2>increment&&increment2>0.0) {
459      if (!noPrinting)
460        printf("cutoff increment increased from %g to %g\n",increment,increment2);
461      increment=increment2;
462    }
463    delete solver;
464    numberChanged=0;
465    delete [] changed;
466    return NULL;
467  }
468}
469#ifdef COIN_HAS_LINK
470/*  Returns OsiSolverInterface (User should delete)
471    On entry numberKnapsack is maximum number of Total entries
472*/
473static OsiSolverInterface * 
474expandKnapsack(CoinModel & model, int * whichColumn, int * knapsackStart, 
475               int * knapsackRow, int &numberKnapsack,
476               CglStored & stored, int logLevel,
477               int fixedPriority, int SOSPriority)
478{
479  int maxTotal = numberKnapsack;
480  // load from coin model
481  OsiSolverLink *si = new OsiSolverLink();
482  OsiSolverInterface * finalModel=NULL;
483  si->setDefaultMeshSize(0.001);
484  // need some relative granularity
485  si->setDefaultBound(100.0);
486  si->setDefaultMeshSize(0.01);
487  si->setDefaultBound(100000.0);
488  si->setIntegerPriority(1000);
489  si->setBiLinearPriority(10000);
490  si->load(model,true,logLevel);
491  // get priorities
492  const int * priorities=model.priorities();
493  int numberColumns = model.numberColumns();
494  if (priorities) {
495    OsiObject ** objects = si->objects();
496    int numberObjects = si->numberObjects();
497    for (int iObj = 0;iObj<numberObjects;iObj++) {
498      int iColumn = objects[iObj]->columnNumber();
499      if (iColumn>=0&&iColumn<numberColumns) {
500        OsiSimpleInteger * obj =
501          dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
502        assert (obj);
503        int iPriority = priorities[iColumn];
504        if (iPriority>0)
505          objects[iObj]->setPriority(iPriority);
506      }
507    }
508    if (fixedPriority>0) {
509      si->setFixedPriority(fixedPriority);
510    }
511    if (SOSPriority<0)
512      SOSPriority=100000;
513  } 
514  CoinModel coinModel=*si->coinModel();
515  assert(coinModel.numberRows()>0);
516  int numberRows = coinModel.numberRows();
517  // Mark variables
518  int * whichKnapsack = new int [numberColumns];
519  int iRow,iColumn;
520  for (iColumn=0;iColumn<numberColumns;iColumn++) 
521    whichKnapsack[iColumn]=-1;
522  int kRow;
523  bool badModel=false;
524  // analyze
525  if (logLevel>1) {
526    for (iRow=0;iRow<numberRows;iRow++) {
527      /* Just obvious one at first
528         positive non unit coefficients
529         all integer
530         positive rowUpper
531         for now - linear (but further down in code may use nonlinear)
532         column bounds should be tight
533      */
534      //double lower = coinModel.getRowLower(iRow);
535      double upper = coinModel.getRowUpper(iRow);
536      if (upper<1.0e10) {
537        CoinModelLink triple=coinModel.firstInRow(iRow);
538        bool possible=true;
539        int n=0;
540        int n1=0;
541        while (triple.column()>=0) {
542          int iColumn = triple.column();
543          const char *  el = coinModel.getElementAsString(iRow,iColumn);
544          if (!strcmp("Numeric",el)) {
545            if (coinModel.columnLower(iColumn)==coinModel.columnUpper(iColumn)) {
546              triple=coinModel.next(triple);
547              continue; // fixed
548            }
549            double value=coinModel.getElement(iRow,iColumn);
550            if (value<0.0) {
551              possible=false;
552            } else {
553              n++;
554              if (value==1.0)
555                n1++;
556              if (coinModel.columnLower(iColumn)<0.0)
557                possible=false;
558              if (!coinModel.isInteger(iColumn))
559                possible=false;
560              if (whichKnapsack[iColumn]>=0)
561                possible=false;
562            }
563          } else {
564            possible=false; // non linear
565          } 
566          triple=coinModel.next(triple);
567        }
568        if (n-n1>1&&possible) {
569          double lower = coinModel.getRowLower(iRow);
570          double upper = coinModel.getRowUpper(iRow);
571          CoinModelLink triple=coinModel.firstInRow(iRow);
572          while (triple.column()>=0) {
573            int iColumn = triple.column();
574            lower -= coinModel.columnLower(iColumn)*triple.value();
575            upper -= coinModel.columnLower(iColumn)*triple.value();
576            triple=coinModel.next(triple);
577          }
578          printf("%d is possible %g <=",iRow,lower);
579          // print
580          triple=coinModel.firstInRow(iRow);
581          while (triple.column()>=0) {
582            int iColumn = triple.column();
583            if (coinModel.columnLower(iColumn)!=coinModel.columnUpper(iColumn))
584              printf(" (%d,el %g up %g)",iColumn,triple.value(),
585                     coinModel.columnUpper(iColumn)-coinModel.columnLower(iColumn));
586            triple=coinModel.next(triple);
587          }
588          printf(" <= %g\n",upper);
589        }
590      }
591    }
592  }
593  numberKnapsack=0;
594  for (kRow=0;kRow<numberRows;kRow++) {
595    iRow=kRow;
596    /* Just obvious one at first
597       positive non unit coefficients
598       all integer
599       positive rowUpper
600       for now - linear (but further down in code may use nonlinear)
601       column bounds should be tight
602    */
603    //double lower = coinModel.getRowLower(iRow);
604    double upper = coinModel.getRowUpper(iRow);
605    if (upper<1.0e10) {
606      CoinModelLink triple=coinModel.firstInRow(iRow);
607      bool possible=true;
608      int n=0;
609      int n1=0;
610      while (triple.column()>=0) {
611        int iColumn = triple.column();
612        const char *  el = coinModel.getElementAsString(iRow,iColumn);
613        if (!strcmp("Numeric",el)) {
614          if (coinModel.columnLower(iColumn)==coinModel.columnUpper(iColumn)) {
615            triple=coinModel.next(triple);
616            continue; // fixed
617          }
618          double value=coinModel.getElement(iRow,iColumn);
619          if (value<0.0) {
620            possible=false;
621          } else {
622            n++;
623            if (value==1.0)
624              n1++;
625            if (coinModel.columnLower(iColumn)<0.0)
626              possible=false;
627            if (!coinModel.isInteger(iColumn))
628              possible=false;
629            if (whichKnapsack[iColumn]>=0)
630              possible=false;
631          }
632        } else {
633          possible=false; // non linear
634        } 
635        triple=coinModel.next(triple);
636      }
637      if (n-n1>1&&possible) {
638        // try
639        CoinModelLink triple=coinModel.firstInRow(iRow);
640        while (triple.column()>=0) {
641          int iColumn = triple.column();
642          if (coinModel.columnLower(iColumn)!=coinModel.columnUpper(iColumn))
643            whichKnapsack[iColumn]=numberKnapsack;
644          triple=coinModel.next(triple);
645        }
646        knapsackRow[numberKnapsack++]=iRow;
647      }
648    }
649  }
650  if (logLevel>0)
651    printf("%d out of %d candidate rows are possible\n",numberKnapsack,numberRows);
652  // Check whether we can get rid of nonlinearities
653  /* mark rows
654     -2 in knapsack and other variables
655     -1 not involved
656     n only in knapsack n
657  */
658  int * markRow = new int [numberRows];
659  for (iRow=0;iRow<numberRows;iRow++) 
660    markRow[iRow]=-1;
661  int canDo=1; // OK and linear
662  for (iColumn=0;iColumn<numberColumns;iColumn++) {
663    CoinModelLink triple=coinModel.firstInColumn(iColumn);
664    int iKnapsack = whichKnapsack[iColumn];
665    bool linear=true;
666    // See if quadratic objective
667    const char * expr = coinModel.getColumnObjectiveAsString(iColumn);
668    if (strcmp(expr,"Numeric")) {
669      linear=false;
670    }
671    while (triple.row()>=0) {
672      int iRow = triple.row();
673      if (iKnapsack>=0) {
674        if (markRow[iRow]==-1) {
675          markRow[iRow]=iKnapsack;
676        } else if (markRow[iRow]!=iKnapsack) {
677          markRow[iRow]=-2;
678        }
679      }
680      const char * expr = coinModel.getElementAsString(iRow,iColumn);
681      if (strcmp(expr,"Numeric")) {
682        linear=false;
683      }
684      triple=coinModel.next(triple);
685    }
686    if (!linear) {
687      if (whichKnapsack[iColumn]<0) {
688        canDo=0;
689        break;
690      } else {
691        canDo=2;
692      }
693    }
694  }
695  int * markKnapsack = NULL;
696  double * coefficient = NULL;
697  double * linear = NULL;
698  int * whichRow = NULL;
699  int * lookupRow = NULL;
700  badModel=(canDo==0);
701  if (numberKnapsack&&canDo) {
702    /* double check - OK if
703       no nonlinear
704       nonlinear only on columns in knapsack
705       nonlinear only on columns in knapsack * ONE other - same for all in knapsack
706       AND that is only row connected to knapsack
707       (theoretically could split knapsack if two other and small numbers)
708       also ONE could be ONE expression - not just a variable
709    */
710    int iKnapsack;
711    markKnapsack = new int [numberKnapsack];
712    coefficient = new double [numberKnapsack];
713    linear = new double [numberColumns];
714    for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) 
715      markKnapsack[iKnapsack]=-1;
716    if (canDo==2) {
717      for (iRow=-1;iRow<numberRows;iRow++) {
718        int numberOdd;
719        CoinPackedMatrix * row = coinModel.quadraticRow(iRow,linear,numberOdd);
720        if (row) {
721          // see if valid
722          const double * element = row->getElements();
723          const int * column = row->getIndices();
724          const CoinBigIndex * columnStart = row->getVectorStarts();
725          const int * columnLength = row->getVectorLengths();
726          int numberLook = row->getNumCols();
727          for (int i=0;i<numberLook;i++) {
728            int iKnapsack=whichKnapsack[i];
729            if (iKnapsack<0) {
730              // might be able to swap - but for now can't have knapsack in
731              for (int j=columnStart[i];j<columnStart[i]+columnLength[i];j++) {
732                int iColumn = column[j];
733                if (whichKnapsack[iColumn]>=0) {
734                  canDo=0; // no good
735                  badModel=true;
736                  break;
737                }
738              }
739            } else {
740              // OK if in same knapsack - or maybe just one
741              int marked=markKnapsack[iKnapsack];
742              for (int j=columnStart[i];j<columnStart[i]+columnLength[i];j++) {
743                int iColumn = column[j];
744                if (whichKnapsack[iColumn]!=iKnapsack&&whichKnapsack[iColumn]>=0) {
745                  canDo=0; // no good
746                  badModel=true;
747                  break;
748                } else if (marked==-1) {
749                  markKnapsack[iKnapsack]=iColumn;
750                  marked=iColumn;
751                  coefficient[iKnapsack]=element[j];
752                  coinModel.associateElement(coinModel.columnName(iColumn),1.0);
753                } else if (marked!=iColumn) {
754                  badModel=true;
755                  canDo=0; // no good
756                  break;
757                } else {
758                  // could manage with different coefficients - but for now ...
759                  assert(coefficient[iKnapsack]==element[j]);
760                }
761              }
762            }
763          }
764          delete row;
765        }
766      }
767    }
768    if (canDo) {
769      // for any rows which are cuts
770      whichRow = new int [numberRows];
771      lookupRow = new int [numberRows];
772      bool someNonlinear=false;
773      double maxCoefficient=1.0;
774      for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
775        if (markKnapsack[iKnapsack]>=0) {
776          someNonlinear=true;
777          int iColumn = markKnapsack[iKnapsack];
778          maxCoefficient = CoinMax(maxCoefficient,fabs(coefficient[iKnapsack]*coinModel.columnUpper(iColumn)));
779        }
780      }
781      if (someNonlinear) {
782        // associate all columns to stop possible error messages
783        for (iColumn=0;iColumn<numberColumns;iColumn++) {
784          coinModel.associateElement(coinModel.columnName(iColumn),1.0);
785        }
786      }
787      ClpSimplex tempModel;
788      tempModel.loadProblem(coinModel);
789      // Create final model - first without knapsacks
790      int nCol=0;
791      int nRow=0;
792      for (iRow=0;iRow<numberRows;iRow++) {
793        if (markRow[iRow]<0) {
794          lookupRow[iRow]=nRow;
795          whichRow[nRow++]=iRow;
796        } else {
797          lookupRow[iRow]=-1;
798        }
799      }
800      for (iColumn=0;iColumn<numberColumns;iColumn++) {
801        if (whichKnapsack[iColumn]<0)
802          whichColumn[nCol++]=iColumn;
803      }
804      ClpSimplex finalModelX(&tempModel,nRow,whichRow,nCol,whichColumn,false,false,false);
805      OsiClpSolverInterface finalModelY(&finalModelX,true);
806      finalModel = finalModelY.clone();
807      finalModelY.releaseClp();
808      // Put back priorities
809      const int * priorities=model.priorities();
810      if (priorities) {
811        finalModel->findIntegers(false);
812        OsiObject ** objects = finalModel->objects();
813        int numberObjects = finalModel->numberObjects();
814        for (int iObj = 0;iObj<numberObjects;iObj++) {
815          int iColumn = objects[iObj]->columnNumber();
816          if (iColumn>=0&&iColumn<nCol) {
817            OsiSimpleInteger * obj =
818              dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
819            assert (obj);
820            int iPriority = priorities[whichColumn[iColumn]];
821            if (iPriority>0)
822              objects[iObj]->setPriority(iPriority);
823          }
824        }
825      }
826      for (iRow=0;iRow<numberRows;iRow++) {
827        whichRow[iRow]=iRow;
828      }
829      int numberOther=finalModel->getNumCols();
830      int nLargest=0;
831      int nelLargest=0;
832      int nTotal=0;
833      for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
834        iRow = knapsackRow[iKnapsack];
835        int nCreate = maxTotal;
836        int nelCreate=coinModel.expandKnapsack(iRow,nCreate,NULL,NULL,NULL,NULL);
837        if (nelCreate<0)
838          badModel=true;
839        nTotal+=nCreate;
840        nLargest = CoinMax(nLargest,nCreate);
841        nelLargest = CoinMax(nelLargest,nelCreate);
842      }
843      if (nTotal>maxTotal) 
844        badModel=true;
845      if (!badModel) {
846        // Now arrays for building
847        nelLargest = CoinMax(nelLargest,nLargest)+1;
848        double * buildObj = new double [nLargest];
849        double * buildElement = new double [nelLargest];
850        int * buildStart = new int[nLargest+1];
851        int * buildRow = new int[nelLargest];
852        // alow for integers in knapsacks
853        OsiObject ** object = new OsiObject * [numberKnapsack+nTotal];
854        int nSOS=0;
855        int nObj=numberKnapsack;
856        for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
857          knapsackStart[iKnapsack]=finalModel->getNumCols();
858          iRow = knapsackRow[iKnapsack];
859          int nCreate = 10000;
860          coinModel.expandKnapsack(iRow,nCreate,buildObj,buildStart,buildRow,buildElement);
861          // Redo row numbers
862          for (iColumn=0;iColumn<nCreate;iColumn++) {
863            for (int j=buildStart[iColumn];j<buildStart[iColumn+1];j++) {
864              int jRow=buildRow[j];
865              jRow=lookupRow[jRow];
866              assert (jRow>=0&&jRow<nRow);
867              buildRow[j]=jRow;
868            }
869          }
870          finalModel->addCols(nCreate,buildStart,buildRow,buildElement,NULL,NULL,buildObj);
871          int numberFinal=finalModel->getNumCols();
872          for (iColumn=numberOther;iColumn<numberFinal;iColumn++) {
873            if (markKnapsack[iKnapsack]<0) {
874              finalModel->setColUpper(iColumn,maxCoefficient);
875              finalModel->setInteger(iColumn);
876            } else {
877              finalModel->setColUpper(iColumn,maxCoefficient+1.0);
878              finalModel->setInteger(iColumn);
879            }
880            OsiSimpleInteger * sosObject = new OsiSimpleInteger(finalModel,iColumn);
881            sosObject->setPriority(1000000);
882            object[nObj++]=sosObject;
883            buildRow[iColumn-numberOther]=iColumn;
884            buildElement[iColumn-numberOther]=1.0;
885          }
886          if (markKnapsack[iKnapsack]<0) {
887            // convexity row
888            finalModel->addRow(numberFinal-numberOther,buildRow,buildElement,1.0,1.0);
889          } else {
890            int iColumn = markKnapsack[iKnapsack];
891            int n=numberFinal-numberOther;
892            buildRow[n]=iColumn;
893            buildElement[n++]=-fabs(coefficient[iKnapsack]);
894            // convexity row (sort of)
895            finalModel->addRow(n,buildRow,buildElement,0.0,0.0);
896            OsiSOS * sosObject = new OsiSOS(finalModel,n-1,buildRow,NULL,1);
897            sosObject->setPriority(iKnapsack+SOSPriority);
898            // Say not integral even if is (switch off heuristics)
899            sosObject->setIntegerValued(false);
900            object[nSOS++]=sosObject;
901          }
902          numberOther=numberFinal;
903        }
904        finalModel->addObjects(nObj,object);
905        for (iKnapsack=0;iKnapsack<nObj;iKnapsack++) 
906          delete object[iKnapsack];
907        delete [] object;
908        // Can we move any rows to cuts
909        const int * cutMarker = coinModel.cutMarker();
910        if (cutMarker&&0) {
911          printf("AMPL CUTS OFF until global cuts fixed\n");
912          cutMarker=NULL;
913        }
914        if (cutMarker) {
915          // Row copy
916          const CoinPackedMatrix * matrixByRow = finalModel->getMatrixByRow();
917          const double * elementByRow = matrixByRow->getElements();
918          const int * column = matrixByRow->getIndices();
919          const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
920          const int * rowLength = matrixByRow->getVectorLengths();
921         
922          const double * rowLower = finalModel->getRowLower();
923          const double * rowUpper = finalModel->getRowUpper();
924          int nDelete=0;
925          for (iRow=0;iRow<numberRows;iRow++) {
926            if (cutMarker[iRow]&&lookupRow[iRow]>=0) {
927              int jRow=lookupRow[iRow];
928              whichRow[nDelete++]=jRow;
929              int start = rowStart[jRow];
930              stored.addCut(rowLower[jRow],rowUpper[jRow],
931                            rowLength[jRow],column+start,elementByRow+start);
932            }
933          }
934          finalModel->deleteRows(nDelete,whichRow);
935        }
936        knapsackStart[numberKnapsack]=finalModel->getNumCols();
937        delete [] buildObj;
938        delete [] buildElement;
939        delete [] buildStart;
940        delete [] buildRow;
941        finalModel->writeMps("full");
942      }
943    }
944  }
945  delete [] whichKnapsack;
946  delete [] markRow;
947  delete [] markKnapsack;
948  delete [] coefficient;
949  delete [] linear;
950  delete [] whichRow;
951  delete [] lookupRow;
952  delete si;
953  si=NULL;
954  if (!badModel) {
955    return finalModel;
956  } else {
957    delete finalModel;
958    return NULL;
959  }
960}
961// Fills in original solution (coinModel length)
962static void 
963afterKnapsack(const CoinModel & coinModel2, const int * whichColumn, const int * knapsackStart, 
964                   const int * knapsackRow, int numberKnapsack,
965                   const double * knapsackSolution, double * solution, int logLevel)
966{
967  CoinModel coinModel = coinModel2;
968  int numberColumns = coinModel.numberColumns();
969  int iColumn;
970  // associate all columns to stop possible error messages
971  for (iColumn=0;iColumn<numberColumns;iColumn++) {
972    coinModel.associateElement(coinModel.columnName(iColumn),1.0);
973  }
974  CoinZeroN(solution,numberColumns);
975  int nCol=knapsackStart[0];
976  for (iColumn=0;iColumn<nCol;iColumn++) {
977    int jColumn = whichColumn[iColumn];
978    solution[jColumn]=knapsackSolution[iColumn];
979  }
980  int * buildRow = new int [numberColumns]; // wild overkill
981  double * buildElement = new double [numberColumns];
982  int iKnapsack;
983  for (iKnapsack=0;iKnapsack<numberKnapsack;iKnapsack++) {
984    int k=-1;
985    double value=0.0;
986    for (iColumn=knapsackStart[iKnapsack];iColumn<knapsackStart[iKnapsack+1];iColumn++) {
987      if (knapsackSolution[iColumn]>1.0e-5) {
988        if (k>=0) {
989          printf("Two nonzero values for knapsack %d at (%d,%g) and (%d,%g)\n",iKnapsack,
990                 k,knapsackSolution[k],iColumn,knapsackSolution[iColumn]);
991          abort();
992        }
993        k=iColumn;
994        value=floor(knapsackSolution[iColumn]+0.5);
995        assert (fabs(value-knapsackSolution[iColumn])<1.0e-5);
996      }
997    }
998    if (k>=0) {
999      int iRow = knapsackRow[iKnapsack];
1000      int nCreate = 10000;
1001      int nel=coinModel.expandKnapsack(iRow,nCreate,NULL,NULL,buildRow,buildElement,k-knapsackStart[iKnapsack]);
1002      assert (nel);
1003      if (logLevel>0) 
1004        printf("expanded column %d in knapsack %d has %d nonzero entries:\n",
1005               k-knapsackStart[iKnapsack],iKnapsack,nel);
1006      for (int i=0;i<nel;i++) {
1007        int jColumn = buildRow[i];
1008        double value = buildElement[i];
1009        if (logLevel>0) 
1010          printf("%d - original %d has value %g\n",i,jColumn,value);
1011        solution[jColumn]=value;
1012      }
1013    }
1014  }
1015  delete [] buildRow;
1016  delete [] buildElement;
1017#if 0
1018  for (iColumn=0;iColumn<numberColumns;iColumn++) {
1019    if (solution[iColumn]>1.0e-5&&coinModel.isInteger(iColumn))
1020      printf("%d %g\n",iColumn,solution[iColumn]);
1021  }
1022#endif
1023}
1024#endif
1025static int outDupRow(OsiSolverInterface * solver) 
1026{
1027  CglDuplicateRow dupCuts(solver);
1028  CglTreeInfo info;
1029  info.level = 0;
1030  info.pass = 0;
1031  int numberRows = solver->getNumRows();
1032  info.formulation_rows = numberRows;
1033  info.inTree = false;
1034  info.strengthenRow= NULL;
1035  info.pass = 0;
1036  OsiCuts cs;
1037  dupCuts.generateCuts(*solver,cs,info);
1038  const int * duplicate = dupCuts.duplicate();
1039  // Get rid of duplicate rows
1040  int * which = new int[numberRows]; 
1041  int numberDrop=0;
1042  for (int iRow=0;iRow<numberRows;iRow++) {
1043    if (duplicate[iRow]==-2||duplicate[iRow]>=0) 
1044      which[numberDrop++]=iRow;
1045  }
1046  if (numberDrop) {
1047    solver->deleteRows(numberDrop,which);
1048  }
1049  delete [] which;
1050  // see if we have any column cuts
1051  int numberColumnCuts = cs.sizeColCuts() ;
1052  const double * columnLower = solver->getColLower();
1053  const double * columnUpper = solver->getColUpper();
1054  for (int k = 0;k<numberColumnCuts;k++) {
1055    OsiColCut * thisCut = cs.colCutPtr(k) ;
1056    const CoinPackedVector & lbs = thisCut->lbs() ;
1057    const CoinPackedVector & ubs = thisCut->ubs() ;
1058    int j ;
1059    int n ;
1060    const int * which ;
1061    const double * values ;
1062    n = lbs.getNumElements() ;
1063    which = lbs.getIndices() ;
1064    values = lbs.getElements() ;
1065    for (j = 0;j<n;j++) {
1066      int iColumn = which[j] ;
1067      if (values[j]>columnLower[iColumn]) 
1068        solver->setColLower(iColumn,values[j]) ;
1069    }
1070    n = ubs.getNumElements() ;
1071    which = ubs.getIndices() ;
1072    values = ubs.getElements() ;
1073    for (j = 0;j<n;j++) {
1074      int iColumn = which[j] ;
1075      if (values[j]<columnUpper[iColumn]) 
1076        solver->setColUpper(iColumn,values[j]) ;
1077    }
1078  }
1079  return numberDrop;
1080}
1081void checkSOS(CbcModel * babModel, const OsiSolverInterface * solver)
1082{
1083#ifdef COIN_DEVELOP
1084  if (!babModel->ownObjects())
1085    return;
1086  //const double *objective = solver->getObjCoefficients() ;
1087  const double *columnLower = solver->getColLower() ;
1088  const double * columnUpper = solver->getColUpper() ;
1089  const double * solution = solver->getColSolution();
1090  //int numberColumns = solver->getNumCols() ;
1091  //int numberRows = solver->getNumRows();
1092  //double direction = solver->getObjSense();
1093  //int iRow,iColumn;
1094
1095  // Row copy
1096  CoinPackedMatrix matrixByRow(*solver->getMatrixByRow());
1097  //const double * elementByRow = matrixByRow.getElements();
1098  //const int * column = matrixByRow.getIndices();
1099  //const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
1100  const int * rowLength = matrixByRow.getVectorLengths();
1101
1102  // Column copy
1103  CoinPackedMatrix  matrixByCol(*solver->getMatrixByCol());
1104  const double * element = matrixByCol.getElements();
1105  const int * row = matrixByCol.getIndices();
1106  const CoinBigIndex * columnStart = matrixByCol.getVectorStarts();
1107  const int * columnLength = matrixByCol.getVectorLengths();
1108
1109  const double * rowLower = solver->getRowLower();
1110  const double * rowUpper = solver->getRowUpper();
1111  OsiObject ** objects = babModel->objects();
1112  int numberObjects = babModel->numberObjects();
1113  for (int iObj = 0;iObj<numberObjects;iObj++) {
1114    CbcSOS * objSOS =
1115      dynamic_cast <CbcSOS *>(objects[iObj]) ;
1116    if (objSOS) {
1117      int n=objSOS->numberMembers();
1118      const int * which = objSOS->members();
1119      const double * weight = objSOS->weights();
1120      int type = objSOS->sosType();
1121      // convexity row?
1122      int iColumn;
1123      iColumn=which[0];
1124      int j;
1125      int convex=-1;
1126      for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
1127        int iRow = row[j];
1128        double value = element[j];
1129        if (rowLower[iRow]==1.0&&rowUpper[iRow]==1.0&&
1130            value==1.0) {
1131          // possible
1132          if (rowLength[iRow]==n) {
1133            if (convex==-1)
1134              convex=iRow;
1135            else
1136              convex=-2;
1137          }
1138        }
1139      }
1140      printf ("set %d of type %d has %d members - possible convexity row %d\n",
1141              iObj,type,n,convex);
1142      for (int i=0;i<n;i++) {
1143        iColumn = which[i];
1144        int convex2=-1;
1145        for (j=columnStart[iColumn];j<columnStart[iColumn]+columnLength[iColumn];j++) {
1146          int iRow = row[j];
1147          if (iRow==convex) {
1148            double value = element[j];
1149            if (value==1.0) {
1150              convex2=iRow;
1151            }
1152          }
1153        }
1154        if (convex2<0&&convex>=0) {
1155          printf("odd convexity row\n");
1156          convex=-2;
1157        }
1158        printf("col %d has weight %g and value %g, bounds %g %g\n",
1159               iColumn,weight[i],solution[iColumn],columnLower[iColumn],
1160               columnUpper[iColumn]);
1161      }
1162    }
1163  }
1164#endif
1165}
1166int callCbc(const char * input2, OsiClpSolverInterface& solver1) 
1167{
1168  char * input = strdup(input2);
1169  int length = strlen(input);
1170  bool blank = input[0]=='0';
1171  int n=blank ? 0 : 1;
1172  for (int i=0;i<length;i++) {
1173    if (blank) {
1174      // look for next non blank
1175      if (input[i]==' ') {
1176        continue;
1177      } else {
1178        n++;
1179        blank=false;
1180      }
1181    } else {
1182      // look for next blank
1183      if (input[i]!=' ') {
1184        continue;
1185      } else {
1186        blank=true;
1187      }
1188    }
1189  }
1190  char ** argv = new char * [n+2];
1191  argv[0]=strdup("cbc");
1192  int i=0;
1193  while(input[i]==' ')
1194    i++;
1195  for (int j=0;j<n;j++) {
1196    int saveI=i;
1197    for (;i<length;i++) {
1198      // look for next blank
1199      if (input[i]!=' ') {
1200        continue;
1201      } else {
1202        break;
1203      }
1204    }
1205    input[i]='\0';
1206    argv[j+1]=strdup(input+saveI);
1207    while(input[i]==' ')
1208      i++;
1209  }
1210  argv[n+1]=strdup("-quit");
1211  free(input);
1212  totalTime=0.0;
1213  currentBranchModel = NULL;
1214  CbcOrClpRead_mode=1;
1215  CbcOrClpReadCommand=stdin;
1216  noPrinting=false;
1217  int returnCode = CbcMain(n+2,(const char **)argv,solver1,NULL);
1218  for (int k=0;k<n+2;k++)
1219    free(argv[k]);
1220  delete [] argv;
1221  return returnCode;
1222}
1223int callCbc(const char * input2)
1224{
1225  {
1226    OsiClpSolverInterface solver1;
1227     return callCbc(input2,solver1);
1228  }
1229}
1230int callCbc(const std::string input2, OsiClpSolverInterface& solver1) 
1231{
1232  char * input3 = strdup(input2.c_str());
1233  int returnCode=callCbc(input3,solver1);
1234  free(input3);
1235  return returnCode;
1236}
1237int callCbc(const std::string input2) 
1238{
1239  char * input3 = strdup(input2.c_str());
1240  OsiClpSolverInterface solver1;
1241  int returnCode=callCbc(input3,solver1);
1242  free(input3);
1243  return returnCode;
1244}
1245
1246int CbcMain (int argc, const char *argv[],
1247             OsiClpSolverInterface & originalSolver,
1248             CbcModel **miplibSolver)
1249{
1250  /* Note
1251     This is meant as a stand-alone executable to do as much of coin as possible.
1252     It should only have one solver known to it.
1253  */
1254  {
1255    double time1 = CoinCpuTime(),time2;
1256    bool goodModel=(originalSolver.getNumRows()&&originalSolver.getNumCols()) ? true : false;
1257    OsiClpSolverInterface solver1;
1258    if (goodModel)
1259      solver1=originalSolver;
1260    /*int kkkkkk=0;
1261    while (kkkkkk==0) {
1262      time2=CoinCpuTime();
1263      if (time2-time1>10.0) {
1264        printf("time %g\n",time2);
1265        time1=time2;
1266      }
1267      }*/
1268
1269    CoinSighandler_t saveSignal=SIG_DFL;
1270    // register signal handler
1271    saveSignal = signal(SIGINT,signal_handler);
1272    // Set up all non-standard stuff
1273    CbcModel model(solver1);
1274    CbcModel * babModel = NULL;
1275    model.setNumberBeforeTrust(21);
1276    int cutPass=-1234567;
1277    int tunePreProcess=5;
1278    int testOsiParameters=-1;
1279    // 0 normal, 1 from ampl or MIQP etc (2 allows cuts)
1280    int complicatedInteger=0;
1281    OsiSolverInterface * solver = model.solver();
1282    OsiClpSolverInterface * clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
1283    ClpSimplex * lpSolver = clpSolver->getModelPtr();
1284    clpSolver->messageHandler()->setLogLevel(0) ;
1285    model.messageHandler()->setLogLevel(1);
1286    // For priorities etc
1287    int * priorities=NULL;
1288    int * branchDirection=NULL;
1289    double * pseudoDown=NULL;
1290    double * pseudoUp=NULL;
1291    double * solutionIn = NULL;
1292    int * prioritiesIn = NULL;
1293    int numberSOS = 0;
1294    int * sosStart = NULL;
1295    int * sosIndices = NULL;
1296    char * sosType = NULL;
1297    double * sosReference = NULL;
1298    int * cut=NULL;
1299    int * sosPriority=NULL;
1300    CoinModel * coinModel = NULL;
1301#ifdef COIN_HAS_ASL
1302    ampl_info info;
1303    CglStored storedAmpl;
1304    CoinModel saveCoinModel;
1305    int * whichColumn = NULL;
1306    int * knapsackStart=NULL;
1307    int * knapsackRow=NULL;
1308    int numberKnapsack=0;
1309    memset(&info,0,sizeof(info));
1310    if (argc>2&&!strcmp(argv[2],"-AMPL")) {
1311      usingAmpl=true;
1312      // see if log in list
1313      noPrinting=true;
1314      for (int i=1;i<argc;i++) {
1315        if (!strncmp(argv[i],"log",3)) {
1316          char * equals = strchr(argv[i],'=');
1317          if (equals&&atoi(equals+1)>0) {
1318            noPrinting=false;
1319            info.logLevel=atoi(equals+1);
1320            // mark so won't be overWritten
1321            info.numberRows=-1234567;
1322            break;
1323          }
1324        }
1325      }
1326      int returnCode = readAmpl(&info,argc,const_cast<char **>(argv),(void **) (& coinModel));
1327      if (returnCode)
1328        return returnCode;
1329      CbcOrClpRead_mode=2; // so will start with parameters
1330      // see if log in list (including environment)
1331      for (int i=1;i<info.numberArguments;i++) {
1332        if (!strcmp(info.arguments[i],"log")) {
1333          if (i<info.numberArguments-1&&atoi(info.arguments[i+1])>0)
1334            noPrinting=false;
1335          break;
1336        }
1337      }
1338      if (noPrinting) {
1339        model.messageHandler()->setLogLevel(0);
1340        setCbcOrClpPrinting(false);
1341      }
1342      if (!noPrinting)
1343        printf("%d rows, %d columns and %d elements\n",
1344               info.numberRows,info.numberColumns,info.numberElements);
1345#ifdef COIN_HAS_LINK
1346      if (!coinModel) {
1347#endif
1348      solver->loadProblem(info.numberColumns,info.numberRows,info.starts,
1349                          info.rows,info.elements,
1350                          info.columnLower,info.columnUpper,info.objective,
1351                          info.rowLower,info.rowUpper);
1352      // take off cuts if ampl wants that
1353      if (info.cut&&0) {
1354        printf("AMPL CUTS OFF until global cuts fixed\n");
1355        info.cut=NULL;
1356      }
1357      if (info.cut) {
1358        int numberRows = info.numberRows;
1359        int * whichRow = new int [numberRows];
1360        // Row copy
1361        const CoinPackedMatrix * matrixByRow = solver->getMatrixByRow();
1362        const double * elementByRow = matrixByRow->getElements();
1363        const int * column = matrixByRow->getIndices();
1364        const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
1365        const int * rowLength = matrixByRow->getVectorLengths();
1366       
1367        const double * rowLower = solver->getRowLower();
1368        const double * rowUpper = solver->getRowUpper();
1369        int nDelete=0;
1370        for (int iRow=0;iRow<numberRows;iRow++) {
1371          if (info.cut[iRow]) {
1372            whichRow[nDelete++]=iRow;
1373            int start = rowStart[iRow];
1374            storedAmpl.addCut(rowLower[iRow],rowUpper[iRow],
1375                          rowLength[iRow],column+start,elementByRow+start);
1376          }
1377        }
1378        solver->deleteRows(nDelete,whichRow);
1379        delete [] whichRow;
1380      }
1381#ifdef COIN_HAS_LINK
1382      } else {
1383        // save
1384        saveCoinModel = *coinModel;
1385        // load from coin model
1386        OsiSolverLink solver1;
1387        OsiSolverInterface * solver2 = solver1.clone();
1388        model.assignSolver(solver2,true);
1389        OsiSolverLink * si =
1390          dynamic_cast<OsiSolverLink *>(model.solver()) ;
1391        assert (si != NULL);
1392        si->setDefaultMeshSize(0.001);
1393        // need some relative granularity
1394        si->setDefaultBound(100.0);
1395        si->setDefaultMeshSize(0.01);
1396        si->setDefaultBound(100000.0);
1397        si->setIntegerPriority(1000);
1398        si->setBiLinearPriority(10000);
1399        CoinModel * model2 = (CoinModel *) coinModel;
1400        si->load(*model2,true,info.logLevel);
1401        // redo
1402        solver = model.solver();
1403        clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
1404        lpSolver = clpSolver->getModelPtr();
1405        clpSolver->messageHandler()->setLogLevel(0) ;
1406        testOsiParameters=0;
1407        complicatedInteger=1;
1408      }
1409#endif
1410      // If we had a solution use it
1411      if (info.primalSolution) {
1412        solver->setColSolution(info.primalSolution);
1413      }
1414      // status
1415      if (info.rowStatus) {
1416        unsigned char * statusArray = lpSolver->statusArray();
1417        int i;
1418        for (i=0;i<info.numberColumns;i++)
1419          statusArray[i]=(char)info.columnStatus[i];
1420        statusArray+=info.numberColumns;
1421        for (i=0;i<info.numberRows;i++)
1422          statusArray[i]=(char)info.rowStatus[i];
1423        CoinWarmStartBasis * basis = lpSolver->getBasis();
1424        solver->setWarmStart(basis);
1425        delete basis;
1426      }
1427      freeArrays1(&info);
1428      // modify objective if necessary
1429      solver->setObjSense(info.direction);
1430      solver->setDblParam(OsiObjOffset,info.offset);
1431      // Set integer variables
1432      for (int i=info.numberColumns-info.numberBinary-info.numberIntegers;
1433           i<info.numberColumns;i++)
1434        solver->setInteger(i);
1435      goodModel=true;
1436      // change argc etc
1437      argc = info.numberArguments;
1438      argv = const_cast<const char **>(info.arguments);
1439    }
1440#endif   
1441    // default action on import
1442    int allowImportErrors=0;
1443    int keepImportNames=1;
1444    int doIdiot=-1;
1445    int outputFormat=2;
1446    int slpValue=-1;
1447    int cppValue=-1;
1448    int printOptions=0;
1449    int printMode=0;
1450    int presolveOptions=0;
1451    int substitution=3;
1452    int dualize=0;
1453    int doCrash=0;
1454    int doVector=0;
1455    int doSprint=-1;
1456    int doScaling=1;
1457    // set reasonable defaults
1458    int preSolve=5;
1459    int preProcess=1;
1460    bool useStrategy=false;
1461    bool preSolveFile=false;
1462   
1463    double djFix=1.0e100;
1464    double gapRatio=1.0e100;
1465    double tightenFactor=0.0;
1466    lpSolver->setPerturbation(50);
1467    lpSolver->messageHandler()->setPrefix(false);
1468    const char dirsep =  CoinFindDirSeparator();
1469    std::string directory = (dirsep == '/' ? "./" : ".\\");
1470    std::string defaultDirectory = directory;
1471    std::string importFile ="";
1472    std::string exportFile ="default.mps";
1473    std::string importBasisFile ="";
1474    std::string importPriorityFile ="";
1475    std::string debugFile="";
1476    std::string printMask="";
1477    double * debugValues = NULL;
1478    int numberDebugValues = -1;
1479    int basisHasValues=0;
1480    std::string exportBasisFile ="default.bas";
1481    std::string saveFile ="default.prob";
1482    std::string restoreFile ="default.prob";
1483    std::string solutionFile ="stdout";
1484    std::string solutionSaveFile ="solution.file";
1485#define CBCMAXPARAMETERS 200
1486    CbcOrClpParam parameters[CBCMAXPARAMETERS];
1487    int numberParameters ;
1488    establishParams(numberParameters,parameters) ;
1489    parameters[whichParam(BASISIN,numberParameters,parameters)].setStringValue(importBasisFile);
1490    parameters[whichParam(PRIORITYIN,numberParameters,parameters)].setStringValue(importPriorityFile);
1491    parameters[whichParam(BASISOUT,numberParameters,parameters)].setStringValue(exportBasisFile);
1492    parameters[whichParam(DEBUG,numberParameters,parameters)].setStringValue(debugFile);
1493    parameters[whichParam(PRINTMASK,numberParameters,parameters)].setStringValue(printMask);
1494    parameters[whichParam(DIRECTORY,numberParameters,parameters)].setStringValue(directory);
1495    parameters[whichParam(DUALBOUND,numberParameters,parameters)].setDoubleValue(lpSolver->dualBound());
1496    parameters[whichParam(DUALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->dualTolerance());
1497    parameters[whichParam(EXPORT,numberParameters,parameters)].setStringValue(exportFile);
1498    parameters[whichParam(IDIOT,numberParameters,parameters)].setIntValue(doIdiot);
1499    parameters[whichParam(IMPORT,numberParameters,parameters)].setStringValue(importFile);
1500    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].setDoubleValue(1.0e-8);
1501    int slog = whichParam(SOLVERLOGLEVEL,numberParameters,parameters);
1502    int log = whichParam(LOGLEVEL,numberParameters,parameters);
1503    parameters[slog].setIntValue(0);
1504    parameters[log].setIntValue(1);
1505    parameters[whichParam(MAXFACTOR,numberParameters,parameters)].setIntValue(lpSolver->factorizationFrequency());
1506    parameters[whichParam(MAXITERATION,numberParameters,parameters)].setIntValue(lpSolver->maximumIterations());
1507    parameters[whichParam(OUTPUTFORMAT,numberParameters,parameters)].setIntValue(outputFormat);
1508    parameters[whichParam(PRESOLVEPASS,numberParameters,parameters)].setIntValue(preSolve);
1509    parameters[whichParam(PERTVALUE,numberParameters,parameters)].setIntValue(lpSolver->perturbation());
1510    parameters[whichParam(PRIMALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->primalTolerance());
1511    parameters[whichParam(PRIMALWEIGHT,numberParameters,parameters)].setDoubleValue(lpSolver->infeasibilityCost());
1512    parameters[whichParam(RESTORE,numberParameters,parameters)].setStringValue(restoreFile);
1513    parameters[whichParam(SAVE,numberParameters,parameters)].setStringValue(saveFile);
1514    //parameters[whichParam(TIMELIMIT,numberParameters,parameters)].setDoubleValue(1.0e8);
1515    parameters[whichParam(TIMELIMIT_BAB,numberParameters,parameters)].setDoubleValue(1.0e8);
1516    parameters[whichParam(SOLUTION,numberParameters,parameters)].setStringValue(solutionFile);
1517    parameters[whichParam(SAVESOL,numberParameters,parameters)].setStringValue(solutionSaveFile);
1518    parameters[whichParam(SPRINT,numberParameters,parameters)].setIntValue(doSprint);
1519    parameters[whichParam(SUBSTITUTION,numberParameters,parameters)].setIntValue(substitution);
1520    parameters[whichParam(DUALIZE,numberParameters,parameters)].setIntValue(dualize);
1521    model.setNumberBeforeTrust(5);
1522    parameters[whichParam(NUMBERBEFORE,numberParameters,parameters)].setIntValue(5);
1523    parameters[whichParam(MAXNODES,numberParameters,parameters)].setIntValue(model.getMaximumNodes());
1524    model.setNumberStrong(5);
1525    parameters[whichParam(STRONGBRANCHING,numberParameters,parameters)].setIntValue(model.numberStrong());
1526    parameters[whichParam(INFEASIBILITYWEIGHT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcInfeasibilityWeight));
1527    parameters[whichParam(INTEGERTOLERANCE,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcIntegerTolerance));
1528    double normalIncrement=model.getCutoffIncrement();;
1529    parameters[whichParam(INCREMENT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcCutoffIncrement));
1530    parameters[whichParam(TESTOSI,numberParameters,parameters)].setIntValue(testOsiParameters);
1531    parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].setIntValue(1003);
1532    if (testOsiParameters>=0) {
1533      // trying nonlinear - switch off some stuff
1534      preProcess=0;
1535    }
1536    // Set up likely cut generators and defaults
1537    parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("on");
1538    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(128|64|1);
1539    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(1);
1540    parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].setIntValue(-1);
1541    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].setIntValue(100);
1542    parameters[whichParam(CUTSSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
1543    parameters[whichParam(HEURISTICSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
1544    parameters[whichParam(NODESTRATEGY,numberParameters,parameters)].setCurrentOption("fewest");
1545    int nodeStrategy=0;
1546    int doSOS=1;
1547    int verbose=0;
1548    CglGomory gomoryGen;
1549    // try larger limit
1550    gomoryGen.setLimitAtRoot(512);
1551    gomoryGen.setLimit(50);
1552    // set default action (0=off,1=on,2=root)
1553    int gomoryAction=3;
1554    parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1555
1556    CglProbing probingGen;
1557    probingGen.setUsingObjective(1);
1558    probingGen.setMaxPass(3);
1559    probingGen.setMaxPassRoot(3);
1560    // Number of unsatisfied variables to look at
1561    probingGen.setMaxProbe(10);
1562    probingGen.setMaxProbeRoot(50);
1563    // How far to follow the consequences
1564    probingGen.setMaxLook(10);
1565    probingGen.setMaxLookRoot(50);
1566    probingGen.setMaxLookRoot(10);
1567    // Only look at rows with fewer than this number of elements
1568    probingGen.setMaxElements(200);
1569    probingGen.setRowCuts(3);
1570    // set default action (0=off,1=on,2=root)
1571    int probingAction=1;
1572    parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1573
1574    CglKnapsackCover knapsackGen;
1575    //knapsackGen.switchOnExpensive();
1576    // set default action (0=off,1=on,2=root)
1577    int knapsackAction=3;
1578    parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1579
1580    CglRedSplit redsplitGen;
1581    //redsplitGen.setLimit(100);
1582    // set default action (0=off,1=on,2=root)
1583    // Off as seems to give some bad cuts
1584    int redsplitAction=0;
1585    parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption("off");
1586
1587    CglClique cliqueGen(false,true);
1588    cliqueGen.setStarCliqueReport(false);
1589    cliqueGen.setRowCliqueReport(false);
1590    cliqueGen.setMinViolation(0.1);
1591    // set default action (0=off,1=on,2=root)
1592    int cliqueAction=3;
1593    parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1594
1595    CglMixedIntegerRounding2 mixedGen;
1596    // set default action (0=off,1=on,2=root)
1597    int mixedAction=3;
1598    parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1599
1600    CglFlowCover flowGen;
1601    // set default action (0=off,1=on,2=root)
1602    int flowAction=3;
1603    parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1604
1605    CglTwomir twomirGen;
1606    twomirGen.setMaxElements(250);
1607    // set default action (0=off,1=on,2=root)
1608    int twomirAction=2;
1609    parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption("root");
1610    CglLandP landpGen;
1611    // set default action (0=off,1=on,2=root)
1612    int landpAction=0;
1613    parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption("off");
1614    // Stored cuts
1615    bool storedCuts = false;
1616
1617    bool useRounding=true;
1618    parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption("on");
1619    bool useFpump=true;
1620    parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption("on");
1621    bool useGreedy=true;
1622    parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption("on");
1623    bool useCombine=true;
1624    parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
1625    bool useLocalTree=false;
1626    bool useRINS=false;
1627    parameters[whichParam(RINS,numberParameters,parameters)].setCurrentOption("off");
1628    parameters[whichParam(COSTSTRATEGY,numberParameters,parameters)].setCurrentOption("off");
1629    int useCosts=0;
1630    // don't use input solution
1631    int useSolution=0;
1632   
1633    // total number of commands read
1634    int numberGoodCommands=0;
1635    // Set false if user does anything advanced
1636    bool defaultSettings=true;
1637
1638    // Hidden stuff for barrier
1639    int choleskyType = 0;
1640    int gamma=0;
1641    int scaleBarrier=0;
1642    int doKKT=0;
1643    int crossover=2; // do crossover unless quadratic
1644    // For names
1645    int lengthName = 0;
1646    std::vector<std::string> rowNames;
1647    std::vector<std::string> columnNames;
1648   
1649    std::string field;
1650    if (!noPrinting) {
1651      std::cout<<"Coin Cbc and Clp Solver version "<<CBCVERSION
1652               <<", build "<<__DATE__<<std::endl;
1653      // Print command line
1654      if (argc>1) {
1655        printf("command line - ");
1656        for (int i=0;i<argc;i++)
1657          printf("%s ",argv[i]);
1658        printf("\n");
1659      }
1660    }
1661    while (1) {
1662      // next command
1663      field=CoinReadGetCommand(argc,argv);
1664      // adjust field if has odd trailing characters
1665      char temp [200];
1666      strcpy(temp,field.c_str());
1667      int length = strlen(temp);
1668      for (int k=length-1;k>=0;k--) {
1669        if (temp[k]<' ')
1670          length--;
1671        else
1672          break;
1673      }
1674      temp[length]='\0';
1675      field=temp;
1676      // exit if null or similar
1677      if (!field.length()) {
1678        if (numberGoodCommands==1&&goodModel) {
1679          // we just had file name - do branch and bound
1680          field="branch";
1681        } else if (!numberGoodCommands) {
1682          // let's give the sucker a hint
1683          std::cout
1684            <<"CoinSolver takes input from arguments ( - switches to stdin)"
1685            <<std::endl
1686            <<"Enter ? for list of commands or help"<<std::endl;
1687          field="-";
1688        } else {
1689          break;
1690        }
1691      }
1692     
1693      // see if ? at end
1694      int numberQuery=0;
1695      if (field!="?"&&field!="???") {
1696        int length = field.length();
1697        int i;
1698        for (i=length-1;i>0;i--) {
1699          if (field[i]=='?') 
1700            numberQuery++;
1701          else
1702            break;
1703        }
1704        field=field.substr(0,length-numberQuery);
1705      }
1706      // find out if valid command
1707      int iParam;
1708      int numberMatches=0;
1709      int firstMatch=-1;
1710      for ( iParam=0; iParam<numberParameters; iParam++ ) {
1711        int match = parameters[iParam].matches(field);
1712        if (match==1) {
1713          numberMatches = 1;
1714          firstMatch=iParam;
1715          break;
1716        } else {
1717          if (match&&firstMatch<0)
1718            firstMatch=iParam;
1719          numberMatches += match>>1;
1720        }
1721      }
1722      if (iParam<numberParameters&&!numberQuery) {
1723        // found
1724        CbcOrClpParam found = parameters[iParam];
1725        CbcOrClpParameterType type = found.type();
1726        int valid;
1727        numberGoodCommands++;
1728        if (type==BAB&&goodModel) {
1729          // check if any integers
1730#ifdef COIN_HAS_ASL
1731          if (info.numberSos&&doSOS&&usingAmpl) {
1732            // SOS
1733            numberSOS = info.numberSos;
1734          }
1735#endif
1736          if (!lpSolver->integerInformation()&&!numberSOS&&
1737              !clpSolver->numberSOS())
1738            type=DUALSIMPLEX;
1739        }
1740        if (type==GENERALQUERY) {
1741          bool evenHidden=false;
1742          if ((verbose&8)!=0) {
1743            // even hidden
1744            evenHidden = true;
1745            verbose &= ~8;
1746          }
1747#ifdef COIN_HAS_ASL
1748          if (verbose<4&&usingAmpl)
1749            verbose +=4;
1750#endif
1751          if (verbose<4) {
1752            std::cout<<"In argument list keywords have leading - "
1753              ", -stdin or just - switches to stdin"<<std::endl;
1754            std::cout<<"One command per line (and no -)"<<std::endl;
1755            std::cout<<"abcd? gives list of possibilities, if only one + explanation"<<std::endl;
1756            std::cout<<"abcd?? adds explanation, if only one fuller help"<<std::endl;
1757            std::cout<<"abcd without value (where expected) gives current value"<<std::endl;
1758            std::cout<<"abcd value sets value"<<std::endl;
1759            std::cout<<"Commands are:"<<std::endl;
1760          } else {
1761            std::cout<<"Cbc options are set within AMPL with commands like:"<<std::endl<<std::endl;
1762            std::cout<<"         option cbc_options \"cuts=root log=2 feas=on slog=1\""<<std::endl<<std::endl;
1763            std::cout<<"only maximize, dual, primal, help and quit are recognized without ="<<std::endl;
1764          }
1765          int maxAcross=5;
1766          if ((verbose%4)!=0)
1767            maxAcross=1;
1768          int limits[]={1,51,101,151,201,251,301,351,401};
1769          std::vector<std::string> types;
1770          types.push_back("Double parameters:");
1771          types.push_back("Branch and Cut double parameters:");
1772          types.push_back("Integer parameters:");
1773          types.push_back("Branch and Cut integer parameters:");
1774          types.push_back("Keyword parameters:");
1775          types.push_back("Branch and Cut keyword parameters:");
1776          types.push_back("Actions or string parameters:");
1777          types.push_back("Branch and Cut actions:");
1778          int iType;
1779          for (iType=0;iType<8;iType++) {
1780            int across=0;
1781            if ((verbose%4)!=0)
1782              std::cout<<std::endl;
1783            std::cout<<types[iType]<<std::endl;
1784            if ((verbose&2)!=0)
1785              std::cout<<std::endl;
1786            for ( iParam=0; iParam<numberParameters; iParam++ ) {
1787              int type = parameters[iParam].type();
1788              if ((parameters[iParam].displayThis()||evenHidden)&&
1789                  type>=limits[iType]
1790                  &&type<limits[iType+1]) {
1791                // but skip if not useful for ampl (and in ampl mode)
1792                if (verbose>=4&&(parameters[iParam].whereUsed()&4)==0)
1793                  continue;
1794                if (!across) {
1795                  if ((verbose&2)==0) 
1796                    std::cout<<"  ";
1797                  else
1798                    std::cout<<"Command ";
1799                }
1800                std::cout<<parameters[iParam].matchName()<<"  ";
1801                across++;
1802                if (across==maxAcross) {
1803                  across=0;
1804                  if ((verbose%4)!=0) {
1805                    // put out description as well
1806                    if ((verbose&1)!=0) 
1807                      std::cout<<parameters[iParam].shortHelp();
1808                    std::cout<<std::endl;
1809                    if ((verbose&2)!=0) {
1810                      std::cout<<"---- description"<<std::endl;
1811                      parameters[iParam].printLongHelp();
1812                      std::cout<<"----"<<std::endl<<std::endl;
1813                    }
1814                  } else {
1815                    std::cout<<std::endl;
1816                  }
1817                }
1818              }
1819            }
1820            if (across)
1821              std::cout<<std::endl;
1822          }
1823        } else if (type==FULLGENERALQUERY) {
1824          std::cout<<"Full list of commands is:"<<std::endl;
1825          int maxAcross=5;
1826          int limits[]={1,51,101,151,201,251,301,351,401};
1827          std::vector<std::string> types;
1828          types.push_back("Double parameters:");
1829          types.push_back("Branch and Cut double parameters:");
1830          types.push_back("Integer parameters:");
1831          types.push_back("Branch and Cut integer parameters:");
1832          types.push_back("Keyword parameters:");
1833          types.push_back("Branch and Cut keyword parameters:");
1834          types.push_back("Actions or string parameters:");
1835          types.push_back("Branch and Cut actions:");
1836          int iType;
1837          for (iType=0;iType<8;iType++) {
1838            int across=0;
1839            std::cout<<types[iType]<<"  ";
1840            for ( iParam=0; iParam<numberParameters; iParam++ ) {
1841              int type = parameters[iParam].type();
1842              if (type>=limits[iType]
1843                  &&type<limits[iType+1]) {
1844                if (!across)
1845                  std::cout<<"  ";
1846                std::cout<<parameters[iParam].matchName()<<"  ";
1847                across++;
1848                if (across==maxAcross) {
1849                  std::cout<<std::endl;
1850                  across=0;
1851                }
1852              }
1853            }
1854            if (across)
1855              std::cout<<std::endl;
1856          }
1857        } else if (type<101) {
1858          // get next field as double
1859          double value = CoinReadGetDoubleField(argc,argv,&valid);
1860          if (!valid) {
1861            if (type<51) {
1862              parameters[iParam].setDoubleParameter(lpSolver,value);
1863            } else if (type<81) {
1864              parameters[iParam].setDoubleParameter(model,value);
1865            } else {
1866              parameters[iParam].setDoubleParameter(lpSolver,value);
1867              switch(type) {
1868              case DJFIX:
1869                djFix=value;
1870                preSolve=5;
1871                defaultSettings=false; // user knows what she is doing
1872                break;
1873              case GAPRATIO:
1874                gapRatio=value;
1875                break;
1876              case TIGHTENFACTOR:
1877                tightenFactor=value;
1878                if(!complicatedInteger)
1879                  defaultSettings=false; // user knows what she is doing
1880                break;
1881              default:
1882                abort();
1883              }
1884            }
1885          } else if (valid==1) {
1886            abort();
1887          } else {
1888            std::cout<<parameters[iParam].name()<<" has value "<<
1889              parameters[iParam].doubleValue()<<std::endl;
1890          }
1891        } else if (type<201) {
1892          // get next field as int
1893          int value = CoinReadGetIntField(argc,argv,&valid);
1894          if (!valid) {
1895            if (type<151) {
1896              if (parameters[iParam].type()==PRESOLVEPASS)
1897                preSolve = value;
1898              else if (parameters[iParam].type()==IDIOT)
1899                doIdiot = value;
1900              else if (parameters[iParam].type()==SPRINT)
1901                doSprint = value;
1902              else if (parameters[iParam].type()==OUTPUTFORMAT)
1903                outputFormat = value;
1904              else if (parameters[iParam].type()==SLPVALUE)
1905                slpValue = value;
1906              else if (parameters[iParam].type()==CPP)
1907                cppValue = value;
1908              else if (parameters[iParam].type()==PRESOLVEOPTIONS)
1909                presolveOptions = value;
1910              else if (parameters[iParam].type()==PRINTOPTIONS)
1911                printOptions = value;
1912              else if (parameters[iParam].type()==SUBSTITUTION)
1913                substitution = value;
1914              else if (parameters[iParam].type()==DUALIZE)
1915                dualize = value;
1916              else if (parameters[iParam].type()==CUTPASS)
1917                cutPass = value;
1918              else if (parameters[iParam].type()==PROCESSTUNE)
1919                tunePreProcess = value;
1920              else if (parameters[iParam].type()==VERBOSE)
1921                verbose = value;
1922              else if (parameters[iParam].type()==FPUMPITS)
1923                { useFpump = true;parameters[iParam].setIntValue(value);}
1924              parameters[iParam].setIntParameter(lpSolver,value);
1925            } else {
1926              parameters[iParam].setIntParameter(model,value);
1927            }
1928          } else if (valid==1) {
1929            abort();
1930          } else {
1931            std::cout<<parameters[iParam].name()<<" has value "<<
1932              parameters[iParam].intValue()<<std::endl;
1933          }
1934        } else if (type<301) {
1935          // one of several strings
1936          std::string value = CoinReadGetString(argc,argv);
1937          int action = parameters[iParam].parameterOption(value);
1938          if (action<0) {
1939            if (value!="EOL") {
1940              // no match
1941              parameters[iParam].printOptions();
1942            } else {
1943              // print current value
1944              std::cout<<parameters[iParam].name()<<" has value "<<
1945                parameters[iParam].currentOption()<<std::endl;
1946            }
1947          } else {
1948            parameters[iParam].setCurrentOption(action,!noPrinting);
1949            // for now hard wired
1950            switch (type) {
1951            case DIRECTION:
1952              if (action==0)
1953                lpSolver->setOptimizationDirection(1);
1954              else if (action==1)
1955                lpSolver->setOptimizationDirection(-1);
1956              else
1957                lpSolver->setOptimizationDirection(0);
1958              break;
1959            case DUALPIVOT:
1960              if (action==0) {
1961                ClpDualRowSteepest steep(3);
1962                lpSolver->setDualRowPivotAlgorithm(steep);
1963              } else if (action==1) {
1964                //ClpDualRowDantzig dantzig;
1965                ClpDualRowSteepest dantzig(5);
1966                lpSolver->setDualRowPivotAlgorithm(dantzig);
1967              } else if (action==2) {
1968                // partial steep
1969                ClpDualRowSteepest steep(2);
1970                lpSolver->setDualRowPivotAlgorithm(steep);
1971              } else {
1972                ClpDualRowSteepest steep;
1973                lpSolver->setDualRowPivotAlgorithm(steep);
1974              }
1975              break;
1976            case PRIMALPIVOT:
1977              if (action==0) {
1978                ClpPrimalColumnSteepest steep(3);
1979                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1980              } else if (action==1) {
1981                ClpPrimalColumnSteepest steep(0);
1982                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1983              } else if (action==2) {
1984                ClpPrimalColumnDantzig dantzig;
1985                lpSolver->setPrimalColumnPivotAlgorithm(dantzig);
1986              } else if (action==3) {
1987                ClpPrimalColumnSteepest steep(2);
1988                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1989              } else if (action==4) {
1990                ClpPrimalColumnSteepest steep(1);
1991                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1992              } else if (action==5) {
1993                ClpPrimalColumnSteepest steep(4);
1994                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1995              } else if (action==6) {
1996                ClpPrimalColumnSteepest steep(10);
1997                lpSolver->setPrimalColumnPivotAlgorithm(steep);
1998              }
1999              break;
2000            case SCALING:
2001              lpSolver->scaling(action);
2002              solver->setHintParam(OsiDoScale,action!=0,OsiHintTry);
2003              doScaling = 1-action;
2004              break;
2005            case AUTOSCALE:
2006              lpSolver->setAutomaticScaling(action!=0);
2007              break;
2008            case SPARSEFACTOR:
2009              lpSolver->setSparseFactorization((1-action)!=0);
2010              break;
2011            case BIASLU:
2012              lpSolver->factorization()->setBiasLU(action);
2013              break;
2014            case PERTURBATION:
2015              if (action==0)
2016                lpSolver->setPerturbation(50);
2017              else
2018                lpSolver->setPerturbation(100);
2019              break;
2020            case ERRORSALLOWED:
2021              allowImportErrors = action;
2022              break;
2023            case INTPRINT:
2024              printMode=action;
2025              break;
2026              //case ALGORITHM:
2027              //algorithm  = action;
2028              //defaultSettings=false; // user knows what she is doing
2029              //abort();
2030              //break;
2031            case KEEPNAMES:
2032              keepImportNames = 1-action;
2033              break;
2034            case PRESOLVE:
2035              if (action==0)
2036                preSolve = 5;
2037              else if (action==1)
2038                preSolve=0;
2039              else if (action==2)
2040                preSolve=10;
2041              else
2042                preSolveFile=true;
2043              break;
2044            case PFI:
2045              lpSolver->factorization()->setForrestTomlin(action==0);
2046              break;
2047            case CRASH:
2048              doCrash=action;
2049              break;
2050            case VECTOR:
2051              doVector=action;
2052              break;
2053            case MESSAGES:
2054              lpSolver->messageHandler()->setPrefix(action!=0);
2055              break;
2056            case CHOLESKY:
2057              choleskyType = action;
2058              break;
2059            case GAMMA:
2060              gamma=action;
2061              break;
2062            case BARRIERSCALE:
2063              scaleBarrier=action;
2064              break;
2065            case KKT:
2066              doKKT=action;
2067              break;
2068            case CROSSOVER:
2069              crossover=action;
2070              break;
2071            case SOS:
2072              doSOS=action;
2073              break;
2074            case GOMORYCUTS:
2075              defaultSettings=false; // user knows what she is doing
2076              gomoryAction = action;
2077              break;
2078            case PROBINGCUTS:
2079              defaultSettings=false; // user knows what she is doing
2080              probingAction = action;
2081              break;
2082            case KNAPSACKCUTS:
2083              defaultSettings=false; // user knows what she is doing
2084              knapsackAction = action;
2085              break;
2086            case REDSPLITCUTS:
2087              defaultSettings=false; // user knows what she is doing
2088              redsplitAction = action;
2089              break;
2090            case CLIQUECUTS:
2091              defaultSettings=false; // user knows what she is doing
2092              cliqueAction = action;
2093              break;
2094            case FLOWCUTS:
2095              defaultSettings=false; // user knows what she is doing
2096              flowAction = action;
2097              break;
2098            case MIXEDCUTS:
2099              defaultSettings=false; // user knows what she is doing
2100              mixedAction = action;
2101              break;
2102            case TWOMIRCUTS:
2103              defaultSettings=false; // user knows what she is doing
2104              twomirAction = action;
2105              break;
2106            case LANDPCUTS:
2107              defaultSettings=false; // user knows what she is doing
2108              landpAction = action;
2109              break;
2110            case ROUNDING:
2111              defaultSettings=false; // user knows what she is doing
2112              useRounding = action;
2113              break;
2114            case FPUMP:
2115              defaultSettings=false; // user knows what she is doing
2116              useFpump=action;
2117              break;
2118            case RINS:
2119              useRINS=action;
2120              break;
2121            case CUTSSTRATEGY:
2122              gomoryAction = action;
2123              probingAction = action;
2124              knapsackAction = action;
2125              cliqueAction = action;
2126              flowAction = action;
2127              mixedAction = action;
2128              twomirAction = action;
2129              //landpAction = action;
2130              parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption(action);
2131              parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption(action);
2132              parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption(action);
2133              if (!action) {
2134                redsplitAction = action;
2135                parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption(action);
2136              }
2137              parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption(action);
2138              parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption(action);
2139              parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption(action);
2140              parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption(action);
2141              if (!action) {
2142                landpAction = action;
2143                parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption(action);
2144              }
2145              break;
2146            case HEURISTICSTRATEGY:
2147              useRounding = action;
2148              useGreedy = action;
2149              useCombine = action;
2150              //useLocalTree = action;
2151              useFpump=action;
2152              parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption(action);
2153              parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption(action);
2154              parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption(action);
2155              //parameters[whichParam(LOCALTREE,numberParameters,parameters)].setCurrentOption(action);
2156              parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption(action);
2157              break;
2158            case GREEDY:
2159              defaultSettings=false; // user knows what she is doing
2160              useGreedy = action;
2161              break;
2162            case COMBINE:
2163              defaultSettings=false; // user knows what she is doing
2164              useCombine = action;
2165              break;
2166            case LOCALTREE:
2167              defaultSettings=false; // user knows what she is doing
2168              useLocalTree = action;
2169              break;
2170            case COSTSTRATEGY:
2171              useCosts=action;
2172              break;
2173            case NODESTRATEGY:
2174              nodeStrategy=action;
2175              break;
2176            case PREPROCESS:
2177              preProcess = action;
2178              break;
2179            case USESOLUTION:
2180              useSolution = action;
2181              break;
2182            default:
2183              abort();
2184            }
2185          }
2186        } else {
2187          // action
2188          if (type==EXIT) {
2189#ifdef COIN_HAS_ASL
2190            if(usingAmpl) {
2191              if (info.numberIntegers||info.numberBinary) {
2192                // integer
2193              } else {
2194                // linear
2195              }
2196              writeAmpl(&info);
2197              freeArrays2(&info);
2198              freeArgs(&info);
2199            }
2200#endif
2201            break; // stop all
2202          }
2203          switch (type) {
2204          case DUALSIMPLEX:
2205          case PRIMALSIMPLEX:
2206          case SOLVECONTINUOUS:
2207          case BARRIER:
2208            if (goodModel) {
2209              double objScale = 
2210                parameters[whichParam(OBJSCALE2,numberParameters,parameters)].doubleValue();
2211              if (objScale!=1.0) {
2212                int iColumn;
2213                int numberColumns=lpSolver->numberColumns();
2214                double * dualColumnSolution = 
2215                  lpSolver->dualColumnSolution();
2216                ClpObjective * obj = lpSolver->objectiveAsObject();
2217                assert(dynamic_cast<ClpLinearObjective *> (obj));
2218                double offset;
2219                double * objective = obj->gradient(NULL,NULL,offset,true);
2220                for (iColumn=0;iColumn<numberColumns;iColumn++) {
2221                  dualColumnSolution[iColumn] *= objScale;
2222                  objective[iColumn] *= objScale;;
2223                }
2224                int iRow;
2225                int numberRows=lpSolver->numberRows();
2226                double * dualRowSolution = 
2227                  lpSolver->dualRowSolution();
2228                for (iRow=0;iRow<numberRows;iRow++) 
2229                  dualRowSolution[iRow] *= objScale;
2230                lpSolver->setObjectiveOffset(objScale*lpSolver->objectiveOffset());
2231              }
2232              ClpSolve::SolveType method;
2233              ClpSolve::PresolveType presolveType;
2234              ClpSimplex * model2 = lpSolver;
2235              if (dualize) {
2236                model2 = ((ClpSimplexOther *) model2)->dualOfModel();
2237                printf("Dual of model has %d rows and %d columns\n",
2238                       model2->numberRows(),model2->numberColumns());
2239                model2->setOptimizationDirection(1.0);
2240              }
2241              if (noPrinting)
2242                lpSolver->setLogLevel(0);
2243              ClpSolve solveOptions;
2244              solveOptions.setPresolveActions(presolveOptions);
2245              solveOptions.setSubstitution(substitution);
2246              if (preSolve!=5&&preSolve) {
2247                presolveType=ClpSolve::presolveNumber;
2248                if (preSolve<0) {
2249                  preSolve = - preSolve;
2250                  if (preSolve<=100) {
2251                    presolveType=ClpSolve::presolveNumber;
2252                    printf("Doing %d presolve passes - picking up non-costed slacks\n",
2253                           preSolve);
2254                    solveOptions.setDoSingletonColumn(true);
2255                  } else {
2256                    preSolve -=100;
2257                    presolveType=ClpSolve::presolveNumberCost;
2258                    printf("Doing %d presolve passes - picking up costed slacks\n",
2259                           preSolve);
2260                  }
2261                } 
2262              } else if (preSolve) {
2263                presolveType=ClpSolve::presolveOn;
2264              } else {
2265                presolveType=ClpSolve::presolveOff;
2266              }
2267              solveOptions.setPresolveType(presolveType,preSolve);
2268              if (type==DUALSIMPLEX||type==SOLVECONTINUOUS) {
2269                method=ClpSolve::useDual;
2270              } else if (type==PRIMALSIMPLEX) {
2271                method=ClpSolve::usePrimalorSprint;
2272              } else {
2273                method = ClpSolve::useBarrier;
2274                if (crossover==1) {
2275                  method=ClpSolve::useBarrierNoCross;
2276                } else if (crossover==2) {
2277                  ClpObjective * obj = lpSolver->objectiveAsObject();
2278                  if (obj->type()>1) {
2279                    method=ClpSolve::useBarrierNoCross;
2280                    presolveType=ClpSolve::presolveOff;
2281                    solveOptions.setPresolveType(presolveType,preSolve);
2282                  } 
2283                }
2284              }
2285              solveOptions.setSolveType(method);
2286              if(preSolveFile)
2287                presolveOptions |= 0x40000000;
2288              solveOptions.setSpecialOption(4,presolveOptions);
2289              solveOptions.setSpecialOption(5,printOptions);
2290              if (doVector) {
2291                ClpMatrixBase * matrix = lpSolver->clpMatrix();
2292                if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
2293                  ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
2294                  clpMatrix->makeSpecialColumnCopy();
2295                }
2296              }
2297              if (method==ClpSolve::useDual) {
2298                // dual
2299                if (doCrash)
2300                  solveOptions.setSpecialOption(0,1,doCrash); // crash
2301                else if (doIdiot)
2302                  solveOptions.setSpecialOption(0,2,doIdiot); // possible idiot
2303              } else if (method==ClpSolve::usePrimalorSprint) {
2304                // primal
2305                // if slp turn everything off
2306                if (slpValue>0) {
2307                  doCrash=false;
2308                  doSprint=0;
2309                  doIdiot=-1;
2310                  solveOptions.setSpecialOption(1,10,slpValue); // slp
2311                  method=ClpSolve::usePrimal;
2312                }
2313                if (doCrash) {
2314                  solveOptions.setSpecialOption(1,1,doCrash); // crash
2315                } else if (doSprint>0) {
2316                  // sprint overrides idiot
2317                  solveOptions.setSpecialOption(1,3,doSprint); // sprint
2318                } else if (doIdiot>0) {
2319                  solveOptions.setSpecialOption(1,2,doIdiot); // idiot
2320                } else if (slpValue<=0) {
2321                  if (doIdiot==0) {
2322                    if (doSprint==0)
2323                      solveOptions.setSpecialOption(1,4); // all slack
2324                    else
2325                      solveOptions.setSpecialOption(1,9); // all slack or sprint
2326                  } else {
2327                    if (doSprint==0)
2328                      solveOptions.setSpecialOption(1,8); // all slack or idiot
2329                    else
2330                      solveOptions.setSpecialOption(1,7); // initiative
2331                  }
2332                }
2333                if (basisHasValues==-1)
2334                  solveOptions.setSpecialOption(1,11); // switch off values
2335              } else if (method==ClpSolve::useBarrier||method==ClpSolve::useBarrierNoCross) {
2336                int barrierOptions = choleskyType;
2337                if (scaleBarrier)
2338                  barrierOptions |= 8;
2339                if (doKKT)
2340                  barrierOptions |= 16;
2341                if (gamma)
2342                  barrierOptions |= 32*gamma;
2343                if (crossover==3) 
2344                  barrierOptions |= 256; // try presolve in crossover
2345                solveOptions.setSpecialOption(4,barrierOptions);
2346              }
2347#ifdef COIN_HAS_LINK
2348              OsiSolverInterface * coinSolver = model.solver();
2349              OsiSolverLink * linkSolver = dynamic_cast< OsiSolverLink*> (coinSolver);
2350              if (!linkSolver) {
2351                model2->initialSolve(solveOptions);
2352              } else {
2353                // special solver
2354                double * solution = linkSolver->nonlinearSLP(slpValue,1.0e-5);
2355                if (solution) {
2356                  memcpy(model2->primalColumnSolution(),solution,
2357                         CoinMin(model2->numberColumns(),linkSolver->coinModel()->numberColumns())*sizeof(double));
2358                  delete [] solution;
2359                } else {
2360                  printf("No nonlinear solution\n");
2361                }
2362              }
2363#else
2364              model2->initialSolve(solveOptions);
2365#endif
2366              basisHasValues=1;
2367              if (dualize) {
2368                int returnCode=((ClpSimplexOther *) lpSolver)->restoreFromDual(model2);
2369                delete model2;
2370                if (returnCode)
2371                  lpSolver->primal(1);
2372                model2=lpSolver;
2373              }
2374#ifdef COIN_HAS_ASL
2375              if (usingAmpl) {
2376                double value = model2->getObjValue()*model2->getObjSense();
2377                char buf[300];
2378                int pos=0;
2379                int iStat = model2->status();
2380                if (iStat==0) {
2381                  pos += sprintf(buf+pos,"optimal," );
2382                } else if (iStat==1) {
2383                  // infeasible
2384                  pos += sprintf(buf+pos,"infeasible,");
2385                } else if (iStat==2) {
2386                  // unbounded
2387                  pos += sprintf(buf+pos,"unbounded,");
2388                } else if (iStat==3) {
2389                  pos += sprintf(buf+pos,"stopped on iterations or time,");
2390                } else if (iStat==4) {
2391                  iStat = 7;
2392                  pos += sprintf(buf+pos,"stopped on difficulties,");
2393                } else if (iStat==5) {
2394                  iStat = 3;
2395                  pos += sprintf(buf+pos,"stopped on ctrl-c,");
2396                } else {
2397                  pos += sprintf(buf+pos,"status unknown,");
2398                  iStat=6;
2399                }
2400                info.problemStatus=iStat;
2401                info.objValue = value;
2402                pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
2403                               value);
2404                sprintf(buf+pos,"\n%d iterations",
2405                        model2->getIterationCount());
2406                free(info.primalSolution);
2407                int numberColumns=model2->numberColumns();
2408                info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
2409                CoinCopyN(model2->primalColumnSolution(),numberColumns,info.primalSolution);
2410                int numberRows = model2->numberRows();
2411                free(info.dualSolution);
2412                info.dualSolution = (double *) malloc(numberRows*sizeof(double));
2413                CoinCopyN(model2->dualRowSolution(),numberRows,info.dualSolution);
2414                CoinWarmStartBasis * basis = model2->getBasis();
2415                free(info.rowStatus);
2416                info.rowStatus = (int *) malloc(numberRows*sizeof(int));
2417                free(info.columnStatus);
2418                info.columnStatus = (int *) malloc(numberColumns*sizeof(int));
2419                // Put basis in
2420                int i;
2421                // free,basic,ub,lb are 0,1,2,3
2422                for (i=0;i<numberRows;i++) {
2423                  CoinWarmStartBasis::Status status = basis->getArtifStatus(i);
2424                  info.rowStatus[i]=status;
2425                }
2426                for (i=0;i<numberColumns;i++) {
2427                  CoinWarmStartBasis::Status status = basis->getStructStatus(i);
2428                  info.columnStatus[i]=status;
2429                }
2430                // put buffer into info
2431                strcpy(info.buffer,buf);
2432                delete basis;
2433              }
2434#endif
2435            } else {
2436              std::cout<<"** Current model not valid"<<std::endl;
2437            }
2438            break;
2439          case STATISTICS:
2440            if (goodModel) {
2441              // If presolve on look at presolved
2442              bool deleteModel2=false;
2443              ClpSimplex * model2 = lpSolver;
2444              if (preSolve) {
2445                ClpPresolve pinfo;
2446                int presolveOptions2 = presolveOptions&~0x40000000;
2447                if ((presolveOptions2&0xffff)!=0)
2448                  pinfo.setPresolveActions(presolveOptions2);
2449                pinfo.setSubstitution(substitution);
2450                if ((printOptions&1)!=0)
2451                  pinfo.statistics();
2452                double presolveTolerance = 
2453                  parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
2454                model2 = 
2455                  pinfo.presolvedModel(*lpSolver,presolveTolerance,
2456                                       true,preSolve);
2457                if (model2) {
2458                  printf("Statistics for presolved model\n");
2459                  deleteModel2=true;
2460                } else {
2461                  printf("Presolved model looks infeasible - will use unpresolved\n");
2462                  model2 = lpSolver;
2463                }
2464              } else {
2465                printf("Statistics for unpresolved model\n");
2466                model2 =  lpSolver;
2467              }
2468              statistics(lpSolver,model2);
2469              if (deleteModel2)
2470                delete model2;
2471            } else {
2472              std::cout<<"** Current model not valid"<<std::endl;
2473            }
2474            break;
2475          case TIGHTEN:
2476            if (goodModel) {
2477              int numberInfeasibilities = lpSolver->tightenPrimalBounds();
2478              if (numberInfeasibilities)
2479                std::cout<<"** Analysis indicates model infeasible"<<std::endl;
2480            } else {
2481              std::cout<<"** Current model not valid"<<std::endl;
2482            }
2483            break;
2484          case PLUSMINUS:
2485            if (goodModel) {
2486              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
2487              ClpPackedMatrix* clpMatrix =
2488                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
2489              if (clpMatrix) {
2490                ClpPlusMinusOneMatrix * newMatrix = new ClpPlusMinusOneMatrix(*(clpMatrix->matrix()));
2491                if (newMatrix->getIndices()) {
2492                  lpSolver->replaceMatrix(newMatrix);
2493                  delete saveMatrix;
2494                  std::cout<<"Matrix converted to +- one matrix"<<std::endl;
2495                } else {
2496                  std::cout<<"Matrix can not be converted to +- 1 matrix"<<std::endl;
2497                }
2498              } else {
2499                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
2500              }
2501            } else {
2502              std::cout<<"** Current model not valid"<<std::endl;
2503            }
2504            break;
2505          case OUTDUPROWS:
2506            if (goodModel) {
2507              int numberRows = clpSolver->getNumRows();
2508              //int nOut = outDupRow(clpSolver);
2509              CglDuplicateRow dupcuts(clpSolver);
2510              storedCuts = dupcuts.outDuplicates(clpSolver)!=0;
2511              int nOut = numberRows-clpSolver->getNumRows();
2512              if (nOut&&!noPrinting)
2513                printf("%d rows eliminated\n",nOut);
2514            } else {
2515              std::cout<<"** Current model not valid"<<std::endl;
2516            }
2517            break;
2518          case NETWORK:
2519            if (goodModel) {
2520              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
2521              ClpPackedMatrix* clpMatrix =
2522                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
2523              if (clpMatrix) {
2524                ClpNetworkMatrix * newMatrix = new ClpNetworkMatrix(*(clpMatrix->matrix()));
2525                if (newMatrix->getIndices()) {
2526                  lpSolver->replaceMatrix(newMatrix);
2527                  delete saveMatrix;
2528                  std::cout<<"Matrix converted to network matrix"<<std::endl;
2529                } else {
2530                  std::cout<<"Matrix can not be converted to network matrix"<<std::endl;
2531                }
2532              } else {
2533                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
2534              }
2535            } else {
2536              std::cout<<"** Current model not valid"<<std::endl;
2537            }
2538            break;
2539          case MIPLIB:
2540            // User can set options - main difference is lack of model and CglPreProcess
2541            goodModel=true;
2542/*
2543  Run branch-and-cut. First set a few options -- node comparison, scaling. If
2544  the solver is Clp, consider running some presolve code (not yet converted
2545  this to generic OSI) with branch-and-cut. If presolve is disabled, or the
2546  solver is not Clp, simply run branch-and-cut. Print elapsed time at the end.
2547*/
2548          case BAB: // branchAndBound
2549          case STRENGTHEN:
2550            if (goodModel) {
2551              bool miplib = type==MIPLIB;
2552              int logLevel = parameters[slog].intValue();
2553              // Reduce printout
2554              if (logLevel<=1)
2555                model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
2556#ifdef COIN_HAS_LINK
2557              {
2558                OsiSolverInterface * solver = model.solver();
2559                OsiClpSolverInterface * si =
2560                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2561                assert (si != NULL);
2562                // See if quadratic
2563                if (!complicatedInteger) {
2564                  ClpSimplex * lpSolver = si->getModelPtr();
2565                  ClpQuadraticObjective * obj = (dynamic_cast< ClpQuadraticObjective*>(lpSolver->objectiveAsObject()));
2566                  if (obj) {
2567                    preProcess=0;
2568                    parameters[whichParam(TESTOSI,numberParameters,parameters)].setIntValue(0);
2569                    // create coin model
2570                    coinModel = lpSolver->createCoinModel();
2571                    assert (coinModel);
2572                    // load from coin model
2573                    OsiSolverLink solver1;
2574                    OsiSolverInterface * solver2 = solver1.clone();
2575                    model.assignSolver(solver2,true);
2576                    OsiSolverLink * si =
2577                      dynamic_cast<OsiSolverLink *>(model.solver()) ;
2578                    assert (si != NULL);
2579                    si->setDefaultMeshSize(0.001);
2580                    // need some relative granularity
2581                    si->setDefaultBound(100.0);
2582                    si->setDefaultMeshSize(0.01);
2583                    si->setDefaultBound(1000.0);
2584                    si->setIntegerPriority(1000);
2585                    si->setBiLinearPriority(10000);
2586                    si->setSpecialOptions2(2+4+8);
2587                    CoinModel * model2 = (CoinModel *) coinModel;
2588                    si->load(*model2,true,info.logLevel);
2589                    // redo
2590                    solver = model.solver();
2591                    clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
2592                    lpSolver = clpSolver->getModelPtr();
2593                    clpSolver->messageHandler()->setLogLevel(0) ;
2594                    testOsiParameters=0;
2595                    complicatedInteger=2;  // allow cuts
2596                    OsiSolverInterface * coinSolver = model.solver();
2597                    OsiSolverLink * linkSolver = dynamic_cast< OsiSolverLink*> (coinSolver);
2598                    if (linkSolver->quadraticModel()) {
2599                      ClpSimplex * qp = linkSolver->quadraticModel();
2600                      //linkSolver->nonlinearSLP(CoinMax(slpValue,10),1.0e-5);
2601                      qp->nonlinearSLP(CoinMax(slpValue,20),1.0e-5);
2602                      qp->primal(1);
2603                      OsiSolverLinearizedQuadratic solver2(qp);
2604                      const double * solution=NULL;
2605                      // Reduce printout
2606                      solver2.setHintParam(OsiDoReducePrint,true,OsiHintTry);
2607                      CbcModel model2(solver2);
2608                      // Now do requested saves and modifications
2609                      CbcModel * cbcModel = & model2;
2610                      OsiSolverInterface * osiModel = model2.solver();
2611                      OsiClpSolverInterface * osiclpModel = dynamic_cast< OsiClpSolverInterface*> (osiModel);
2612                      ClpSimplex * clpModel = osiclpModel->getModelPtr();
2613                     
2614                      // Set changed values
2615                     
2616                      CglProbing probing;
2617                      probing.setMaxProbe(10);
2618                      probing.setMaxLook(10);
2619                      probing.setMaxElements(200);
2620                      probing.setMaxProbeRoot(50);
2621                      probing.setMaxLookRoot(10);
2622                      probing.setRowCuts(3);
2623                      probing.setUsingObjective(true);
2624                      cbcModel->addCutGenerator(&probing,-1,"Probing",true,false,false,-100,-1,-1);
2625                      cbcModel->cutGenerator(0)->setTiming(true);
2626                     
2627                      CglGomory gomory;
2628                      gomory.setLimitAtRoot(512);
2629                      cbcModel->addCutGenerator(&gomory,-98,"Gomory",true,false,false,-100,-1,-1);
2630                      cbcModel->cutGenerator(1)->setTiming(true);
2631                     
2632                      CglKnapsackCover knapsackCover;
2633                      cbcModel->addCutGenerator(&knapsackCover,-98,"KnapsackCover",true,false,false,-100,-1,-1);
2634                      cbcModel->cutGenerator(2)->setTiming(true);
2635                     
2636                      CglRedSplit redSplit;
2637                      cbcModel->addCutGenerator(&redSplit,-99,"RedSplit",true,false,false,-100,-1,-1);
2638                      cbcModel->cutGenerator(3)->setTiming(true);
2639                     
2640                      CglClique clique;
2641                      clique.setStarCliqueReport(false);
2642                      clique.setRowCliqueReport(false);
2643                      clique.setMinViolation(0.1);
2644                      cbcModel->addCutGenerator(&clique,-98,"Clique",true,false,false,-100,-1,-1);
2645                      cbcModel->cutGenerator(4)->setTiming(true);
2646                     
2647                      CglMixedIntegerRounding2 mixedIntegerRounding2;
2648                      cbcModel->addCutGenerator(&mixedIntegerRounding2,-98,"MixedIntegerRounding2",true,false,false,-100,-1,-1);
2649                      cbcModel->cutGenerator(5)->setTiming(true);
2650                     
2651                      CglFlowCover flowCover;
2652                      cbcModel->addCutGenerator(&flowCover,-98,"FlowCover",true,false,false,-100,-1,-1);
2653                      cbcModel->cutGenerator(6)->setTiming(true);
2654                     
2655                      CglTwomir twomir;
2656                      twomir.setMaxElements(250);
2657                      cbcModel->addCutGenerator(&twomir,-99,"Twomir",true,false,false,-100,-1,-1);
2658                      cbcModel->cutGenerator(7)->setTiming(true);
2659                     
2660                      CbcHeuristicFPump heuristicFPump(*cbcModel);
2661                      heuristicFPump.setWhen(13);
2662                      heuristicFPump.setMaximumPasses(20);
2663                      heuristicFPump.setMaximumRetries(7);
2664                      heuristicFPump.setAbsoluteIncrement(4332.64);
2665                      cbcModel->addHeuristic(&heuristicFPump);
2666                      heuristicFPump.setInitialWeight(1);
2667                     
2668                      CbcRounding rounding(*cbcModel);
2669                      cbcModel->addHeuristic(&rounding);
2670                     
2671                      CbcHeuristicLocal heuristicLocal(*cbcModel);
2672                      heuristicLocal.setSearchType(1);
2673                      cbcModel->addHeuristic(&heuristicLocal);
2674                     
2675                      CbcHeuristicGreedyCover heuristicGreedyCover(*cbcModel);
2676                      cbcModel->addHeuristic(&heuristicGreedyCover);
2677                     
2678                      CbcHeuristicGreedyEquality heuristicGreedyEquality(*cbcModel);
2679                      cbcModel->addHeuristic(&heuristicGreedyEquality);
2680                     
2681                      CbcCompareDefault compare;
2682                      cbcModel->setNodeComparison(compare);
2683                      cbcModel->setNumberBeforeTrust(5);
2684                      cbcModel->setSpecialOptions(2);
2685                      cbcModel->messageHandler()->setLogLevel(1);
2686                      cbcModel->setMaximumCutPassesAtRoot(-100);
2687                      cbcModel->setMaximumCutPasses(1);
2688                      cbcModel->setMinimumDrop(0.05);
2689                      // For branchAndBound this may help
2690                      clpModel->defaultFactorizationFrequency();
2691                      clpModel->setDualBound(1.0001e+08);
2692                      clpModel->setPerturbation(50);
2693                      osiclpModel->setSpecialOptions(193);
2694                      osiclpModel->messageHandler()->setLogLevel(0);
2695                      osiclpModel->setIntParam(OsiMaxNumIterationHotStart,100);
2696                      osiclpModel->setHintParam(OsiDoReducePrint,true,OsiHintTry);
2697                      // You can save some time by switching off message building
2698                      // clpModel->messagesPointer()->setDetailMessages(100,10000,(int *) NULL);
2699                     
2700                      // Solve
2701                     
2702                      cbcModel->initialSolve();
2703                      if (clpModel->tightenPrimalBounds()!=0) {
2704                        std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
2705                        exit(1);
2706                      }
2707                      clpModel->dual();  // clean up
2708                      cbcModel->initialSolve();
2709                      cbcModel->branchAndBound();
2710                      OsiSolverLinearizedQuadratic * solver3 = dynamic_cast<OsiSolverLinearizedQuadratic *> (model2.solver());
2711                      assert (solver3);
2712                      solution = solver3->bestSolution();
2713                      double bestObjectiveValue = solver3->bestObjectiveValue();
2714                      linkSolver->setBestObjectiveValue(bestObjectiveValue);
2715                      linkSolver->setBestSolution(solution,solver3->getNumCols());
2716                      CbcHeuristicDynamic3 dynamic(model);
2717                      model.addHeuristic(&dynamic);
2718                      // if convex
2719                      if ((linkSolver->specialOptions2()&4)!=0) {
2720                        int numberColumns = coinModel->numberColumns();
2721                        // add OA cut
2722                        double offset;
2723                        double * gradient = new double [numberColumns+1];
2724                        memcpy(gradient,qp->objectiveAsObject()->gradient(qp,solution,offset,true,2),
2725                               numberColumns*sizeof(double));
2726                        double rhs = 0.0;
2727                        int * column = new int[numberColumns+1];
2728                        int n=0;
2729                        for (int i=0;i<numberColumns;i++) {
2730                          double value = gradient[i];
2731                          if (fabs(value)>1.0e-12) {
2732                            gradient[n]=value;
2733                            rhs += value*solution[i];
2734                            column[n++]=i;
2735                          }
2736                        }
2737                        gradient[n]=-1.0;
2738                        column[n++]=numberColumns;
2739                        storedAmpl.addCut(-COIN_DBL_MAX,offset+1.0e-7,n,column,gradient);
2740                        delete [] gradient;
2741                        delete [] column;
2742                      }
2743                      // could do three way branching round a) continuous b) best solution
2744                      printf("obj %g\n",bestObjectiveValue);
2745                      linkSolver->initialSolve();
2746                    }
2747                  }
2748                }
2749                si->setSpecialOptions(0x40000000);
2750              }
2751#endif
2752              if (!miplib) {
2753                if (!preSolve) {
2754                  model.solver()->setHintParam(OsiDoPresolveInInitial,false,OsiHintTry);
2755                  model.solver()->setHintParam(OsiDoPresolveInResolve,false,OsiHintTry);
2756                }
2757                model.initialSolve();
2758                OsiSolverInterface * solver = model.solver();
2759                OsiClpSolverInterface * si =
2760                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2761                ClpSimplex * clpSolver = si->getModelPtr();
2762                clpSolver->setSpecialOptions(clpSolver->specialOptions()|0x01000000); // say is Cbc (and in branch and bound)
2763                if (!complicatedInteger&&clpSolver->tightenPrimalBounds()!=0) {
2764                  std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
2765                  exit(1);
2766                }
2767                if (clpSolver->dualBound()==1.0e10) {
2768                  // user did not set - so modify
2769                  // get largest scaled away from bound
2770                  double largest=1.0e-12;
2771                  int numberRows = clpSolver->numberRows();
2772                  const double * rowPrimal = clpSolver->primalRowSolution();
2773                  const double * rowLower = clpSolver->rowLower();
2774                  const double * rowUpper = clpSolver->rowUpper();
2775                  const double * rowScale = clpSolver->rowScale();
2776                  int iRow;
2777                  for (iRow=0;iRow<numberRows;iRow++) {
2778                    double value = rowPrimal[iRow];
2779                    double above = value-rowLower[iRow];
2780                    double below = rowUpper[iRow]-value;
2781                    if (rowScale) {
2782                      double multiplier = rowScale[iRow];
2783                      above *= multiplier;
2784                      below *= multiplier;
2785                    }
2786                    if (above<1.0e12)
2787                      largest = CoinMax(largest,above);
2788                    if (below<1.0e12)
2789                      largest = CoinMax(largest,below);
2790                  }
2791                 
2792                  int numberColumns = clpSolver->numberColumns();
2793                  const double * columnPrimal = clpSolver->primalColumnSolution();
2794                  const double * columnLower = clpSolver->columnLower();
2795                  const double * columnUpper = clpSolver->columnUpper();
2796                  const double * columnScale = clpSolver->columnScale();
2797                  int iColumn;
2798                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2799                    double value = columnPrimal[iColumn];
2800                    double above = value-columnLower[iColumn];
2801                    double below = columnUpper[iColumn]-value;
2802                    if (columnScale) {
2803                      double multiplier = 1.0/columnScale[iColumn];
2804                      above *= multiplier;
2805                      below *= multiplier;
2806                    }
2807                    if (above<1.0e12)
2808                      largest = CoinMax(largest,above);
2809                    if (below<1.0e12)
2810                      largest = CoinMax(largest,below);
2811                  }
2812                  //if (!noPrinting)
2813                  //std::cout<<"Largest (scaled) away from bound "<<largest<<std::endl;
2814                  clpSolver->setDualBound(CoinMax(1.0001e8,CoinMin(1000.0*largest,1.00001e10)));
2815                }
2816                clpSolver->dual();  // clean up
2817              }
2818              // If user made settings then use them
2819              if (!defaultSettings) {
2820                OsiSolverInterface * solver = model.solver();
2821                if (!doScaling)
2822                  solver->setHintParam(OsiDoScale,false,OsiHintTry);
2823                OsiClpSolverInterface * si =
2824                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2825                assert (si != NULL);
2826                // get clp itself
2827                ClpSimplex * modelC = si->getModelPtr();
2828                //if (modelC->tightenPrimalBounds()!=0) {
2829                //std::cout<<"Problem is infeasible!"<<std::endl;
2830                //break;
2831                //}
2832                // bounds based on continuous
2833                if (tightenFactor&&!complicatedInteger) {
2834                  if (modelC->tightenPrimalBounds(tightenFactor)!=0) {
2835                    std::cout<<"Problem is infeasible!"<<std::endl;
2836                    break;
2837                  }
2838                }
2839                if (djFix<1.0e20) {
2840                  // do some fixing
2841                  int numberColumns = modelC->numberColumns();
2842                  int i;
2843                  const char * type = modelC->integerInformation();
2844                  double * lower = modelC->columnLower();
2845                  double * upper = modelC->columnUpper();
2846                  double * solution = modelC->primalColumnSolution();
2847                  double * dj = modelC->dualColumnSolution();
2848                  int numberFixed=0;
2849                  for (i=0;i<numberColumns;i++) {
2850                    if (type[i]) {
2851                      double value = solution[i];
2852                      if (value<lower[i]+1.0e-5&&dj[i]>djFix) {
2853                        solution[i]=lower[i];
2854                        upper[i]=lower[i];
2855                        numberFixed++;
2856                      } else if (value>upper[i]-1.0e-5&&dj[i]<-djFix) {
2857                        solution[i]=upper[i];
2858                        lower[i]=upper[i];
2859                        numberFixed++;
2860                      }
2861                    }
2862                  }
2863                  printf("%d columns fixed\n",numberFixed);
2864                }
2865              }
2866              // See if we want preprocessing
2867              OsiSolverInterface * saveSolver=NULL;
2868              CglPreProcess process;
2869              delete babModel;
2870              babModel = new CbcModel(model);
2871              OsiSolverInterface * solver3 = clpSolver->clone();
2872              babModel->assignSolver(solver3);
2873              OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2874              int numberChanged=0;
2875              if (clpSolver2->messageHandler()->logLevel())
2876                clpSolver2->messageHandler()->setLogLevel(1);
2877              if (logLevel>-1)
2878                clpSolver2->messageHandler()->setLogLevel(logLevel);
2879              lpSolver = clpSolver2->getModelPtr();
2880              if (lpSolver->factorizationFrequency()==200&&!miplib) {
2881                // User did not touch preset
2882                int numberRows = lpSolver->numberRows();
2883                const int cutoff1=10000;
2884                const int cutoff2=100000;
2885                const int base=75;
2886                const int freq0 = 50;
2887                const int freq1=200;
2888                const int freq2=400;
2889                const int maximum=1000;
2890                int frequency;
2891                if (numberRows<cutoff1)
2892                  frequency=base+numberRows/freq0;
2893                else if (numberRows<cutoff2)
2894                  frequency=base+cutoff1/freq0 + (numberRows-cutoff1)/freq1;
2895                else
2896                  frequency=base+cutoff1/freq0 + (cutoff2-cutoff1)/freq1 + (numberRows-cutoff2)/freq2;
2897                lpSolver->setFactorizationFrequency(CoinMin(maximum,frequency));
2898              }
2899              time2 = CoinCpuTime();
2900              totalTime += time2-time1;
2901              time1 = time2;
2902              double timeLeft = babModel->getMaximumSeconds();
2903              int numberOriginalColumns = babModel->solver()->getNumCols();
2904              if (preProcess==7) {
2905                // use strategy instead
2906                preProcess=0;
2907                useStrategy=true;
2908#ifdef COIN_HAS_LINK
2909                // empty out any cuts
2910                if (storedAmpl.sizeRowCuts()) {
2911                  printf("Emptying ampl stored cuts as internal preprocessing\n");
2912                  CglStored temp;
2913                  storedAmpl=temp;
2914                }
2915#endif
2916              }
2917              if (preProcess&&type==BAB) {
2918                // See if sos from mps file
2919                if (numberSOS==0&&clpSolver->numberSOS()&&doSOS) {
2920                  // SOS
2921                  numberSOS = clpSolver->numberSOS();
2922                  const CoinSet * setInfo = clpSolver->setInfo();
2923                  sosStart = new int [numberSOS+1];
2924                  sosType = new char [numberSOS];
2925                  int i;
2926                  int nTotal=0;
2927                  sosStart[0]=0;
2928                  for ( i=0;i<numberSOS;i++) {
2929                    int type = setInfo[i].setType();
2930                    int n=setInfo[i].numberEntries();
2931                    sosType[i]=type;
2932                    nTotal += n;
2933                    sosStart[i+1] = nTotal;
2934                  }
2935                  sosIndices = new int[nTotal];
2936                  sosReference = new double [nTotal];
2937                  for (i=0;i<numberSOS;i++) {
2938                    int n=setInfo[i].numberEntries();
2939                    const int * which = setInfo[i].which();
2940                    const double * weights = setInfo[i].weights();
2941                    int base = sosStart[i];
2942                    for (int j=0;j<n;j++) {
2943                      int k=which[j];
2944                      sosIndices[j+base]=k;
2945                      sosReference[j+base] = weights ? weights[j] : (double) j;
2946                    }
2947                  }
2948                }
2949                saveSolver=babModel->solver()->clone();
2950                /* Do not try and produce equality cliques and
2951                   do up to 10 passes */
2952                OsiSolverInterface * solver2;
2953                {
2954                  // Tell solver we are in Branch and Cut
2955                  saveSolver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo) ;
2956                  // Default set of cut generators
2957                  CglProbing generator1;
2958                  generator1.setUsingObjective(1);
2959                  generator1.setMaxPass(3);
2960                  generator1.setMaxProbeRoot(saveSolver->getNumCols());
2961                  generator1.setMaxElements(100);
2962                  generator1.setMaxLookRoot(50);
2963                  generator1.setRowCuts(3);
2964                  // Add in generators
2965                  process.addCutGenerator(&generator1);
2966                  int translate[]={9999,0,0,-1,2,3,-2};
2967                  process.messageHandler()->setLogLevel(babModel->logLevel());
2968#ifdef COIN_HAS_ASL
2969                  if (info.numberSos&&doSOS&&usingAmpl) {
2970                    // SOS
2971                    numberSOS = info.numberSos;
2972                    sosStart = info.sosStart;
2973                    sosIndices = info.sosIndices;
2974                  }
2975#endif
2976                  if (numberSOS&&doSOS) {
2977                    // SOS
2978                    int numberColumns = saveSolver->getNumCols();
2979                    char * prohibited = new char[numberColumns];
2980                    memset(prohibited,0,numberColumns);
2981                    int n=sosStart[numberSOS];
2982                    for (int i=0;i<n;i++) {
2983                      int iColumn = sosIndices[i];
2984                      prohibited[iColumn]=1;
2985                    }
2986                    process.passInProhibited(prohibited,numberColumns);
2987                    delete [] prohibited;
2988                  }
2989                  int numberPasses = 10;
2990                  if (tunePreProcess>=1000000) {
2991                    numberPasses = (tunePreProcess/1000000)-1;
2992                    tunePreProcess = tunePreProcess % 1000000;
2993                  } else if (tunePreProcess>=1000) {
2994                    numberPasses = (tunePreProcess/1000)-1;
2995                    tunePreProcess = tunePreProcess % 1000;
2996                  }
2997                  if (doSprint>0) {
2998                    // Sprint for primal solves
2999                    ClpSolve::SolveType method = ClpSolve::usePrimalorSprint;
3000                    ClpSolve::PresolveType presolveType = ClpSolve::presolveOff;
3001                    int numberPasses = 5;
3002                    int options[] = {0,3,0,0,0,0};
3003                    int extraInfo[] = {-1,20,-1,-1,-1,-1};
3004                    extraInfo[1]=doSprint;
3005                    int independentOptions[] = {0,0,3};
3006                    ClpSolve clpSolve(method,presolveType,numberPasses,
3007                                      options,extraInfo,independentOptions);
3008                    // say use in OsiClp
3009                    clpSolve.setSpecialOption(6,1);
3010                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (saveSolver);
3011                    osiclp->setSolveOptions(clpSolve);
3012                    osiclp->setHintParam(OsiDoDualInResolve,false);
3013                    // switch off row copy
3014                    osiclp->getModelPtr()->setSpecialOptions(osiclp->getModelPtr()->specialOptions()|256);
3015                    osiclp->getModelPtr()->setInfeasibilityCost(1.0e11);
3016                  }
3017                  solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],numberPasses,
3018                                                         tunePreProcess);
3019                  // Tell solver we are not in Branch and Cut
3020                  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
3021                  if (solver2)
3022                    solver2->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
3023                }
3024#ifdef COIN_HAS_ASL
3025                if (!solver2&&usingAmpl) {
3026                  // infeasible
3027                  info.problemStatus=1;
3028                  info.objValue = 1.0e100;
3029                  sprintf(info.buffer,"infeasible/unbounded by pre-processing");
3030                  info.primalSolution=NULL;
3031                  info.dualSolution=NULL;
3032                  break;
3033                }
3034#endif
3035                if (!noPrinting) {
3036                  if (!solver2) {
3037                    printf("Pre-processing says infeasible or unbounded\n");
3038                    break;
3039                  } else {
3040                    printf("processed model has %d rows, %d columns and %d elements\n",
3041                           solver2->getNumRows(),solver2->getNumCols(),solver2->getNumElements());
3042                  }
3043                }
3044                //solver2->resolve();
3045                if (preProcess==2) {
3046                  OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (solver2);
3047                  ClpSimplex * lpSolver = clpSolver2->getModelPtr();
3048                  lpSolver->writeMps("presolved.mps",0,1,lpSolver->optimizationDirection());
3049                  printf("Preprocessed model (minimization) on presolved.mps\n");
3050                }
3051                // we have to keep solver2 so pass clone
3052                solver2 = solver2->clone();
3053                babModel->assignSolver(solver2);
3054                babModel->initialSolve();
3055                babModel->setMaximumSeconds(timeLeft-(CoinCpuTime()-time1));
3056              }
3057              // now tighten bounds
3058              if (!miplib) {
3059                OsiClpSolverInterface * si =
3060                  dynamic_cast<OsiClpSolverInterface *>(babModel->solver()) ;
3061                assert (si != NULL);
3062                // get clp itself
3063                ClpSimplex * modelC = si->getModelPtr();
3064                if (noPrinting)
3065                  modelC->setLogLevel(0);
3066                if (!complicatedInteger&&modelC->tightenPrimalBounds()!=0) {
3067                  std::cout<<"Problem is infeasible!"<<std::endl;
3068                  break;
3069                }
3070                modelC->dual();
3071              }
3072#if 0
3073              numberDebugValues=599;
3074              debugValues = new double[numberDebugValues];
3075              CoinZeroN(debugValues,numberDebugValues);
3076              debugValues[3]=1.0;
3077              debugValues[6]=25.0;
3078              debugValues[9]=4.0;
3079              debugValues[26]=4.0;
3080              debugValues[27]=6.0;
3081              debugValues[35]=8.0;
3082              debugValues[53]=21.0;
3083              debugValues[56]=4.0;
3084#endif
3085              if (debugValues) {
3086                // for debug
3087                std::string problemName ;
3088                babModel->solver()->getStrParam(OsiProbName,problemName) ;
3089                babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
3090                twomirGen.probname_=strdup(problemName.c_str());
3091                // checking seems odd
3092                //redsplitGen.set_given_optsol(babModel->solver()->getRowCutDebuggerAlways()->optimalSolution(),
3093                //                         babModel->getNumCols());
3094              }
3095              int testOsiOptions = parameters[whichParam(TESTOSI,numberParameters,parameters)].intValue();
3096#ifdef COIN_HAS_LINK
3097              // If linked then see if expansion wanted
3098              {
3099                OsiSolverLink * solver3 = dynamic_cast<OsiSolverLink *> (babModel->solver());
3100                if (solver3) {
3101                  int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000;
3102                  if (options) {
3103                    /*
3104                      1 - force mini branch and bound
3105                      2 - set priorities high on continuous
3106                      4 - try adding OA cuts
3107                      8 - try doing quadratic linearization
3108                      16 - try expanding knapsacks
3109                    */
3110                    if ((options&16)) {
3111                      int numberColumns = saveCoinModel.numberColumns();
3112                      int numberRows = saveCoinModel.numberRows();
3113                      whichColumn = new int[numberColumns];
3114                      knapsackStart=new int[numberRows+1];
3115                      knapsackRow=new int[numberRows];
3116                      numberKnapsack=10000;
3117                      int extra1 = parameters[whichParam(EXTRA1,numberParameters,parameters)].intValue();
3118                      int extra2 = parameters[whichParam(EXTRA2,numberParameters,parameters)].intValue();
3119                      int logLevel = parameters[log].intValue();
3120                      OsiSolverInterface * solver = expandKnapsack(saveCoinModel,whichColumn,knapsackStart,
3121                                                                   knapsackRow,numberKnapsack,
3122                                                                   storedAmpl,logLevel,extra1,extra2);
3123                      if (solver) {
3124                        clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
3125                        assert (clpSolver);
3126                        lpSolver = clpSolver->getModelPtr();
3127                        babModel->assignSolver(solver);
3128                        testOsiOptions=0;
3129                        // Priorities already done
3130                        free(info.priorities);
3131                        info.priorities=NULL;
3132                      } else {
3133                        numberKnapsack=0;
3134                        delete [] whichColumn;
3135                        delete [] knapsackStart;
3136                        delete [] knapsackRow;
3137                        whichColumn=NULL;
3138                        knapsackStart=NULL;
3139                        knapsackRow=NULL;
3140                      }
3141                    }
3142                  }
3143                }
3144              }
3145#endif
3146              if (useCosts&&testOsiOptions<0) {
3147                int numberColumns = babModel->getNumCols();
3148                int * sort = new int[numberColumns];
3149                double * dsort = new double[numberColumns];
3150                int * priority = new int [numberColumns];
3151                const double * objective = babModel->getObjCoefficients();
3152                const double * lower = babModel->getColLower() ;
3153                const double * upper = babModel->getColUpper() ;
3154                const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
3155                const int * columnLength = matrix->getVectorLengths();
3156                int iColumn;
3157                int n=0;
3158                for (iColumn=0;iColumn<numberColumns;iColumn++) {
3159                  if (babModel->isInteger(iColumn)) {
3160                    sort[n]=n;
3161                    if (useCosts==1)
3162                      dsort[n++]=-fabs(objective[iColumn]);
3163                    else if (useCosts==2)
3164                      dsort[n++]=iColumn;
3165                    else if (useCosts==3)
3166                      dsort[n++]=upper[iColumn]-lower[iColumn];
3167                    else if (useCosts==4)
3168                      dsort[n++]=-(upper[iColumn]-lower[iColumn]);
3169                    else if (useCosts==5)
3170                      dsort[n++]=-columnLength[iColumn];
3171                  }
3172                }
3173                CoinSort_2(dsort,dsort+n,sort);
3174                int level=0;
3175                double last = -1.0e100;
3176                for (int i=0;i<n;i++) {
3177                  int iPut=sort[i];
3178                  if (dsort[i]!=last) {
3179                    level++;
3180                    last=dsort[i];
3181                  }
3182                  priority[iPut]=level;
3183                }
3184                babModel->passInPriorities( priority,false);
3185                delete [] priority;
3186                delete [] sort;
3187                delete [] dsort;
3188              }
3189              // FPump done first as it only works if no solution
3190              CbcHeuristicFPump heuristic4(*babModel);
3191              if (useFpump) {
3192                heuristic4.setMaximumPasses(parameters[whichParam(FPUMPITS,numberParameters,parameters)].intValue());
3193                int pumpTune=parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].intValue();
3194                if (pumpTune>0) {
3195                  /*
3196                    >=10000000 for using obj
3197                    >=1000000 use as accumulate switch
3198                    >=1000 use index+1 as number of large loops
3199                    >=100 use 0.05 objvalue as increment
3200                    >=10 use +0.1 objvalue for cutoff (add)
3201                    1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
3202                    4 and static continuous, 5 as 3 but no internal integers
3203                    6 as 3 but all slack basis!
3204                  */
3205                  int logLevel = parameters[log].intValue();
3206                  double value = babModel->solver()->getObjSense()*babModel->solver()->getObjValue();
3207                  int w = pumpTune/10;
3208                  int c = w % 10;
3209                  w /= 10;
3210                  int i = w % 10;
3211                  w /= 10;
3212                  int r = w;
3213                  int accumulate = r/1000;
3214                  r -= 1000*accumulate;
3215                  if (accumulate>=10) {
3216                    int which = accumulate/10;
3217                    accumulate -= 10*which;
3218                    which--;
3219                    // weights and factors
3220                    double weight[]={0.1,0.1,0.5,0.5,1.0,1.0,5.0,5.0};
3221                    double factor[] = {0.1,0.5,0.1,0.5,0.1,0.5,0.1,0.5};
3222                    heuristic4.setInitialWeight(weight[which]);
3223                    heuristic4.setWeightFactor(factor[which]);
3224                  }
3225                  // fake cutoff
3226                  if (logLevel>1)
3227                    printf("Setting ");
3228                  if (c) {
3229                    double cutoff;
3230                    babModel->solver()->getDblParam(OsiDualObjectiveLimit,cutoff);
3231                    cutoff = CoinMin(cutoff,value + 0.1*fabs(value)*c);
3232                    heuristic4.setFakeCutoff(cutoff);
3233                    if (logLevel>1)
3234                      printf("fake cutoff of %g ",cutoff);
3235                  }
3236                  if (i||r) {
3237                    // also set increment
3238                    heuristic4.setAbsoluteIncrement((0.01*i+0.005)*(fabs(value)+1.0e-12));
3239                    heuristic4.setAccumulate(accumulate);
3240                    heuristic4.setMaximumRetries(r+1);
3241                    if (logLevel>1) {
3242                      if (i) 
3243                        printf("increment of %g ",heuristic4.absoluteIncrement());
3244                      if (accumulate) 
3245                        printf("accumulate of %d ",accumulate);
3246                      printf("%d retries ",r+2);
3247                    }
3248                  }
3249                  pumpTune = pumpTune%100;
3250                  if (logLevel>1)
3251                    printf("and setting when to %d\n",pumpTune+10);
3252                  if (pumpTune==6)
3253                    pumpTune =13;
3254                  heuristic4.setWhen(pumpTune+10);
3255                }
3256                babModel->addHeuristic(&heuristic4);
3257              }
3258              if (!miplib) {
3259                CbcRounding heuristic1(*babModel);
3260                if (useRounding)
3261                  babModel->addHeuristic(&heuristic1) ;
3262                CbcHeuristicLocal heuristic2(*babModel);
3263                heuristic2.setSearchType(1);
3264                if (useCombine)
3265                  babModel->addHeuristic(&heuristic2);
3266                CbcHeuristicGreedyCover heuristic3(*babModel);
3267                CbcHeuristicGreedyEquality heuristic3a(*babModel);
3268                if (useGreedy) {
3269                  babModel->addHeuristic(&heuristic3);
3270                  babModel->addHeuristic(&heuristic3a);
3271                }
3272                if (useLocalTree) {
3273                  CbcTreeLocal localTree(babModel,NULL,10,0,0,10000,2000);
3274                  babModel->passInTreeHandler(localTree);
3275                }
3276              }
3277              CbcHeuristicRINS heuristic5(*babModel);
3278              if (useRINS)
3279                babModel->addHeuristic(&heuristic5) ;
3280              if (type==MIPLIB) {
3281                if (babModel->numberStrong()==5&&babModel->numberBeforeTrust()==5) 
3282                  babModel->setNumberBeforeTrust(50);
3283              }
3284              // add cut generators if wanted
3285              int switches[20];
3286              int numberGenerators=0;
3287              int translate[]={-100,-1,-99,-98,1,1,1,1};
3288              if (probingAction) {
3289                if (probingAction==5||probingAction==7)
3290                  probingGen.setRowCuts(-3); // strengthening etc just at root
3291                if (probingAction==6||probingAction==7) {
3292                  // Number of unsatisfied variables to look at
3293                  probingGen.setMaxProbe(1000);
3294                  probingGen.setMaxProbeRoot(1000);
3295                  // How far to follow the consequences
3296                  probingGen.setMaxLook(50);
3297                  probingGen.setMaxLookRoot(50);
3298                }
3299                babModel->addCutGenerator(&probingGen,translate[probingAction],"Probing");
3300                switches[numberGenerators++]=0;
3301              }
3302              if (gomoryAction&&(complicatedInteger!=1||gomoryAction==1)) {
3303                babModel->addCutGenerator(&gomoryGen,translate[gomoryAction],"Gomory");
3304                switches[numberGenerators++]=-1;
3305              }
3306              if (knapsackAction) {
3307                babModel->addCutGenerator(&knapsackGen,translate[knapsackAction],"Knapsack");
3308                switches[numberGenerators++]=0;
3309              }
3310              if (redsplitAction&&!complicatedInteger) {
3311                babModel->addCutGenerator(&redsplitGen,translate[redsplitAction],"Reduce-and-split");
3312                switches[numberGenerators++]=1;
3313              }
3314              if (cliqueAction) {
3315                babModel->addCutGenerator(&cliqueGen,translate[cliqueAction],"Clique");
3316                switches[numberGenerators++]=0;
3317              }
3318              if (mixedAction) {
3319                babModel->addCutGenerator(&mixedGen,translate[mixedAction],"MixedIntegerRounding2");
3320                switches[numberGenerators++]=-1;
3321              }
3322              if (flowAction) {
3323                babModel->addCutGenerator(&flowGen,translate[flowAction],"FlowCover");
3324                switches[numberGenerators++]=1;
3325              }
3326              if (twomirAction&&!complicatedInteger) {
3327                babModel->addCutGenerator(&twomirGen,translate[twomirAction],"TwoMirCuts");
3328                switches[numberGenerators++]=1;
3329              }
3330              if (landpAction) {
3331                babModel->addCutGenerator(&landpGen,translate[landpAction],"LiftAndProject");
3332                switches[numberGenerators++]=1;
3333              }
3334              if (storedCuts) 
3335                babModel->setSpecialOptions(babModel->specialOptions()|64);
3336              // Say we want timings
3337              numberGenerators = babModel->numberCutGenerators();
3338              int iGenerator;
3339              int cutDepth=
3340                parameters[whichParam(CUTDEPTH,numberParameters,parameters)].intValue();
3341              for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
3342                CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
3343                int howOften = generator->howOften();
3344                if (howOften==-98||howOften==-99) 
3345                  generator->setSwitchOffIfLessThan(switches[iGenerator]);
3346                generator->setTiming(true);
3347                if (cutDepth>=0)
3348                  generator->setWhatDepth(cutDepth) ;
3349              }
3350              // Could tune more
3351              if (!miplib) {
3352                babModel->setMinimumDrop(min(5.0e-2,
3353                                             fabs(babModel->getMinimizationObjValue())*1.0e-3+1.0e-4));
3354                if (cutPass==-1234567) {
3355                  if (babModel->getNumCols()<500)
3356                    babModel->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
3357                  else if (babModel->getNumCols()<5000)
3358                    babModel->setMaximumCutPassesAtRoot(100); // use minimum drop
3359                  else
3360                    babModel->setMaximumCutPassesAtRoot(20);
3361                } else {
3362                  babModel->setMaximumCutPassesAtRoot(cutPass);
3363                }
3364                babModel->setMaximumCutPasses(1);
3365              }
3366              // Do more strong branching if small
3367              //if (babModel->getNumCols()<5000)
3368              //babModel->setNumberStrong(20);
3369              // Switch off strong branching if wanted
3370              //if (babModel->getNumCols()>10*babModel->getNumRows())
3371              //babModel->setNumberStrong(0);
3372              if (!noPrinting) {
3373                int iLevel = parameters[log].intValue();
3374                if (iLevel<0) {
3375                  babModel->setPrintingMode(1);
3376                  iLevel = -iLevel;
3377                }
3378                babModel->messageHandler()->setLogLevel(iLevel);
3379                if (babModel->getNumCols()>2000||babModel->getNumRows()>1500||
3380                    babModel->messageHandler()->logLevel()>1)
3381                  babModel->setPrintFrequency(100);
3382              }
3383             
3384              babModel->solver()->setIntParam(OsiMaxNumIterationHotStart,
3385                    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].intValue());
3386              OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3387              // go faster stripes
3388              if (osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
3389                osiclp->setupForRepeatedUse(2,parameters[slog].intValue());
3390              } else {
3391                osiclp->setupForRepeatedUse(0,parameters[slog].intValue());
3392              }
3393              double increment=babModel->getCutoffIncrement();;
3394              int * changed = NULL;
3395              if (!miplib&&increment==normalIncrement)
3396                changed=analyze( osiclp,numberChanged,increment,false);
3397              if (debugValues) {
3398                if (numberDebugValues==babModel->getNumCols()) {
3399                  // for debug
3400                  babModel->solver()->activateRowCutDebugger(debugValues) ;
3401                } else {
3402                  printf("debug file has incorrect number of columns\n");
3403                }
3404              }
3405              babModel->setCutoffIncrement(CoinMax(babModel->getCutoffIncrement(),increment));
3406              // Turn this off if you get problems
3407              // Used to be automatically set
3408              int mipOptions = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()%10000;
3409              if (mipOptions!=(1))
3410                printf("mip options %d\n",mipOptions);
3411              osiclp->setSpecialOptions(mipOptions);
3412              if (gapRatio < 1.0e100) {
3413                double value = babModel->solver()->getObjValue() ;
3414                double value2 = gapRatio*(1.0e-5+fabs(value)) ;
3415                babModel->setAllowableGap(value2) ;
3416                std::cout << "Continuous " << value
3417                          << ", so allowable gap set to "
3418                          << value2 << std::endl ;
3419              }
3420              // probably faster to use a basis to get integer solutions
3421              babModel->setSpecialOptions(babModel->specialOptions()|2);
3422              currentBranchModel = babModel;
3423              OsiSolverInterface * strengthenedModel=NULL;
3424              if (type==BAB||type==MIPLIB) {
3425                int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
3426                if (moreMipOptions>=0) {
3427                  printf("more mip options %d\n",moreMipOptions);
3428                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3429                  if (moreMipOptions==10000) {
3430                    // test memory saving
3431                    moreMipOptions -= 10000;
3432                    ClpSimplex * lpSolver = osiclp->getModelPtr();
3433                    lpSolver->setPersistenceFlag(1);
3434                    // switch off row copy if few rows
3435                    if (lpSolver->numberRows()<150)
3436                      lpSolver->setSpecialOptions(lpSolver->specialOptions()|256);
3437                  }
3438                  if (((moreMipOptions+1)%1000000)!=0)
3439                    babModel->setSearchStrategy(moreMipOptions%1000000);
3440                  // go faster stripes
3441                  if( moreMipOptions >=999999) {
3442                    if (osiclp) {
3443                      int save = osiclp->specialOptions();
3444                      osiclp->setupForRepeatedUse(2,0);
3445                      osiclp->setSpecialOptions(save|osiclp->specialOptions());
3446                    }
3447                  } 
3448                }
3449              }
3450              if (type==BAB) {
3451#ifdef COIN_HAS_ASL
3452                if (usingAmpl) {
3453                  priorities=info.priorities;
3454                  branchDirection=info.branchDirection;
3455                  pseudoDown=info.pseudoDown;
3456                  pseudoUp=info.pseudoUp;
3457                  solutionIn=info.primalSolution;
3458                  prioritiesIn = info.priorities;
3459                  if (info.numberSos&&doSOS) {
3460                    // SOS
3461                    numberSOS = info.numberSos;
3462                    sosStart = info.sosStart;
3463                    sosIndices = info.sosIndices;
3464                    sosType = info.sosType;
3465                    sosReference = info.sosReference;
3466                    sosPriority = info.sosPriority;
3467                  }
3468                }
3469#endif               
3470                const int * originalColumns = preProcess ? process.originalColumns() : NULL;
3471                if (solutionIn&&useSolution) {
3472                  if (preProcess) {
3473                    int numberColumns = babModel->getNumCols();
3474                    // extend arrays in case SOS
3475                    int n = originalColumns[numberColumns-1]+1;
3476                    int nSmaller = CoinMin(n,numberOriginalColumns);
3477                    double * solutionIn2 = new double [n];
3478                    int * prioritiesIn2 = new int[n];
3479                    int i;
3480                    for (i=0;i<nSmaller;i++) {
3481                      solutionIn2[i]=solutionIn[i];
3482                      prioritiesIn2[i]=prioritiesIn[i];
3483                    }
3484                    for (;i<n;i++) {
3485                      solutionIn2[i]=0.0;
3486                      prioritiesIn2[i]=1000000;
3487                    }
3488                    int iLast=-1;
3489                    for (i=0;i<numberColumns;i++) {
3490                      int iColumn = originalColumns[i];
3491                      assert (iColumn>iLast);
3492                      iLast=iColumn;
3493                      solutionIn2[i]=solutionIn2[iColumn];
3494                      if (prioritiesIn)
3495                        prioritiesIn2[i]=prioritiesIn2[iColumn];
3496                    }
3497                    babModel->setHotstartSolution(solutionIn2,prioritiesIn2);
3498                    delete [] solutionIn2;
3499                    delete [] prioritiesIn2;
3500                  } else {
3501                    babModel->setHotstartSolution(solutionIn,prioritiesIn);
3502                  }
3503                }
3504                OsiSolverInterface * testOsiSolver= (testOsiOptions>=0) ? babModel->solver() : NULL;
3505                if (!testOsiSolver) {
3506                  // *************************************************************
3507                  // CbcObjects
3508                  if (preProcess&&process.numberSOS()) {
3509                    int numberSOS = process.numberSOS();
3510                    int numberIntegers = babModel->numberIntegers();
3511                    /* model may not have created objects
3512                       If none then create
3513                    */
3514                    if (!numberIntegers||!babModel->numberObjects()) {
3515                      int type = (pseudoUp) ? 1 : 0;
3516                      babModel->findIntegers(true,type);
3517                      numberIntegers = babModel->numberIntegers();
3518                    }
3519                    OsiObject ** oldObjects = babModel->objects();
3520                    // Do sets and priorities
3521                    OsiObject ** objects = new OsiObject * [numberSOS];
3522                    // set old objects to have low priority
3523                    int numberOldObjects = babModel->numberObjects();
3524                    int numberColumns = babModel->getNumCols();
3525                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
3526                      oldObjects[iObj]->setPriority(numberColumns+1);
3527                      int iColumn = oldObjects[iObj]->columnNumber();
3528                      assert (iColumn>=0);
3529                      if (iColumn>=numberOriginalColumns)
3530                        continue;
3531                      if (originalColumns)
3532                        iColumn = originalColumns[iColumn];
3533                      if (branchDirection) {
3534                        CbcSimpleInteger * obj =
3535                          dynamic_cast <CbcSimpleInteger *>(oldObjects[iObj]) ;
3536                        if (obj) { 
3537                          obj->setPreferredWay(branchDirection[iColumn]);
3538                        } else {
3539                          CbcObject * obj =
3540                            dynamic_cast <CbcObject *>(oldObjects[iObj]) ;
3541                          assert (obj);
3542                          obj->setPreferredWay(branchDirection[iColumn]);
3543                        }
3544                      }
3545                      if (pseudoUp) {
3546                        CbcSimpleIntegerPseudoCost * obj1a =
3547                          dynamic_cast <CbcSimpleIntegerPseudoCost *>(oldObjects[iObj]) ;
3548                        assert (obj1a);
3549                        if (pseudoDown[iColumn]>0.0)
3550                          obj1a->setDownPseudoCost(pseudoDown[iColumn]);
3551                        if (pseudoUp[iColumn]>0.0)
3552                          obj1a->setUpPseudoCost(pseudoUp[iColumn]);
3553                      }
3554                    }
3555                    const int * starts = process.startSOS();
3556                    const int * which = process.whichSOS();
3557                    const int * type = process.typeSOS();
3558                    const double * weight = process.weightSOS();
3559                    int iSOS;
3560                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
3561                      int iStart = starts[iSOS];
3562                      int n=starts[iSOS+1]-iStart;
3563                      objects[iSOS] = new CbcSOS(babModel,n,which+iStart,weight+iStart,
3564                                                 iSOS,type[iSOS]);
3565                      // branch on long sets first
3566                      objects[iSOS]->setPriority(numberColumns-n);
3567                    }
3568                    babModel->addObjects(numberSOS,objects);
3569                    for (iSOS=0;iSOS<numberSOS;iSOS++)
3570                      delete objects[iSOS];
3571                    delete [] objects;
3572                  } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
3573                    // do anyway for priorities etc
3574                    int numberIntegers = babModel->numberIntegers();
3575                    /* model may not have created objects
3576                       If none then create
3577                    */
3578                    if (!numberIntegers||!babModel->numberObjects()) {
3579                      int type = (pseudoUp) ? 1 : 0;
3580                      babModel->findIntegers(true,type);
3581                    }
3582                    if (numberSOS) {
3583                      // Do sets and priorities
3584                      OsiObject ** objects = new OsiObject * [numberSOS];
3585                      int iSOS;
3586                      if (originalColumns) {
3587                        // redo sequence numbers
3588                        int numberColumns = babModel->getNumCols();
3589                        int nOld = originalColumns[numberColumns-1]+1;
3590                        int * back = new int[nOld];
3591                        int i;
3592                        for (i=0;i<nOld;i++)
3593                          back[i]=-1;
3594                        for (i=0;i<numberColumns;i++)
3595                          back[originalColumns[i]]=i;
3596                        // Really need better checks
3597                        int nMissing=0;
3598                        int n=sosStart[numberSOS];
3599                        for (i=0;i<n;i++) {
3600                          int iColumn = sosIndices[i];
3601                          int jColumn = back[iColumn];
3602                          if (jColumn>=0) 
3603                            sosIndices[i] = jColumn;
3604                          else 
3605                            nMissing++;
3606                        }
3607                        delete [] back;
3608                        if (nMissing)
3609                          printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
3610                      }
3611                      for (iSOS =0;iSOS<numberSOS;iSOS++) {
3612                        int iStart = sosStart[iSOS];
3613                        int n=sosStart[iSOS+1]-iStart;
3614                        objects[iSOS] = new CbcSOS(babModel,n,sosIndices+iStart,sosReference+iStart,
3615                                                   iSOS,sosType[iSOS]);
3616                        if (sosPriority)
3617                          objects[iSOS]->setPriority(sosPriority[iSOS]);
3618                        else if (!prioritiesIn)
3619                          objects[iSOS]->setPriority(10);  // rather than 1000
3620                      }
3621                      // delete any existing SOS objects
3622                      int numberObjects=babModel->numberObjects();
3623                      OsiObject ** oldObjects=babModel->objects();
3624                      int nNew=0;
3625                      for (int i=0;i<numberObjects;i++) {
3626                        OsiObject * objThis = oldObjects[i];
3627                        CbcSOS * obj1 =
3628                          dynamic_cast <CbcSOS *>(objThis) ;
3629                        OsiSOS * obj2 =
3630                          dynamic_cast <OsiSOS *>(objThis) ;
3631                        if (!obj1&&!obj2) {
3632                          oldObjects[nNew++]=objThis;
3633                        } else {
3634                          delete objThis;
3635                        }
3636                      }
3637                      babModel->setNumberObjects(nNew);
3638                      babModel->addObjects(numberSOS,objects);
3639                      for (iSOS=0;iSOS<numberSOS;iSOS++)
3640                        delete objects[iSOS];
3641                      delete [] objects;
3642                    }
3643                  }
3644                  OsiObject ** objects = babModel->objects();
3645                  int numberObjects = babModel->numberObjects();
3646                  for (int iObj = 0;iObj<numberObjects;iObj++) {
3647                    // skip sos
3648                    CbcSOS * objSOS =
3649                      dynamic_cast <CbcSOS *>(objects[iObj]) ;
3650                    if (objSOS)
3651                      continue;
3652                    int iColumn = objects[iObj]->columnNumber();
3653                    assert (iColumn>=0);
3654                    if (originalColumns)
3655                      iColumn = originalColumns[iColumn];
3656                    if (branchDirection) {
3657                      CbcSimpleInteger * obj =
3658                        dynamic_cast <CbcSimpleInteger *>(objects[iObj]) ;
3659                      if (obj) { 
3660                        obj->setPreferredWay(branchDirection[iColumn]);
3661                      } else {
3662                        CbcObject * obj =
3663                          dynamic_cast <CbcObject *>(objects[iObj]) ;
3664                        assert (obj);
3665                        obj->setPreferredWay(branchDirection[iColumn]);
3666                      }
3667                    }
3668                    if (priorities) {
3669                      int iPriority = priorities[iColumn];
3670                      if (iPriority>0)
3671                        objects[iObj]->setPriority(iPriority);
3672                    }
3673                    if (pseudoUp&&pseudoUp[iColumn]) {
3674                      CbcSimpleIntegerPseudoCost * obj1a =
3675                        dynamic_cast <CbcSimpleIntegerPseudoCost *>(objects[iObj]) ;
3676                      assert (obj1a);
3677                      if (pseudoDown[iColumn]>0.0)
3678                        obj1a->setDownPseudoCost(pseudoDown[iColumn]);
3679                      if (pseudoUp[iColumn]>0.0)
3680                        obj1a->setUpPseudoCost(pseudoUp[iColumn]);
3681                    }
3682                  }
3683                  // *************************************************************
3684                } else {
3685                  // *************************************************************
3686                  // OsiObjects
3687                  // Find if none
3688                  int numberIntegers = testOsiSolver->getNumIntegers();
3689                  /* model may not have created objects
3690                     If none then create
3691                  */
3692                  if (!numberIntegers||!testOsiSolver->numberObjects()) {
3693                    //int type = (pseudoUp) ? 1 : 0;
3694                    testOsiSolver->findIntegers(false);
3695                    numberIntegers = testOsiSolver->getNumIntegers();
3696                  }
3697                  if (preProcess&&process.numberSOS()) {
3698                    int numberSOS = process.numberSOS();
3699                    OsiObject ** oldObjects = testOsiSolver->objects();
3700                    // Do sets and priorities
3701                    OsiObject ** objects = new OsiObject * [numberSOS];
3702                    // set old objects to have low priority
3703                    int numberOldObjects = testOsiSolver->numberObjects();
3704                    int numberColumns = testOsiSolver->getNumCols();
3705                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
3706                      oldObjects[iObj]->setPriority(numberColumns+1);
3707                      int iColumn = oldObjects[iObj]->columnNumber();
3708                      assert (iColumn>=0);
3709                      if (iColumn>=numberOriginalColumns)
3710                        continue;
3711                      if (originalColumns)
3712                        iColumn = originalColumns[iColumn];
3713                      if (branchDirection) {
3714                        OsiSimpleInteger * obj =
3715                          dynamic_cast <OsiSimpleInteger *>(oldObjects[iObj]) ;
3716                        if (obj) { 
3717                          obj->setPreferredWay(branchDirection[iColumn]);
3718                        } else {
3719                          OsiObject2 * obj =
3720                            dynamic_cast <OsiObject2 *>(oldObjects[iObj]) ;
3721                          if (obj)
3722                            obj->setPreferredWay(branchDirection[iColumn]);
3723                        }
3724                      }
3725                      if (pseudoUp) {
3726                        abort();
3727                      }
3728                    }
3729                    const int * starts = process.startSOS();
3730                    const int * which = process.whichSOS();
3731                    const int * type = process.typeSOS();
3732                    const double * weight = process.weightSOS();
3733                    int iSOS;
3734                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
3735                      int iStart = starts[iSOS];
3736                      int n=starts[iSOS+1]-iStart;
3737                      objects[iSOS] = new OsiSOS(testOsiSolver,n,which+iStart,weight+iStart,
3738                                                 type[iSOS]);
3739                      // branch on long sets first
3740                      objects[iSOS]->setPriority(numberColumns-n);
3741                    }
3742                    testOsiSolver->addObjects(numberSOS,objects);
3743                    for (iSOS=0;iSOS<numberSOS;iSOS++)
3744                      delete objects[iSOS];
3745                    delete [] objects;
3746                  } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
3747                    if (numberSOS) {
3748                      // Do sets and priorities
3749                      OsiObject ** objects = new OsiObject * [numberSOS];
3750                      int iSOS;
3751                      if (originalColumns) {
3752                        // redo sequence numbers
3753                        int numberColumns = testOsiSolver->getNumCols();
3754                        int nOld = originalColumns[numberColumns-1]+1;
3755                        int * back = new int[nOld];
3756                        int i;
3757                        for (i=0;i<nOld;i++)
3758                          back[i]=-1;
3759                        for (i=0;i<numberColumns;i++)
3760                          back[originalColumns[i]]=i;
3761                        // Really need better checks
3762                        int nMissing=0;
3763                        int n=sosStart[numberSOS];
3764                        for (i=0;i<n;i++) {
3765                          int iColumn = sosIndices[i];
3766                          int jColumn = back[iColumn];
3767                          if (jColumn>=0) 
3768                            sosIndices[i] = jColumn;
3769                          else 
3770                            nMissing++;
3771                        }
3772                        delete [] back;
3773                        if (nMissing)
3774                          printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
3775                      }
3776                      for (iSOS =0;iSOS<numberSOS;iSOS++) {
3777                        int iStart = sosStart[iSOS];
3778                        int n=sosStart[iSOS+1]-iStart;
3779                        objects[iSOS] = new OsiSOS(testOsiSolver,n,sosIndices+iStart,sosReference+iStart,
3780                                                   sosType[iSOS]);
3781                        if (sosPriority)
3782                          objects[iSOS]->setPriority(sosPriority[iSOS]);
3783                        else if (!prioritiesIn)
3784                          objects[iSOS]->setPriority(10);  // rather than 1000
3785                      }
3786                      // delete any existing SOS objects
3787                      int numberObjects=testOsiSolver->numberObjects();
3788                      OsiObject ** oldObjects=testOsiSolver->objects();
3789                      int nNew=0;
3790                      for (int i=0;i<numberObjects;i++) {
3791                        OsiObject * objThis = oldObjects[i];
3792                        OsiSOS * obj1 =
3793                          dynamic_cast <OsiSOS *>(objThis) ;
3794                        OsiSOS * obj2 =
3795                          dynamic_cast <OsiSOS *>(objThis) ;
3796                        if (!obj1&&!obj2) {
3797                          oldObjects[nNew++]=objThis;
3798                        } else {
3799                          delete objThis;
3800                        }
3801                      }
3802                      testOsiSolver->setNumberObjects(nNew);
3803                      testOsiSolver->addObjects(numberSOS,objects);
3804                      for (iSOS=0;iSOS<numberSOS;iSOS++)
3805                        delete objects[iSOS];
3806                      delete [] objects;
3807                    }
3808                  }
3809                  OsiObject ** objects = testOsiSolver->objects();
3810                  int numberObjects = testOsiSolver->numberObjects();
3811                  int logLevel = parameters[log].intValue();
3812                  for (int iObj = 0;iObj<numberObjects;iObj++) {
3813                    // skip sos
3814                    OsiSOS * objSOS =
3815                      dynamic_cast <OsiSOS *>(objects[iObj]) ;
3816                    if (objSOS) {
3817                      if (logLevel>2)
3818                        printf("Set %d is SOS - priority %d\n",iObj,objSOS->priority());
3819                      continue;
3820                    }
3821                    int iColumn = objects[iObj]->columnNumber();
3822                    if (iColumn>=0) {
3823                      if (originalColumns)
3824                        iColumn = originalColumns[iColumn];
3825                      if (branchDirection) {
3826                        OsiSimpleInteger * obj =
3827                          dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
3828                        if (obj) { 
3829                          obj->setPreferredWay(branchDirection[iColumn]);
3830                        } else {
3831                          OsiObject2 * obj =
3832                            dynamic_cast <OsiObject2 *>(objects[iObj]) ;
3833                          if (obj)
3834                            obj->setPreferredWay(branchDirection[iColumn]);
3835                        }
3836                      }
3837                    }
3838                    if (priorities) {
3839                      int iPriority = priorities[iColumn];
3840                      if (iPriority>0)
3841                        objects[iObj]->setPriority(iPriority);
3842                    }
3843                    if (logLevel>2)
3844                      printf("Obj %d is int? - priority %d\n",iObj,objects[iObj]->priority());
3845                    if (pseudoUp&&pseudoUp[iColumn]) {
3846                      abort();
3847                    }
3848                  }
3849                  // *************************************************************
3850                }
3851                int statistics = (printOptions>0) ? printOptions: 0;
3852#ifdef COIN_HAS_ASL
3853                if (!usingAmpl) {
3854#endif
3855                  free(priorities);
3856                  priorities=NULL;
3857                  free(branchDirection);
3858                  branchDirection=NULL;
3859                  free(pseudoDown);
3860                  pseudoDown=NULL;
3861                  free(pseudoUp);
3862                  pseudoUp=NULL;
3863                  free(solutionIn);
3864                  solutionIn=NULL;
3865                  free(prioritiesIn);
3866                  prioritiesIn=NULL;
3867                  free(sosStart);
3868                  sosStart=NULL;
3869                  free(sosIndices);
3870                  sosIndices=NULL;
3871                  free(sosType);
3872                  sosType=NULL;
3873                  free(sosReference);
3874                  sosReference=NULL;
3875                  free(cut);
3876                  cut=NULL;
3877                  free(sosPriority);
3878                  sosPriority=NULL;
3879#ifdef COIN_HAS_ASL
3880                }
3881#endif               
3882                if (nodeStrategy) {
3883                  // change default
3884                  if (nodeStrategy>1) {
3885                    // up or down
3886                    int way = ((nodeStrategy%1)==1) ? -1 : +1;
3887                    babModel->setPreferredWay(way);
3888#if 0
3889                    OsiObject ** objects = babModel->objects();
3890                    int numberObjects = babModel->numberObjects();
3891                    for (int iObj = 0;iObj<numberObjects;iObj++) {
3892                      CbcObject * obj =
3893                        dynamic_cast <CbcObject *>(objects[iObj]) ;
3894                      assert (obj);
3895                      obj->setPreferredWay(way);
3896                    }
3897#endif
3898                  }
3899                  if (nodeStrategy==1||nodeStrategy>3) {
3900                    // depth
3901                    CbcCompareDefault compare;
3902                    compare.setWeight(-3.0);
3903                    babModel->setNodeComparison(compare);
3904                  }
3905                }
3906                if (cppValue>=0) {
3907                  int prepro = useStrategy ? -1 : preProcess;
3908                  // generate code
3909                  FILE * fp = fopen("user_driver.cpp","w");
3910                  if (fp) {
3911                    // generate enough to do BAB
3912                    babModel->generateCpp(fp,1);
3913                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3914                    // Make general so do factorization
3915                    int factor = osiclp->getModelPtr()->factorizationFrequency();
3916                    osiclp->getModelPtr()->setFactorizationFrequency(200);
3917                    osiclp->generateCpp(fp);
3918                    osiclp->getModelPtr()->setFactorizationFrequency(factor);
3919                    //solveOptions.generateCpp(fp);
3920                    fclose(fp);
3921                    // now call generate code
3922                    generateCode(babModel,"user_driver.cpp",cppValue,prepro);
3923                  } else {
3924                    std::cout<<"Unable to open file user_driver.cpp"<<std::endl;
3925                  }
3926                }
3927                if (!babModel->numberStrong())
3928                  babModel->setNumberBeforeTrust(0);
3929                if (useStrategy) {
3930                  CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
3931                  strategy.setupPreProcessing(1);
3932                  babModel->setStrategy(strategy);
3933                }
3934                if (testOsiOptions>=0) {
3935                  printf("Testing OsiObject options %d\n",testOsiOptions);
3936                  if (!numberSOS) {
3937                    babModel->solver()->findIntegersAndSOS(false);
3938#ifdef COIN_HAS_LINK
3939                    // If linked then pass in model
3940                    OsiSolverLink * solver3 = dynamic_cast<OsiSolverLink *> (babModel->solver());
3941                    if (solver3) {
3942                      int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000;
3943                      CglStored stored;
3944                      if (options) {
3945                        printf("nlp options %d\n",options);
3946                        /*
3947                          1 - force mini branch and bound
3948                          2 - set priorities high on continuous
3949                          4 - try adding OA cuts
3950                          8 - try doing quadratic linearization
3951                          16 - try expanding knapsacks
3952                          32 - OA cuts strictly concave
3953                        */
3954                        if ((options&2)) {
3955                          solver3->setBiLinearPriorities(10,tightenFactor > 0.0 ? tightenFactor : 1.0);
3956                        } else if (tightenFactor>0.0) {
3957                          // set grid size for all continuous bi-linear
3958                          solver3->setMeshSizes(tightenFactor);
3959                        }
3960                        if ((options&4)) {
3961                          solver3->setSpecialOptions2(solver3->specialOptions2()|(8+4));
3962                          // say convex
3963                          solver3->sayConvex((options&32)==0);
3964                        }
3965                        if ((options&1)!=0&&slpValue>0)
3966                          solver3->setFixedPriority(slpValue);
3967                        double cutoff=COIN_DBL_MAX;
3968                        if ((options&8))
3969                          cutoff=solver3->linearizedBAB(&stored);
3970                        if (cutoff<babModel->getCutoff())
3971                          babModel->setCutoff(cutoff);
3972                      }
3973                      solver3->setCbcModel(babModel);
3974                      babModel->addCutGenerator(&stored,1,"Stored");
3975                      CglTemporary temp;
3976                      babModel->addCutGenerator(&temp,1,"OnceOnly");
3977                      //choose.setNumberBeforeTrusted(2000);
3978                      //choose.setNumberStrong(20);
3979                    }
3980#endif
3981                  } else {
3982                    // move across
3983                    babModel->deleteObjects(false);
3984                    //babModel->addObjects(babModel->solver()->numberObjects(),babModel->solver()->objects());
3985                  }
3986                  CbcBranchDefaultDecision decision;
3987                  if (babModel->numberStrong()) {
3988                    OsiChooseStrong choose(babModel->solver());
3989                    choose.setNumberBeforeTrusted(babModel->numberBeforeTrust());
3990                    choose.setNumberStrong(babModel->numberStrong());
3991                    choose.setShadowPriceMode(testOsiOptions);
3992                    decision.setChooseMethod(choose);
3993                  } else {
3994                    OsiChooseVariable choose(babModel->solver());
3995                    decision.setChooseMethod(choose);
3996                  }
3997                  babModel->setBranchingMethod(decision);
3998                  if (useCosts&&testOsiOptions>=0) {
3999                    int numberColumns = babModel->getNumCols();
4000                    int * sort = new int[numberColumns];
4001                    double * dsort = new double[numberColumns];
4002                    int * priority = new int [numberColumns];
4003                    const double * objective = babModel->getObjCoefficients();
4004                    const double * lower = babModel->getColLower() ;
4005                    const double * upper = babModel->getColUpper() ;
4006                    const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
4007                    const int * columnLength = matrix->getVectorLengths();
4008                    int iColumn;
4009                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4010                      sort[iColumn]=iColumn;
4011                      if (useCosts==1)
4012                        dsort[iColumn]=-fabs(objective[iColumn]);
4013                      else if (useCosts==2)
4014                        dsort[iColumn]=iColumn;
4015                      else if (useCosts==3)
4016                        dsort[iColumn]=upper[iColumn]-lower[iColumn];
4017                      else if (useCosts==4)
4018                        dsort[iColumn]=-(upper[iColumn]-lower[iColumn]);
4019                      else if (useCosts==5)
4020                        dsort[iColumn]=-columnLength[iColumn];
4021                    }
4022                    CoinSort_2(dsort,dsort+numberColumns,sort);
4023                    int level=0;
4024                    double last = -1.0e100;
4025                    for (int i=0;i<numberColumns;i++) {
4026                      int iPut=sort[i];
4027                      if (dsort[i]!=last) {
4028                        level++;
4029                        last=dsort[i];
4030                      }
4031                      priority[iPut]=level;
4032                    }
4033                    OsiObject ** objects = babModel->objects();
4034                    int numberObjects = babModel->numberObjects();
4035                    for (int iObj = 0;iObj<numberObjects;iObj++) {
4036                      OsiObject * obj = objects[iObj] ;
4037                      int iColumn = obj->columnNumber();
4038                      if (iColumn>=0)
4039                        obj->setPriority(priority[iColumn]);
4040                    }
4041                    delete [] priority;
4042                    delete [] sort;
4043                    delete [] dsort;
4044                  }
4045                }
4046                checkSOS(babModel, babModel->solver());
4047                if (doSprint>0) {
4048                  // Sprint for primal solves
4049                  ClpSolve::SolveType method = ClpSolve::usePrimalorSprint;
4050                  ClpSolve::PresolveType presolveType = ClpSolve::presolveOff;
4051                  int numberPasses = 5;
4052                  int options[] = {0,3,0,0,0,0};
4053                  int extraInfo[] = {-1,20,-1,-1,-1,-1};
4054                  extraInfo[1]=doSprint;
4055                  int independentOptions[] = {0,0,3};
4056                  ClpSolve clpSolve(method,presolveType,numberPasses,
4057                                    options,extraInfo,independentOptions);
4058                  // say use in OsiClp
4059                  clpSolve.setSpecialOption(6,1);
4060                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4061                  osiclp->setSolveOptions(clpSolve);
4062                  osiclp->setHintParam(OsiDoDualInResolve,false);
4063                  // switch off row copy
4064                  osiclp->getModelPtr()->setSpecialOptions(osiclp->getModelPtr()->specialOptions()|256);
4065                  osiclp->getModelPtr()->setInfeasibilityCost(1.0e11);
4066                }
4067#ifdef COIN_HAS_LINK
4068                if (storedAmpl.sizeRowCuts()) {
4069                  if (preProcess) {
4070                    const int * originalColumns = process.originalColumns();
4071                    int numberColumns = babModel->getNumCols();
4072                    int * newColumn = new int[numberOriginalColumns];
4073                    int i;
4074                    for (i=0;i<numberOriginalColumns;i++) 
4075                      newColumn[i]=-1;
4076                    for (i=0;i<numberColumns;i++) {
4077                      int iColumn = originalColumns[i];
4078                      newColumn[iColumn]=i;
4079                    }
4080                    int * buildColumn = new int[numberColumns];
4081                    // Build up valid cuts
4082                    int nBad=0;
4083                    int nCuts = storedAmpl.sizeRowCuts();
4084                    CglStored newCuts;
4085                    for (i=0;i<nCuts;i++) {
4086                      const OsiRowCut * cut = storedAmpl.rowCutPointer(i);
4087                      double lb = cut->lb();
4088                      double ub = cut->ub();
4089                      int n=cut->row().getNumElements();
4090                      const int * column = cut->row().getIndices();
4091                      const double * element = cut->row().getElements();
4092                      bool bad=false;
4093                      for (int i=0;i<n;i++) {
4094                        int iColumn = column[i];
4095                        iColumn = newColumn[iColumn];
4096                        if (iColumn>=0) {
4097                          buildColumn[i]=iColumn;
4098                        } else {
4099                          bad=true;
4100                          break;
4101                        }
4102                      }
4103                      if (!bad) {
4104                        newCuts.addCut(lb,ub,n,buildColumn,element);
4105                      } else {
4106                        nBad++;
4107                      }
4108                    }
4109                    storedAmpl=newCuts;
4110                    if (nBad)
4111                      printf("%d cuts dropped\n",nBad);
4112                    delete [] newColumn;
4113                    delete [] buildColumn;
4114                  }
4115                  if (storedAmpl.sizeRowCuts()) {
4116                    //babModel->addCutGenerator(&storedAmpl,1,"AmplStored");
4117                    int numberRowCuts = storedAmpl.sizeRowCuts();
4118                    for (int i=0;i<numberRowCuts;i++) {
4119                      const OsiRowCut * rowCutPointer = storedAmpl.rowCutPointer(i);
4120                      babModel->makeGlobalCut(rowCutPointer);
4121                    }
4122                  }
4123                }
4124#endif
4125#ifdef CLP_MALLOC_STATISTICS
4126                malloc_stats();
4127                malloc_stats2();
4128#endif
4129                if (outputFormat==5) {
4130                  osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4131                  lpSolver = osiclp->getModelPtr();
4132                  lpSolver->setPersistenceFlag(1);
4133                }
4134#ifdef COIN_HAS_ASL
4135                // add in lotsizing
4136                if (usingAmpl&&info.special) {
4137                  int numberColumns = babModel->getNumCols();
4138                  int i;
4139                  int n=0;
4140                  if (preProcess) {
4141                    const int * originalColumns = process.originalColumns();
4142                    for (i=0;i<numberColumns;i++) {
4143                      int iColumn = originalColumns[i];
4144                      assert (iColumn>=i);
4145                      int iType = info.special[iColumn];
4146                      if (iType) {
4147                        assert (iType==1);
4148                        n++;
4149                      }
4150                      info.special[i]=iType;
4151                    }
4152                  }
4153                  if (n) {
4154                    int numberIntegers=0;
4155                    int numberOldObjects=0;
4156                    OsiObject ** oldObjects=NULL;
4157                    const double * lower = babModel->solver()->getColLower();
4158                    const double * upper = babModel->solver()->getColUpper();
4159                    if (testOsiOptions<0) {
4160                      // *************************************************************
4161                      // CbcObjects
4162                      numberIntegers = babModel->numberIntegers();
4163                      /* model may not have created objects
4164                         If none then create
4165                      */
4166                      if (!numberIntegers||!babModel->numberObjects()) {
4167                        int type = (pseudoUp) ? 1 : 0;
4168                        babModel->findIntegers(true,type);
4169                        numberIntegers = babModel->numberIntegers();
4170                      }
4171                      oldObjects = babModel->objects();
4172                      numberOldObjects = babModel->numberObjects();
4173                    } else {
4174                      numberIntegers = testOsiSolver->getNumIntegers();
4175                      if (!numberIntegers||!testOsiSolver->numberObjects()) {
4176                        /* model may not have created objects
4177                           If none then create
4178                        */
4179                        testOsiSolver->findIntegers(false);
4180                        numberIntegers = testOsiSolver->getNumIntegers();
4181                      }
4182                      oldObjects = testOsiSolver->objects();
4183                      numberOldObjects = testOsiSolver->numberObjects();
4184                    }
4185                    OsiObject ** objects = new OsiObject * [n];
4186                    n=0;
4187                    // set new objects to have one lower priority
4188                    double ranges[] = {-COIN_DBL_MAX,-1.0,1.0,COIN_DBL_MAX};
4189                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
4190                      int iColumn = oldObjects[iObj]->columnNumber();
4191                      if (iColumn>=0&&info.special[iColumn]) {
4192                        if (lower[iColumn]<=-1.0&&upper[iColumn]>=0.0) {
4193                          ranges[0]=lower[iColumn];
4194                          ranges[3]=upper[iColumn];
4195                          int priority = oldObjects[iObj]->priority();
4196                          if (testOsiOptions<0) {
4197                            objects[n] = new CbcLotsize(babModel,iColumn,2,ranges,true);
4198                          } else {
4199                            objects[n] = new OsiLotsize(testOsiSolver,iColumn,2,ranges,true);
4200                          }
4201                          objects[n++]->setPriority (priority-1);
4202                        }
4203                      }
4204                    }
4205                    if (testOsiOptions<0) {
4206                      babModel->addObjects(n,objects);
4207                    } else {
4208                      testOsiSolver->addObjects(n,objects);
4209                    }
4210                    for (i=0;i<n;i++)
4211                      delete objects[i];
4212                    delete [] objects;
4213                  }
4214                }
4215#endif
4216                babModel->branchAndBound(statistics);
4217#ifdef CLP_MALLOC_STATISTICS
4218                malloc_stats();
4219                malloc_stats2();
4220#endif
4221                checkSOS(babModel, babModel->solver());
4222              } else if (type==MIPLIB) {
4223                CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
4224                // Set up pre-processing
4225                int translate2[]={9999,1,1,3,2,4,5};
4226                if (preProcess)
4227                  strategy.setupPreProcessing(translate2[preProcess]);
4228                babModel->setStrategy(strategy);
4229                if (outputFormat==5) {
4230                  osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4231                  lpSolver = osiclp->getModelPtr();
4232                  lpSolver->setPersistenceFlag(1);
4233                }
4234                if (testOsiOptions>=0) {
4235                  printf("Testing OsiObject options %d\n",testOsiOptions);
4236                  CbcBranchDefaultDecision decision;
4237                  OsiChooseStrong choose(babModel->solver());
4238                  choose.setNumberBeforeTrusted(babModel->numberBeforeTrust());
4239                  choose.setNumberStrong(babModel->numberStrong());
4240                  choose.setShadowPriceMode(testOsiOptions);
4241                  //babModel->deleteObjects(false);
4242                  decision.setChooseMethod(choose);
4243                  babModel->setBranchingMethod(decision);
4244                }
4245                *miplibSolver = babModel;
4246                return 777;
4247              } else {
4248                strengthenedModel = babModel->strengthenedModel();
4249              }
4250              currentBranchModel = NULL;
4251              osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4252              if (debugFile=="createAfterPre"&&babModel->bestSolution()) {
4253                lpSolver = osiclp->getModelPtr();
4254                //move best solution (should be there -- but ..)
4255                int n = lpSolver->getNumCols();
4256                memcpy(lpSolver->primalColumnSolution(),babModel->bestSolution(),n*sizeof(double));
4257                saveSolution(osiclp->getModelPtr(),"debug.file");
4258              }
4259              if (!noPrinting) {
4260                // Print more statistics
4261                std::cout<<"Cuts at root node changed objective from "<<babModel->getContinuousObjective()
4262                         <<" to "<<babModel->rootObjectiveAfterCuts()<<std::endl;
4263               
4264                numberGenerators = babModel->numberCutGenerators();
4265                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
4266                  CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
4267                  std::cout<<generator->cutGeneratorName()<<" was tried "
4268                           <<generator->numberTimesEntered()<<" times and created "
4269                           <<generator->numberCutsInTotal()<<" cuts of which "
4270                           <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
4271                  if (generator->timing())
4272                    std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
4273                  else
4274                    std::cout<<std::endl;
4275                }
4276              }
4277              time2 = CoinCpuTime();
4278              totalTime += time2-time1;
4279              // For best solution
4280              double * bestSolution = NULL;
4281              if (babModel->getMinimizationObjValue()<1.0e50&&type==BAB) {
4282                // post process
4283                if (preProcess) {
4284                  int n = saveSolver->getNumCols();
4285                  bestSolution = new double [n];
4286                  process.postProcess(*babModel->solver());
4287                  // Solution now back in saveSolver
4288                  babModel->assignSolver(saveSolver);
4289                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
4290                } else {
4291                  int n = babModel->solver()->getNumCols();
4292                  bestSolution = new double [n];
4293                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
4294                }
4295                // and put back in very original solver
4296                {
4297                  ClpSimplex * original = originalSolver.getModelPtr();
4298                  double * lower = original->columnLower();
4299                  double * upper = original->columnUpper();
4300                  double * solution = original->primalColumnSolution();
4301                  int n = original->numberColumns();
4302                  assert (!n||n==babModel->solver()->getNumCols());
4303                  for (int i=0;i<n;i++) {
4304                    solution[i]=bestSolution[i];
4305                    if (solver1.isInteger(i)) {
4306                      lower[i]=solution[i];
4307                      upper[i]=solution[i];
4308                    }
4309                  }
4310                }
4311                checkSOS(babModel, babModel->solver());
4312              }
4313              if (type==STRENGTHEN&&strengthenedModel)
4314                clpSolver = dynamic_cast< OsiClpSolverInterface*> (strengthenedModel);
4315#ifdef COIN_HAS_ASL
4316              else if (usingAmpl) 
4317                clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4318#endif
4319              lpSolver = clpSolver->getModelPtr();
4320              if (numberChanged) {
4321                for (int i=0;i<numberChanged;i++) {
4322                  int iColumn=changed[i];
4323                  clpSolver->setContinuous(iColumn);
4324                }
4325                delete [] changed;
4326              }
4327              if (type==BAB) {
4328                //move best solution (should be there -- but ..)
4329                int n = lpSolver->getNumCols();
4330                if (bestSolution) {
4331                  memcpy(lpSolver->primalColumnSolution(),bestSolution,n*sizeof(double));
4332                  // now see what that does to row solution
4333                  int numberRows=lpSolver->numberRows();
4334                  double * rowSolution = lpSolver->primalRowSolution();
4335                  memset (rowSolution,0,numberRows*sizeof(double));
4336                  lpSolver->clpMatrix()->times(1.0,bestSolution,rowSolution);
4337                }
4338                if (debugFile=="create"&&bestSolution) {
4339                  saveSolution(lpSolver,"debug.file");
4340                }
4341                delete [] bestSolution;
4342                std::string statusName[]={"Finished","Stopped on ","Difficulties",
4343                                          "","","User ctrl-c"};
4344                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
4345                int iStat = babModel->status();
4346                int iStat2 = babModel->secondaryStatus();
4347                if (!noPrinting)
4348                  std::cout<<"Result - "<<statusName[iStat]<<minor[iStat2]
4349                           <<" objective "<<babModel->getObjValue()<<
4350                    " after "<<babModel->getNodeCount()<<" nodes and "
4351                           <<babModel->getIterationCount()<<
4352                    " iterations - took "<<time2-time1<<" seconds"<<std::endl;
4353#ifdef COIN_HAS_ASL
4354                if (usingAmpl) {
4355                  clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4356                  lpSolver = clpSolver->getModelPtr();
4357                  double value = babModel->getObjValue()*lpSolver->getObjSense();
4358                  char buf[300];
4359                  int pos=0;
4360                  if (iStat==0) {
4361                    if (babModel->getObjValue()<1.0e40) {
4362                      pos += sprintf(buf+pos,"optimal," );
4363                    } else {
4364                      // infeasible
4365                      iStat=1;
4366                      pos += sprintf(buf+pos,"infeasible,");
4367                    }
4368                  } else if (iStat==1) {
4369                    if (iStat2!=6)
4370                      iStat=3;
4371                    else
4372                      iStat=4;
4373                    pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
4374                  } else if (iStat==2) {
4375                    iStat = 7;
4376                    pos += sprintf(buf+pos,"stopped on difficulties,");
4377                  } else if (iStat==5) {
4378                    iStat = 3;
4379                    pos += sprintf(buf+pos,"stopped on ctrl-c,");
4380                  } else {
4381                    pos += sprintf(buf+pos,"status unknown,");
4382                    iStat=6;
4383                  }
4384                  info.problemStatus=iStat;
4385                  info.objValue = value;
4386                  if (babModel->getObjValue()<1.0e40) {
4387                    int precision = ampl_obj_prec();
4388                    if (precision>0)
4389                      pos += sprintf(buf+pos," objective %.*g",precision,
4390                                     value);
4391                    else
4392                      pos += sprintf(buf+pos," objective %g",value);
4393                  }
4394                  sprintf(buf+pos,"\n%d nodes, %d iterations, %g seconds",
4395                          babModel->getNodeCount(),
4396                          babModel->getIterationCount(),
4397                          totalTime);
4398                  if (bestSolution) {
4399                    free(info.primalSolution);
4400                    if (!numberKnapsack) {
4401                      info.primalSolution = (double *) malloc(n*sizeof(double));
4402                      CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
4403                      int numberRows = lpSolver->numberRows();
4404                      free(info.dualSolution);
4405                      info.dualSolution = (double *) malloc(numberRows*sizeof(double));
4406                      CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
4407                    } else {
4408                      // expanded knapsack
4409                      info.dualSolution=NULL;
4410                      int numberColumns = saveCoinModel.numberColumns();
4411                      info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
4412                      // Fills in original solution (coinModel length)
4413                      afterKnapsack(saveCoinModel,  whichColumn,  knapsackStart, 
4414                                    knapsackRow,  numberKnapsack,
4415                                    lpSolver->primalColumnSolution(), info.primalSolution,1);
4416                    }
4417                  } else {
4418                    info.primalSolution=NULL;
4419                    info.dualSolution=NULL;
4420                  }
4421                  // put buffer into info
4422                  strcpy(info.buffer,buf);
4423                }
4424#endif
4425              } else {
4426                std::cout<<"Model strengthened - now has "<<clpSolver->getNumRows()
4427                         <<" rows"<<std::endl;
4428              }
4429              time1 = time2;
4430              delete babModel;
4431              babModel=NULL;
4432            } else {
4433              std::cout << "** Current model not valid" << std::endl ; 
4434            }
4435            break ;
4436          case IMPORT:
4437            {
4438#ifdef COIN_HAS_ASL
4439              if (!usingAmpl) {
4440#endif
4441                free(priorities);
4442                priorities=NULL;
4443                free(branchDirection);
4444                branchDirection=NULL;
4445                free(pseudoDown);
4446                pseudoDown=NULL;
4447                free(pseudoUp);
4448                pseudoUp=NULL;
4449                free(solutionIn);
4450                solutionIn=NULL;
4451                free(prioritiesIn);
4452                prioritiesIn=NULL;
4453                free(sosStart);
4454                sosStart=NULL;
4455                free(sosIndices);
4456                sosIndices=NULL;
4457                free(sosType);
4458                sosType=NULL;
4459                free(sosReference);
4460                sosReference=NULL;
4461                free(cut);
4462                cut=NULL;
4463                free(sosPriority);
4464                sosPriority=NULL;
4465#ifdef COIN_HAS_ASL
4466              }
4467#endif               
4468              delete babModel;
4469              babModel=NULL;
4470              // get next field
4471              field = CoinReadGetString(argc,argv);
4472              if (field=="$") {
4473                field = parameters[iParam].stringValue();
4474              } else if (field=="EOL") {
4475                parameters[iParam].printString();
4476                break;
4477              } else {
4478                parameters[iParam].setStringValue(field);
4479              }
4480              std::string fileName;
4481              bool canOpen=false;
4482              // See if gmpl file
4483              int gmpl=0;
4484              std::string gmplData;
4485              if (field=="-") {
4486                // stdin
4487                canOpen=true;
4488                fileName = "-";
4489              } else {
4490                // See if .lp
4491                {
4492                  const char * c_name = field.c_str();
4493                  int length = strlen(c_name);
4494                  if (length>3&&!strncmp(c_name+length-3,".lp",3))
4495                    gmpl=-1; // .lp
4496                }
4497                bool absolutePath;
4498                if (dirsep=='/') {
4499                  // non Windows (or cygwin)
4500                  absolutePath=(field[0]=='/');
4501                } else {
4502                  //Windows (non cycgwin)
4503                  absolutePath=(field[0]=='\\');
4504                  // but allow for :
4505                  if (strchr(field.c_str(),':'))
4506                    absolutePath=true;
4507                }
4508                if (absolutePath) {
4509                  fileName = field;
4510                } else if (field[0]=='~') {
4511                  char * environVar = getenv("HOME");
4512                  if (environVar) {
4513                    std::string home(environVar);
4514                    field=field.erase(0,1);
4515                    fileName = home+field;
4516                  } else {
4517                    fileName=field;
4518                  }
4519                } else {
4520                  fileName = directory+field;
4521                  // See if gmpl (model & data) - or even lp file
4522                  int length = field.size();
4523                  int percent = field.find('%');
4524                  if (percent<length&&percent>0) {
4525                    gmpl=1;
4526                    fileName = directory+field.substr(0,percent);
4527                    gmplData = directory+field.substr(percent+1);
4528                    if (percent<length-1)
4529                      gmpl=2; // two files
4530                    printf("GMPL model file %s and data file %s\n",
4531                           fileName.c_str(),gmplData.c_str());
4532                  }
4533                }
4534                std::string name=fileName;
4535                if (fileCoinReadable(name)) {
4536                  // can open - lets go for it
4537                  canOpen=true;
4538                  if (gmpl==2) {
4539                    FILE *fp;
4540                    fp=fopen(gmplData.c_str(),"r");
4541                    if (fp) {
4542                      fclose(fp);
4543                    } else {
4544                      canOpen=false;
4545                      std::cout<<"Unable to open file "<<gmplData<<std::endl;
4546                    }
4547                  }
4548                } else {
4549                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4550                }
4551              }
4552              if (canOpen) {
4553                int status;
4554                ClpSimplex * lpSolver = clpSolver->getModelPtr();
4555                if (!gmpl)
4556                  status =clpSolver->readMps(fileName.c_str(),
4557                                                 keepImportNames!=0,
4558                                                 allowImportErrors!=0);
4559                else if (gmpl>0)
4560                  status= lpSolver->readGMPL(fileName.c_str(),
4561                                                  (gmpl==2) ? gmplData.c_str() : NULL,
4562                                                  keepImportNames!=0);
4563                else
4564                  status= lpSolver->readLp(fileName.c_str(),1.0e-12);
4565                if (!status||(status>0&&allowImportErrors)) {
4566                  goodModel=true;
4567                  // sets to all slack (not necessary?)
4568                  lpSolver->createStatus();
4569                  // make sure integer
4570                  int numberColumns = lpSolver->numberColumns();
4571                  for (int i=0;i<numberColumns;i++) {
4572                    if (lpSolver->isInteger(i))
4573                      clpSolver->setInteger(i);
4574                  }
4575                  time2 = CoinCpuTime();
4576                  totalTime += time2-time1;
4577                  time1=time2;
4578                  // Go to canned file if just input file
4579                  if (CbcOrClpRead_mode==2&&argc==2) {
4580                    // only if ends .mps
4581                    char * find = (char *)strstr(fileName.c_str(),".mps");
4582                    if (find&&find[4]=='\0') {
4583                      find[1]='p'; find[2]='a';find[3]='r';
4584                      FILE *fp=fopen(fileName.c_str(),"r");
4585                      if (fp) {
4586                        CbcOrClpReadCommand=fp; // Read from that file
4587                        CbcOrClpRead_mode=-1;
4588                      }
4589                    }
4590                  }
4591                } else {
4592                  // errors
4593                  std::cout<<"There were "<<status<<
4594                    " errors on input"<<std::endl;
4595                }
4596              }
4597            }
4598            break;
4599          case MODELIN:
4600#ifdef COIN_HAS_LINK
4601            {
4602              // get next field
4603              field = CoinReadGetString(argc,argv);
4604              if (field=="$") {
4605                field = parameters[iParam].stringValue();
4606              } else if (field=="EOL") {
4607                parameters[iParam].printString();
4608                break;
4609              } else {
4610                parameters[iParam].setStringValue(field);
4611              }
4612              std::string fileName;
4613              bool canOpen=false;
4614              if (field=="-") {
4615                // stdin
4616                canOpen=true;
4617                fileName = "-";
4618              } else {
4619                bool absolutePath;
4620                if (dirsep=='/') {
4621                  // non Windows (or cygwin)
4622                  absolutePath=(field[0]=='/');
4623                } else {
4624                  //Windows (non cycgwin)
4625                  absolutePath=(field[0]=='\\');
4626                  // but allow for :
4627                  if (strchr(field.c_str(),':'))
4628                    absolutePath=true;
4629                }
4630                if (absolutePath) {
4631                  fileName = field;
4632                } else if (field[0]=='~') {
4633                  char * environVar = getenv("HOME");
4634                  if (environVar) {
4635                    std::string home(environVar);
4636                    field=field.erase(0,1);
4637                    fileName = home+field;
4638                  } else {
4639                    fileName=field;
4640                  }
4641                } else {
4642                  fileName = directory+field;
4643                }
4644                FILE *fp=fopen(fileName.c_str(),"r");
4645                if (fp) {
4646                  // can open - lets go for it
4647                  fclose(fp);
4648                  canOpen=true;
4649                } else {
4650                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4651                }
4652              }
4653              if (canOpen) {
4654                CoinModel coinModel(fileName.c_str(),2);
4655                // load from coin model
4656                OsiSolverLink solver1;
4657                OsiSolverInterface * solver2 = solver1.clone();
4658                model.assignSolver(solver2,true);
4659                OsiSolverLink * si =
4660                  dynamic_cast<OsiSolverLink *>(model.solver()) ;
4661                assert (si != NULL);
4662                si->setDefaultMeshSize(0.001);
4663                // need some relative granularity
4664                si->setDefaultBound(100.0);
4665                si->setDefaultMeshSize(0.01);
4666                si->setDefaultBound(100.0);
4667                si->setIntegerPriority(1000);
4668                si->setBiLinearPriority(10000);
4669                CoinModel * model2 = (CoinModel *) &coinModel;
4670                si->load(*model2);
4671                // redo
4672                solver = model.solver();
4673                clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
4674                lpSolver = clpSolver->getModelPtr();
4675                clpSolver->messageHandler()->setLogLevel(0) ;
4676                testOsiParameters=0;
4677                complicatedInteger=2;
4678              }
4679            }
4680#endif
4681            break;
4682          case EXPORT:
4683            if (goodModel) {
4684              // get next field
4685              field = CoinReadGetString(argc,argv);
4686              if (field=="$") {
4687                field = parameters[iParam].stringValue();
4688              } else if (field=="EOL") {
4689                parameters[iParam].printString();
4690                break;
4691              } else {
4692                parameters[iParam].setStringValue(field);
4693              }
4694              std::string fileName;
4695              bool canOpen=false;
4696              if (field[0]=='/'||field[0]=='\\') {
4697                fileName = field;
4698              } else if (field[0]=='~') {
4699                char * environVar = getenv("HOME");
4700                if (environVar) {
4701                  std::string home(environVar);
4702                  field=field.erase(0,1);
4703                  fileName = home+field;
4704                } else {
4705                  fileName=field;
4706                }
4707              } else {
4708                fileName = directory+field;
4709              }
4710              FILE *fp=fopen(fileName.c_str(),"w");
4711              if (fp) {
4712                // can open - lets go for it
4713                fclose(fp);
4714                canOpen=true;
4715              } else {
4716                std::cout<<"Unable to open file "<<fileName<<std::endl;
4717              }
4718              if (canOpen) {
4719                // If presolve on then save presolved
4720                bool deleteModel2=false;
4721                ClpSimplex * model2 = lpSolver;
4722                if (dualize) {
4723                  model2 = ((ClpSimplexOther *) model2)->dualOfModel();
4724                  printf("Dual of model has %d rows and %d columns\n",
4725                         model2->numberRows(),model2->numberColumns());
4726                  model2->setOptimizationDirection(1.0);
4727                }
4728#ifdef COIN_HAS_ASL
4729                if (info.numberSos&&doSOS&&usingAmpl) {
4730                  // SOS
4731                  numberSOS = info.numberSos;
4732                  sosStart = info.sosStart;
4733                  sosIndices = info.sosIndices;
4734                  sosReference = info.sosReference;
4735                  preSolve=false;
4736                  clpSolver->setSOSData(numberSOS,info.sosType,sosStart,sosIndices,sosReference);
4737                }
4738#endif
4739                if (preSolve) {
4740                  ClpPresolve pinfo;
4741                  int presolveOptions2 = presolveOptions&~0x40000000;
4742                  if ((presolveOptions2&0xffff)!=0)
4743                    pinfo.setPresolveActions(presolveOptions2);
4744                  if ((printOptions&1)!=0)
4745                    pinfo.statistics();
4746                  double presolveTolerance = 
4747                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
4748                  model2 = 
4749                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
4750                                         true,preSolve);
4751                  if (model2) {
4752                    printf("Saving presolved model on %s\n",
4753                           fileName.c_str());
4754                    deleteModel2=true;
4755                  } else {
4756                    printf("Presolved model looks infeasible - saving original on %s\n",
4757                           fileName.c_str());
4758                    deleteModel2=false;
4759                    model2 = lpSolver;
4760
4761                  }
4762                  model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
4763                  if (deleteModel2)
4764                    delete model2;
4765                } else {
4766                  printf("Saving model on %s\n",
4767                           fileName.c_str());
4768                  if (numberSOS) {
4769                    // Convert names
4770                    int iRow;
4771                    int numberRows=model2->numberRows();
4772                    int iColumn;
4773                    int numberColumns=model2->numberColumns();
4774                   
4775                    char ** rowNames = NULL;
4776                    char ** columnNames = NULL;
4777                    if (model2->lengthNames()) {
4778                      rowNames = new char * [numberRows];
4779                      for (iRow=0;iRow<numberRows;iRow++) {
4780                        rowNames[iRow] = 
4781                          strdup(model2->rowName(iRow).c_str());
4782                      }
4783                     
4784                      columnNames = new char * [numberColumns];
4785                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
4786                        columnNames[iColumn] = 
4787                          strdup(model2->columnName(iColumn).c_str());
4788                      }
4789                    }
4790                    clpSolver->writeMpsNative(fileName.c_str(),(const char **) rowNames,(const char **) columnNames,
4791                                              (outputFormat-1)/2,1+((outputFormat-1)&1));
4792                    if (rowNames) {
4793                      for (iRow=0;iRow<numberRows;iRow++) {
4794                        free(rowNames[iRow]);
4795                      }
4796                      delete [] rowNames;
4797                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
4798                        free(columnNames[iColumn]);
4799                      }
4800                      delete [] columnNames;
4801                    }
4802                  } else {
4803                    model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
4804                  }
4805                }
4806                time2 = CoinCpuTime();
4807                totalTime += time2-time1;
4808                time1=time2;
4809              }
4810            } else {
4811              std::cout<<"** Current model not valid"<<std::endl;
4812            }
4813            break;
4814          case BASISIN:
4815            if (goodModel) {
4816              // get next field
4817              field = CoinReadGetString(argc,argv);
4818              if (field=="$") {
4819                field = parameters[iParam].stringValue();
4820              } else if (field=="EOL") {
4821                parameters[iParam].printString();
4822                break;
4823              } else {
4824                parameters[iParam].setStringValue(field);
4825              }
4826              std::string fileName;
4827              bool canOpen=false;
4828              if (field=="-") {
4829                // stdin
4830                canOpen=true;
4831                fileName = "-";
4832              } else {
4833                if (field[0]=='/'||field[0]=='\\') {
4834                  fileName = field;
4835                } else if (field[0]=='~') {
4836                  char * environVar = getenv("HOME");
4837                  if (environVar) {
4838                    std::string home(environVar);
4839                    field=field.erase(0,1);
4840                    fileName = home+field;
4841                  } else {
4842                    fileName=field;
4843                  }
4844                } else {
4845                  fileName = directory+field;
4846                }
4847                FILE *fp=fopen(fileName.c_str(),"r");
4848                if (fp) {
4849                  // can open - lets go for it
4850                  fclose(fp);
4851                  canOpen=true;
4852                } else {
4853                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4854                }
4855              }
4856              if (canOpen) {
4857                int values = lpSolver->readBasis(fileName.c_str());
4858                if (values==0)
4859                  basisHasValues=-1;
4860                else
4861                  basisHasValues=1;
4862              }
4863            } else {
4864              std::cout<<"** Current model not valid"<<std::endl;
4865            }
4866            break;
4867          case PRIORITYIN:
4868            if (goodModel) {
4869              // get next field
4870              field = CoinReadGetString(argc,argv);
4871              if (field=="$") {
4872                field = parameters[iParam].stringValue();
4873              } else if (field=="EOL") {
4874                parameters[iParam].printString();
4875                break;
4876              } else {
4877                parameters[iParam].setStringValue(field);
4878              }
4879              std::string fileName;
4880              if (field[0]=='/'||field[0]=='\\') {
4881                fileName = field;
4882              } else if (field[0]=='~') {
4883                char * environVar = getenv("HOME");
4884                if (environVar) {
4885                  std::string home(environVar);
4886                  field=field.erase(0,1);
4887                  fileName = home+field;
4888                } else {
4889                  fileName=field;
4890                }
4891              } else {
4892                fileName = directory+field;
4893              }
4894              FILE *fp=fopen(fileName.c_str(),"r");
4895              if (fp) {
4896                // can open - lets go for it
4897                std::string headings[]={"name","number","direction","priority","up","down",
4898                                        "solution","priin"};
4899                int got[]={-1,-1,-1,-1,-1,-1,-1,-1};
4900                int order[8];
4901                assert(sizeof(got)==sizeof(order));
4902                int nAcross=0;
4903                char line[1000];
4904                int numberColumns = lpSolver->numberColumns();
4905                if (!fgets(line,1000,fp)) {
4906                  std::cout<<"Odd file "<<fileName<<std::endl;
4907                } else {
4908                  char * pos = line;
4909                  char * put = line;
4910                  while (*pos>=' '&&*pos!='\n') {
4911                    if (*pos!=' '&&*pos!='\t') {
4912                      *put=tolower(*pos);
4913                      put++;
4914                    }
4915                    pos++;
4916                  }
4917                  *put='\0';
4918                  pos=line;
4919                  int i;
4920                  bool good=true;
4921                  while (pos) {
4922                    char * comma = strchr(pos,',');
4923                    if (comma)
4924                      *comma='\0';
4925                    for (i=0;i<(int) (sizeof(got)/sizeof(int));i++) {
4926                      if (headings[i]==pos) {
4927                        if (got[i]<0) {
4928                          order[nAcross]=i;
4929                          got[i]=nAcross++;
4930                        } else {
4931                          // duplicate
4932                          good=false;
4933                        }
4934                        break;
4935                      }
4936                    }
4937                    if (i==(int) (sizeof(got)/sizeof(int)))
4938                      good=false;
4939                    if (comma) {
4940                      *comma=',';
4941                      pos=comma+1;
4942                    } else {
4943                      break;
4944                    }
4945                  }
4946                  if (got[0]<0&&got[1]<0)
4947                    good=false;
4948                  if (got[0]>=0&&got[1]>=0)
4949                    good=false;
4950                  if (got[0]>=0&&!lpSolver->lengthNames())
4951                    good=false;
4952                  if (good) {
4953                    char ** columnNames = new char * [numberColumns];
4954                    pseudoDown= (double *) malloc(numberColumns*sizeof(double));
4955                    pseudoUp = (double *) malloc(numberColumns*sizeof(double));
4956                    branchDirection = (int *) malloc(numberColumns*sizeof(int));
4957                    priorities= (int *) malloc(numberColumns*sizeof(int));
4958                    free(solutionIn);
4959                    solutionIn=NULL;
4960                    free(prioritiesIn);
4961                    prioritiesIn=NULL;
4962                    int iColumn;
4963                    if (got[6]>=0) {
4964                      solutionIn = (double *) malloc(numberColumns*sizeof(double));
4965                      CoinZeroN(solutionIn,numberColumns);
4966                    }
4967                    if (got[7]>=0) {
4968                      prioritiesIn = (int *) malloc(numberColumns*sizeof(int));
4969                      for (iColumn=0;iColumn<numberColumns;iColumn++) 
4970                        prioritiesIn[iColumn]=10000;
4971                    }
4972                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4973                      columnNames[iColumn] = 
4974                        strdup(lpSolver->columnName(iColumn).c_str());
4975                      pseudoDown[iColumn]=0.0;
4976                      pseudoUp[iColumn]=0.0;
4977                      branchDirection[iColumn]=0;
4978                      priorities[iColumn]=0;
4979                    }
4980                    int nBadPseudo=0;
4981                    int nBadDir=0;
4982                    int nBadPri=0;
4983                    int nBadName=0;
4984                    int nBadLine=0;
4985                    int nLine=0;
4986                    while (fgets(line,1000,fp)) {
4987                      nLine++;
4988                      iColumn = -1;
4989                      double up =0.0;
4990                      double down=0.0;
4991                      int pri=0;
4992                      int dir=0;
4993                      double solValue=COIN_DBL_MAX;
4994                      int priValue=1000000;
4995                      char * pos = line;
4996                      char * put = line;
4997                      while (*pos>=' '&&*pos!='\n') {
4998                        if (*pos!=' '&&*pos!='\t') {
4999                          *put=tolower(*pos);
5000                          put++;
5001                        }
5002                        pos++;
5003                      }
5004                      *put='\0';
5005                      pos=line;
5006                      for (int i=0;i<nAcross;i++) {
5007