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

Last change on this file since 601 was 601, checked in by forrest, 13 years ago

cuts?

File size: 229.4 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        if (info.cut) {
1409          printf("warning - cuts with LOS\n");
1410          int numberRows = info.numberRows;
1411          int * whichRow = new int [numberRows];
1412          // Row copy
1413          const CoinPackedMatrix * matrixByRow = solver->getMatrixByRow();
1414          const double * elementByRow = matrixByRow->getElements();
1415          const int * column = matrixByRow->getIndices();
1416          const CoinBigIndex * rowStart = matrixByRow->getVectorStarts();
1417          const int * rowLength = matrixByRow->getVectorLengths();
1418         
1419          const double * rowLower = solver->getRowLower();
1420          const double * rowUpper = solver->getRowUpper();
1421          int nDelete=0;
1422          for (int iRow=0;iRow<numberRows;iRow++) {
1423            if (info.cut[iRow]) {
1424              whichRow[nDelete++]=iRow;
1425              int start = rowStart[iRow];
1426              storedAmpl.addCut(rowLower[iRow],rowUpper[iRow],
1427                                rowLength[iRow],column+start,elementByRow+start);
1428            }
1429          }
1430          solver->deleteRows(nDelete,whichRow);
1431          delete [] whichRow;
1432        }
1433      }
1434#endif
1435      // If we had a solution use it
1436      if (info.primalSolution) {
1437        solver->setColSolution(info.primalSolution);
1438      }
1439      // status
1440      if (info.rowStatus) {
1441        unsigned char * statusArray = lpSolver->statusArray();
1442        int i;
1443        for (i=0;i<info.numberColumns;i++)
1444          statusArray[i]=(char)info.columnStatus[i];
1445        statusArray+=info.numberColumns;
1446        for (i=0;i<info.numberRows;i++)
1447          statusArray[i]=(char)info.rowStatus[i];
1448        CoinWarmStartBasis * basis = lpSolver->getBasis();
1449        solver->setWarmStart(basis);
1450        delete basis;
1451      }
1452      freeArrays1(&info);
1453      // modify objective if necessary
1454      solver->setObjSense(info.direction);
1455      solver->setDblParam(OsiObjOffset,info.offset);
1456      // Set integer variables
1457      for (int i=info.numberColumns-info.numberBinary-info.numberIntegers;
1458           i<info.numberColumns;i++)
1459        solver->setInteger(i);
1460      goodModel=true;
1461      // change argc etc
1462      argc = info.numberArguments;
1463      argv = const_cast<const char **>(info.arguments);
1464    }
1465#endif   
1466    // default action on import
1467    int allowImportErrors=0;
1468    int keepImportNames=1;
1469    int doIdiot=-1;
1470    int outputFormat=2;
1471    int slpValue=-1;
1472    int cppValue=-1;
1473    int printOptions=0;
1474    int printMode=0;
1475    int presolveOptions=0;
1476    int substitution=3;
1477    int dualize=0;
1478    int doCrash=0;
1479    int doVector=0;
1480    int doSprint=-1;
1481    int doScaling=1;
1482    // set reasonable defaults
1483    int preSolve=5;
1484    int preProcess=1;
1485    bool useStrategy=false;
1486    bool preSolveFile=false;
1487   
1488    double djFix=1.0e100;
1489    double gapRatio=1.0e100;
1490    double tightenFactor=0.0;
1491    lpSolver->setPerturbation(50);
1492    lpSolver->messageHandler()->setPrefix(false);
1493    const char dirsep =  CoinFindDirSeparator();
1494    std::string directory = (dirsep == '/' ? "./" : ".\\");
1495    std::string defaultDirectory = directory;
1496    std::string importFile ="";
1497    std::string exportFile ="default.mps";
1498    std::string importBasisFile ="";
1499    std::string importPriorityFile ="";
1500    std::string debugFile="";
1501    std::string printMask="";
1502    double * debugValues = NULL;
1503    int numberDebugValues = -1;
1504    int basisHasValues=0;
1505    std::string exportBasisFile ="default.bas";
1506    std::string saveFile ="default.prob";
1507    std::string restoreFile ="default.prob";
1508    std::string solutionFile ="stdout";
1509    std::string solutionSaveFile ="solution.file";
1510#define CBCMAXPARAMETERS 200
1511    CbcOrClpParam parameters[CBCMAXPARAMETERS];
1512    int numberParameters ;
1513    establishParams(numberParameters,parameters) ;
1514    parameters[whichParam(BASISIN,numberParameters,parameters)].setStringValue(importBasisFile);
1515    parameters[whichParam(PRIORITYIN,numberParameters,parameters)].setStringValue(importPriorityFile);
1516    parameters[whichParam(BASISOUT,numberParameters,parameters)].setStringValue(exportBasisFile);
1517    parameters[whichParam(DEBUG,numberParameters,parameters)].setStringValue(debugFile);
1518    parameters[whichParam(PRINTMASK,numberParameters,parameters)].setStringValue(printMask);
1519    parameters[whichParam(DIRECTORY,numberParameters,parameters)].setStringValue(directory);
1520    parameters[whichParam(DUALBOUND,numberParameters,parameters)].setDoubleValue(lpSolver->dualBound());
1521    parameters[whichParam(DUALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->dualTolerance());
1522    parameters[whichParam(EXPORT,numberParameters,parameters)].setStringValue(exportFile);
1523    parameters[whichParam(IDIOT,numberParameters,parameters)].setIntValue(doIdiot);
1524    parameters[whichParam(IMPORT,numberParameters,parameters)].setStringValue(importFile);
1525    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].setDoubleValue(1.0e-8);
1526    int slog = whichParam(SOLVERLOGLEVEL,numberParameters,parameters);
1527    int log = whichParam(LOGLEVEL,numberParameters,parameters);
1528    parameters[slog].setIntValue(0);
1529    parameters[log].setIntValue(1);
1530    parameters[whichParam(MAXFACTOR,numberParameters,parameters)].setIntValue(lpSolver->factorizationFrequency());
1531    parameters[whichParam(MAXITERATION,numberParameters,parameters)].setIntValue(lpSolver->maximumIterations());
1532    parameters[whichParam(OUTPUTFORMAT,numberParameters,parameters)].setIntValue(outputFormat);
1533    parameters[whichParam(PRESOLVEPASS,numberParameters,parameters)].setIntValue(preSolve);
1534    parameters[whichParam(PERTVALUE,numberParameters,parameters)].setIntValue(lpSolver->perturbation());
1535    parameters[whichParam(PRIMALTOLERANCE,numberParameters,parameters)].setDoubleValue(lpSolver->primalTolerance());
1536    parameters[whichParam(PRIMALWEIGHT,numberParameters,parameters)].setDoubleValue(lpSolver->infeasibilityCost());
1537    parameters[whichParam(RESTORE,numberParameters,parameters)].setStringValue(restoreFile);
1538    parameters[whichParam(SAVE,numberParameters,parameters)].setStringValue(saveFile);
1539    //parameters[whichParam(TIMELIMIT,numberParameters,parameters)].setDoubleValue(1.0e8);
1540    parameters[whichParam(TIMELIMIT_BAB,numberParameters,parameters)].setDoubleValue(1.0e8);
1541    parameters[whichParam(SOLUTION,numberParameters,parameters)].setStringValue(solutionFile);
1542    parameters[whichParam(SAVESOL,numberParameters,parameters)].setStringValue(solutionSaveFile);
1543    parameters[whichParam(SPRINT,numberParameters,parameters)].setIntValue(doSprint);
1544    parameters[whichParam(SUBSTITUTION,numberParameters,parameters)].setIntValue(substitution);
1545    parameters[whichParam(DUALIZE,numberParameters,parameters)].setIntValue(dualize);
1546    model.setNumberBeforeTrust(5);
1547    parameters[whichParam(NUMBERBEFORE,numberParameters,parameters)].setIntValue(5);
1548    parameters[whichParam(MAXNODES,numberParameters,parameters)].setIntValue(model.getMaximumNodes());
1549    model.setNumberStrong(5);
1550    parameters[whichParam(STRONGBRANCHING,numberParameters,parameters)].setIntValue(model.numberStrong());
1551    parameters[whichParam(INFEASIBILITYWEIGHT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcInfeasibilityWeight));
1552    parameters[whichParam(INTEGERTOLERANCE,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcIntegerTolerance));
1553    double normalIncrement=model.getCutoffIncrement();;
1554    parameters[whichParam(INCREMENT,numberParameters,parameters)].setDoubleValue(model.getDblParam(CbcModel::CbcCutoffIncrement));
1555    parameters[whichParam(TESTOSI,numberParameters,parameters)].setIntValue(testOsiParameters);
1556    parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].setIntValue(1003);
1557    if (testOsiParameters>=0) {
1558      // trying nonlinear - switch off some stuff
1559      preProcess=0;
1560    }
1561    // Set up likely cut generators and defaults
1562    parameters[whichParam(PREPROCESS,numberParameters,parameters)].setCurrentOption("on");
1563    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(128|64|1);
1564    parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].setIntValue(1);
1565    parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].setIntValue(-1);
1566    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].setIntValue(100);
1567    parameters[whichParam(CUTSSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
1568    parameters[whichParam(HEURISTICSTRATEGY,numberParameters,parameters)].setCurrentOption("on");
1569    parameters[whichParam(NODESTRATEGY,numberParameters,parameters)].setCurrentOption("fewest");
1570    int nodeStrategy=0;
1571    int doSOS=1;
1572    int verbose=0;
1573    CglGomory gomoryGen;
1574    // try larger limit
1575    gomoryGen.setLimitAtRoot(512);
1576    gomoryGen.setLimit(50);
1577    // set default action (0=off,1=on,2=root)
1578    int gomoryAction=3;
1579    parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1580
1581    CglProbing probingGen;
1582    probingGen.setUsingObjective(1);
1583    probingGen.setMaxPass(3);
1584    probingGen.setMaxPassRoot(3);
1585    // Number of unsatisfied variables to look at
1586    probingGen.setMaxProbe(10);
1587    probingGen.setMaxProbeRoot(50);
1588    // How far to follow the consequences
1589    probingGen.setMaxLook(10);
1590    probingGen.setMaxLookRoot(50);
1591    probingGen.setMaxLookRoot(10);
1592    // Only look at rows with fewer than this number of elements
1593    probingGen.setMaxElements(200);
1594    probingGen.setRowCuts(3);
1595    // set default action (0=off,1=on,2=root)
1596    int probingAction=1;
1597    parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1598
1599    CglKnapsackCover knapsackGen;
1600    //knapsackGen.switchOnExpensive();
1601    // set default action (0=off,1=on,2=root)
1602    int knapsackAction=3;
1603    parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1604
1605    CglRedSplit redsplitGen;
1606    //redsplitGen.setLimit(100);
1607    // set default action (0=off,1=on,2=root)
1608    // Off as seems to give some bad cuts
1609    int redsplitAction=0;
1610    parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption("off");
1611
1612    CglClique cliqueGen(false,true);
1613    cliqueGen.setStarCliqueReport(false);
1614    cliqueGen.setRowCliqueReport(false);
1615    cliqueGen.setMinViolation(0.1);
1616    // set default action (0=off,1=on,2=root)
1617    int cliqueAction=3;
1618    parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1619
1620    CglMixedIntegerRounding2 mixedGen;
1621    // set default action (0=off,1=on,2=root)
1622    int mixedAction=3;
1623    parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1624
1625    CglFlowCover flowGen;
1626    // set default action (0=off,1=on,2=root)
1627    int flowAction=3;
1628    parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption("ifmove");
1629
1630    CglTwomir twomirGen;
1631    twomirGen.setMaxElements(250);
1632    // set default action (0=off,1=on,2=root)
1633    int twomirAction=2;
1634    parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption("root");
1635    CglLandP landpGen;
1636    // set default action (0=off,1=on,2=root)
1637    int landpAction=0;
1638    parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption("off");
1639    // Stored cuts
1640    bool storedCuts = false;
1641
1642    bool useRounding=true;
1643    parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption("on");
1644    bool useFpump=true;
1645    parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption("on");
1646    bool useGreedy=true;
1647    parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption("on");
1648    bool useCombine=true;
1649    parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption("on");
1650    bool useLocalTree=false;
1651    bool useRINS=false;
1652    parameters[whichParam(RINS,numberParameters,parameters)].setCurrentOption("off");
1653    parameters[whichParam(COSTSTRATEGY,numberParameters,parameters)].setCurrentOption("off");
1654    int useCosts=0;
1655    // don't use input solution
1656    int useSolution=0;
1657   
1658    // total number of commands read
1659    int numberGoodCommands=0;
1660    // Set false if user does anything advanced
1661    bool defaultSettings=true;
1662
1663    // Hidden stuff for barrier
1664    int choleskyType = 0;
1665    int gamma=0;
1666    int scaleBarrier=0;
1667    int doKKT=0;
1668    int crossover=2; // do crossover unless quadratic
1669    // For names
1670    int lengthName = 0;
1671    std::vector<std::string> rowNames;
1672    std::vector<std::string> columnNames;
1673   
1674    std::string field;
1675    if (!noPrinting) {
1676      std::cout<<"Coin Cbc and Clp Solver version "<<CBCVERSION
1677               <<", build "<<__DATE__<<std::endl;
1678      // Print command line
1679      if (argc>1) {
1680        printf("command line - ");
1681        for (int i=0;i<argc;i++)
1682          printf("%s ",argv[i]);
1683        printf("\n");
1684      }
1685    }
1686    while (1) {
1687      // next command
1688      field=CoinReadGetCommand(argc,argv);
1689      // adjust field if has odd trailing characters
1690      char temp [200];
1691      strcpy(temp,field.c_str());
1692      int length = strlen(temp);
1693      for (int k=length-1;k>=0;k--) {
1694        if (temp[k]<' ')
1695          length--;
1696        else
1697          break;
1698      }
1699      temp[length]='\0';
1700      field=temp;
1701      // exit if null or similar
1702      if (!field.length()) {
1703        if (numberGoodCommands==1&&goodModel) {
1704          // we just had file name - do branch and bound
1705          field="branch";
1706        } else if (!numberGoodCommands) {
1707          // let's give the sucker a hint
1708          std::cout
1709            <<"CoinSolver takes input from arguments ( - switches to stdin)"
1710            <<std::endl
1711            <<"Enter ? for list of commands or help"<<std::endl;
1712          field="-";
1713        } else {
1714          break;
1715        }
1716      }
1717     
1718      // see if ? at end
1719      int numberQuery=0;
1720      if (field!="?"&&field!="???") {
1721        int length = field.length();
1722        int i;
1723        for (i=length-1;i>0;i--) {
1724          if (field[i]=='?') 
1725            numberQuery++;
1726          else
1727            break;
1728        }
1729        field=field.substr(0,length-numberQuery);
1730      }
1731      // find out if valid command
1732      int iParam;
1733      int numberMatches=0;
1734      int firstMatch=-1;
1735      for ( iParam=0; iParam<numberParameters; iParam++ ) {
1736        int match = parameters[iParam].matches(field);
1737        if (match==1) {
1738          numberMatches = 1;
1739          firstMatch=iParam;
1740          break;
1741        } else {
1742          if (match&&firstMatch<0)
1743            firstMatch=iParam;
1744          numberMatches += match>>1;
1745        }
1746      }
1747      if (iParam<numberParameters&&!numberQuery) {
1748        // found
1749        CbcOrClpParam found = parameters[iParam];
1750        CbcOrClpParameterType type = found.type();
1751        int valid;
1752        numberGoodCommands++;
1753        if (type==BAB&&goodModel) {
1754          // check if any integers
1755#ifdef COIN_HAS_ASL
1756          if (info.numberSos&&doSOS&&usingAmpl) {
1757            // SOS
1758            numberSOS = info.numberSos;
1759          }
1760#endif
1761          if (!lpSolver->integerInformation()&&!numberSOS&&
1762              !clpSolver->numberSOS())
1763            type=DUALSIMPLEX;
1764        }
1765        if (type==GENERALQUERY) {
1766          bool evenHidden=false;
1767          if ((verbose&8)!=0) {
1768            // even hidden
1769            evenHidden = true;
1770            verbose &= ~8;
1771          }
1772#ifdef COIN_HAS_ASL
1773          if (verbose<4&&usingAmpl)
1774            verbose +=4;
1775#endif
1776          if (verbose<4) {
1777            std::cout<<"In argument list keywords have leading - "
1778              ", -stdin or just - switches to stdin"<<std::endl;
1779            std::cout<<"One command per line (and no -)"<<std::endl;
1780            std::cout<<"abcd? gives list of possibilities, if only one + explanation"<<std::endl;
1781            std::cout<<"abcd?? adds explanation, if only one fuller help"<<std::endl;
1782            std::cout<<"abcd without value (where expected) gives current value"<<std::endl;
1783            std::cout<<"abcd value sets value"<<std::endl;
1784            std::cout<<"Commands are:"<<std::endl;
1785          } else {
1786            std::cout<<"Cbc options are set within AMPL with commands like:"<<std::endl<<std::endl;
1787            std::cout<<"         option cbc_options \"cuts=root log=2 feas=on slog=1\""<<std::endl<<std::endl;
1788            std::cout<<"only maximize, dual, primal, help and quit are recognized without ="<<std::endl;
1789          }
1790          int maxAcross=5;
1791          if ((verbose%4)!=0)
1792            maxAcross=1;
1793          int limits[]={1,51,101,151,201,251,301,351,401};
1794          std::vector<std::string> types;
1795          types.push_back("Double parameters:");
1796          types.push_back("Branch and Cut double parameters:");
1797          types.push_back("Integer parameters:");
1798          types.push_back("Branch and Cut integer parameters:");
1799          types.push_back("Keyword parameters:");
1800          types.push_back("Branch and Cut keyword parameters:");
1801          types.push_back("Actions or string parameters:");
1802          types.push_back("Branch and Cut actions:");
1803          int iType;
1804          for (iType=0;iType<8;iType++) {
1805            int across=0;
1806            if ((verbose%4)!=0)
1807              std::cout<<std::endl;
1808            std::cout<<types[iType]<<std::endl;
1809            if ((verbose&2)!=0)
1810              std::cout<<std::endl;
1811            for ( iParam=0; iParam<numberParameters; iParam++ ) {
1812              int type = parameters[iParam].type();
1813              if ((parameters[iParam].displayThis()||evenHidden)&&
1814                  type>=limits[iType]
1815                  &&type<limits[iType+1]) {
1816                // but skip if not useful for ampl (and in ampl mode)
1817                if (verbose>=4&&(parameters[iParam].whereUsed()&4)==0)
1818                  continue;
1819                if (!across) {
1820                  if ((verbose&2)==0) 
1821                    std::cout<<"  ";
1822                  else
1823                    std::cout<<"Command ";
1824                }
1825                std::cout<<parameters[iParam].matchName()<<"  ";
1826                across++;
1827                if (across==maxAcross) {
1828                  across=0;
1829                  if ((verbose%4)!=0) {
1830                    // put out description as well
1831                    if ((verbose&1)!=0) 
1832                      std::cout<<parameters[iParam].shortHelp();
1833                    std::cout<<std::endl;
1834                    if ((verbose&2)!=0) {
1835                      std::cout<<"---- description"<<std::endl;
1836                      parameters[iParam].printLongHelp();
1837                      std::cout<<"----"<<std::endl<<std::endl;
1838                    }
1839                  } else {
1840                    std::cout<<std::endl;
1841                  }
1842                }
1843              }
1844            }
1845            if (across)
1846              std::cout<<std::endl;
1847          }
1848        } else if (type==FULLGENERALQUERY) {
1849          std::cout<<"Full list of commands is:"<<std::endl;
1850          int maxAcross=5;
1851          int limits[]={1,51,101,151,201,251,301,351,401};
1852          std::vector<std::string> types;
1853          types.push_back("Double parameters:");
1854          types.push_back("Branch and Cut double parameters:");
1855          types.push_back("Integer parameters:");
1856          types.push_back("Branch and Cut integer parameters:");
1857          types.push_back("Keyword parameters:");
1858          types.push_back("Branch and Cut keyword parameters:");
1859          types.push_back("Actions or string parameters:");
1860          types.push_back("Branch and Cut actions:");
1861          int iType;
1862          for (iType=0;iType<8;iType++) {
1863            int across=0;
1864            std::cout<<types[iType]<<"  ";
1865            for ( iParam=0; iParam<numberParameters; iParam++ ) {
1866              int type = parameters[iParam].type();
1867              if (type>=limits[iType]
1868                  &&type<limits[iType+1]) {
1869                if (!across)
1870                  std::cout<<"  ";
1871                std::cout<<parameters[iParam].matchName()<<"  ";
1872                across++;
1873                if (across==maxAcross) {
1874                  std::cout<<std::endl;
1875                  across=0;
1876                }
1877              }
1878            }
1879            if (across)
1880              std::cout<<std::endl;
1881          }
1882        } else if (type<101) {
1883          // get next field as double
1884          double value = CoinReadGetDoubleField(argc,argv,&valid);
1885          if (!valid) {
1886            if (type<51) {
1887              parameters[iParam].setDoubleParameter(lpSolver,value);
1888            } else if (type<81) {
1889              parameters[iParam].setDoubleParameter(model,value);
1890            } else {
1891              parameters[iParam].setDoubleParameter(lpSolver,value);
1892              switch(type) {
1893              case DJFIX:
1894                djFix=value;
1895                preSolve=5;
1896                defaultSettings=false; // user knows what she is doing
1897                break;
1898              case GAPRATIO:
1899                gapRatio=value;
1900                break;
1901              case TIGHTENFACTOR:
1902                tightenFactor=value;
1903                if(!complicatedInteger)
1904                  defaultSettings=false; // user knows what she is doing
1905                break;
1906              default:
1907                abort();
1908              }
1909            }
1910          } else if (valid==1) {
1911            abort();
1912          } else {
1913            std::cout<<parameters[iParam].name()<<" has value "<<
1914              parameters[iParam].doubleValue()<<std::endl;
1915          }
1916        } else if (type<201) {
1917          // get next field as int
1918          int value = CoinReadGetIntField(argc,argv,&valid);
1919          if (!valid) {
1920            if (type<151) {
1921              if (parameters[iParam].type()==PRESOLVEPASS)
1922                preSolve = value;
1923              else if (parameters[iParam].type()==IDIOT)
1924                doIdiot = value;
1925              else if (parameters[iParam].type()==SPRINT)
1926                doSprint = value;
1927              else if (parameters[iParam].type()==OUTPUTFORMAT)
1928                outputFormat = value;
1929              else if (parameters[iParam].type()==SLPVALUE)
1930                slpValue = value;
1931              else if (parameters[iParam].type()==CPP)
1932                cppValue = value;
1933              else if (parameters[iParam].type()==PRESOLVEOPTIONS)
1934                presolveOptions = value;
1935              else if (parameters[iParam].type()==PRINTOPTIONS)
1936                printOptions = value;
1937              else if (parameters[iParam].type()==SUBSTITUTION)
1938                substitution = value;
1939              else if (parameters[iParam].type()==DUALIZE)
1940                dualize = value;
1941              else if (parameters[iParam].type()==CUTPASS)
1942                cutPass = value;
1943              else if (parameters[iParam].type()==PROCESSTUNE)
1944                tunePreProcess = value;
1945              else if (parameters[iParam].type()==VERBOSE)
1946                verbose = value;
1947              else if (parameters[iParam].type()==FPUMPITS)
1948                { useFpump = true;parameters[iParam].setIntValue(value);}
1949              parameters[iParam].setIntParameter(lpSolver,value);
1950            } else {
1951              parameters[iParam].setIntParameter(model,value);
1952            }
1953          } else if (valid==1) {
1954            abort();
1955          } else {
1956            std::cout<<parameters[iParam].name()<<" has value "<<
1957              parameters[iParam].intValue()<<std::endl;
1958          }
1959        } else if (type<301) {
1960          // one of several strings
1961          std::string value = CoinReadGetString(argc,argv);
1962          int action = parameters[iParam].parameterOption(value);
1963          if (action<0) {
1964            if (value!="EOL") {
1965              // no match
1966              parameters[iParam].printOptions();
1967            } else {
1968              // print current value
1969              std::cout<<parameters[iParam].name()<<" has value "<<
1970                parameters[iParam].currentOption()<<std::endl;
1971            }
1972          } else {
1973            parameters[iParam].setCurrentOption(action,!noPrinting);
1974            // for now hard wired
1975            switch (type) {
1976            case DIRECTION:
1977              if (action==0)
1978                lpSolver->setOptimizationDirection(1);
1979              else if (action==1)
1980                lpSolver->setOptimizationDirection(-1);
1981              else
1982                lpSolver->setOptimizationDirection(0);
1983              break;
1984            case DUALPIVOT:
1985              if (action==0) {
1986                ClpDualRowSteepest steep(3);
1987                lpSolver->setDualRowPivotAlgorithm(steep);
1988              } else if (action==1) {
1989                //ClpDualRowDantzig dantzig;
1990                ClpDualRowSteepest dantzig(5);
1991                lpSolver->setDualRowPivotAlgorithm(dantzig);
1992              } else if (action==2) {
1993                // partial steep
1994                ClpDualRowSteepest steep(2);
1995                lpSolver->setDualRowPivotAlgorithm(steep);
1996              } else {
1997                ClpDualRowSteepest steep;
1998                lpSolver->setDualRowPivotAlgorithm(steep);
1999              }
2000              break;
2001            case PRIMALPIVOT:
2002              if (action==0) {
2003                ClpPrimalColumnSteepest steep(3);
2004                lpSolver->setPrimalColumnPivotAlgorithm(steep);
2005              } else if (action==1) {
2006                ClpPrimalColumnSteepest steep(0);
2007                lpSolver->setPrimalColumnPivotAlgorithm(steep);
2008              } else if (action==2) {
2009                ClpPrimalColumnDantzig dantzig;
2010                lpSolver->setPrimalColumnPivotAlgorithm(dantzig);
2011              } else if (action==3) {
2012                ClpPrimalColumnSteepest steep(2);
2013                lpSolver->setPrimalColumnPivotAlgorithm(steep);
2014              } else if (action==4) {
2015                ClpPrimalColumnSteepest steep(1);
2016                lpSolver->setPrimalColumnPivotAlgorithm(steep);
2017              } else if (action==5) {
2018                ClpPrimalColumnSteepest steep(4);
2019                lpSolver->setPrimalColumnPivotAlgorithm(steep);
2020              } else if (action==6) {
2021                ClpPrimalColumnSteepest steep(10);
2022                lpSolver->setPrimalColumnPivotAlgorithm(steep);
2023              }
2024              break;
2025            case SCALING:
2026              lpSolver->scaling(action);
2027              solver->setHintParam(OsiDoScale,action!=0,OsiHintTry);
2028              doScaling = 1-action;
2029              break;
2030            case AUTOSCALE:
2031              lpSolver->setAutomaticScaling(action!=0);
2032              break;
2033            case SPARSEFACTOR:
2034              lpSolver->setSparseFactorization((1-action)!=0);
2035              break;
2036            case BIASLU:
2037              lpSolver->factorization()->setBiasLU(action);
2038              break;
2039            case PERTURBATION:
2040              if (action==0)
2041                lpSolver->setPerturbation(50);
2042              else
2043                lpSolver->setPerturbation(100);
2044              break;
2045            case ERRORSALLOWED:
2046              allowImportErrors = action;
2047              break;
2048            case INTPRINT:
2049              printMode=action;
2050              break;
2051              //case ALGORITHM:
2052              //algorithm  = action;
2053              //defaultSettings=false; // user knows what she is doing
2054              //abort();
2055              //break;
2056            case KEEPNAMES:
2057              keepImportNames = 1-action;
2058              break;
2059            case PRESOLVE:
2060              if (action==0)
2061                preSolve = 5;
2062              else if (action==1)
2063                preSolve=0;
2064              else if (action==2)
2065                preSolve=10;
2066              else
2067                preSolveFile=true;
2068              break;
2069            case PFI:
2070              lpSolver->factorization()->setForrestTomlin(action==0);
2071              break;
2072            case CRASH:
2073              doCrash=action;
2074              break;
2075            case VECTOR:
2076              doVector=action;
2077              break;
2078            case MESSAGES:
2079              lpSolver->messageHandler()->setPrefix(action!=0);
2080              break;
2081            case CHOLESKY:
2082              choleskyType = action;
2083              break;
2084            case GAMMA:
2085              gamma=action;
2086              break;
2087            case BARRIERSCALE:
2088              scaleBarrier=action;
2089              break;
2090            case KKT:
2091              doKKT=action;
2092              break;
2093            case CROSSOVER:
2094              crossover=action;
2095              break;
2096            case SOS:
2097              doSOS=action;
2098              break;
2099            case GOMORYCUTS:
2100              defaultSettings=false; // user knows what she is doing
2101              gomoryAction = action;
2102              break;
2103            case PROBINGCUTS:
2104              defaultSettings=false; // user knows what she is doing
2105              probingAction = action;
2106              break;
2107            case KNAPSACKCUTS:
2108              defaultSettings=false; // user knows what she is doing
2109              knapsackAction = action;
2110              break;
2111            case REDSPLITCUTS:
2112              defaultSettings=false; // user knows what she is doing
2113              redsplitAction = action;
2114              break;
2115            case CLIQUECUTS:
2116              defaultSettings=false; // user knows what she is doing
2117              cliqueAction = action;
2118              break;
2119            case FLOWCUTS:
2120              defaultSettings=false; // user knows what she is doing
2121              flowAction = action;
2122              break;
2123            case MIXEDCUTS:
2124              defaultSettings=false; // user knows what she is doing
2125              mixedAction = action;
2126              break;
2127            case TWOMIRCUTS:
2128              defaultSettings=false; // user knows what she is doing
2129              twomirAction = action;
2130              break;
2131            case LANDPCUTS:
2132              defaultSettings=false; // user knows what she is doing
2133              landpAction = action;
2134              break;
2135            case ROUNDING:
2136              defaultSettings=false; // user knows what she is doing
2137              useRounding = action;
2138              break;
2139            case FPUMP:
2140              defaultSettings=false; // user knows what she is doing
2141              useFpump=action;
2142              break;
2143            case RINS:
2144              useRINS=action;
2145              break;
2146            case CUTSSTRATEGY:
2147              gomoryAction = action;
2148              probingAction = action;
2149              knapsackAction = action;
2150              cliqueAction = action;
2151              flowAction = action;
2152              mixedAction = action;
2153              twomirAction = action;
2154              //landpAction = action;
2155              parameters[whichParam(GOMORYCUTS,numberParameters,parameters)].setCurrentOption(action);
2156              parameters[whichParam(PROBINGCUTS,numberParameters,parameters)].setCurrentOption(action);
2157              parameters[whichParam(KNAPSACKCUTS,numberParameters,parameters)].setCurrentOption(action);
2158              if (!action) {
2159                redsplitAction = action;
2160                parameters[whichParam(REDSPLITCUTS,numberParameters,parameters)].setCurrentOption(action);
2161              }
2162              parameters[whichParam(CLIQUECUTS,numberParameters,parameters)].setCurrentOption(action);
2163              parameters[whichParam(FLOWCUTS,numberParameters,parameters)].setCurrentOption(action);
2164              parameters[whichParam(MIXEDCUTS,numberParameters,parameters)].setCurrentOption(action);
2165              parameters[whichParam(TWOMIRCUTS,numberParameters,parameters)].setCurrentOption(action);
2166              if (!action) {
2167                landpAction = action;
2168                parameters[whichParam(LANDPCUTS,numberParameters,parameters)].setCurrentOption(action);
2169              }
2170              break;
2171            case HEURISTICSTRATEGY:
2172              useRounding = action;
2173              useGreedy = action;
2174              useCombine = action;
2175              //useLocalTree = action;
2176              useFpump=action;
2177              parameters[whichParam(ROUNDING,numberParameters,parameters)].setCurrentOption(action);
2178              parameters[whichParam(GREEDY,numberParameters,parameters)].setCurrentOption(action);
2179              parameters[whichParam(COMBINE,numberParameters,parameters)].setCurrentOption(action);
2180              //parameters[whichParam(LOCALTREE,numberParameters,parameters)].setCurrentOption(action);
2181              parameters[whichParam(FPUMP,numberParameters,parameters)].setCurrentOption(action);
2182              break;
2183            case GREEDY:
2184              defaultSettings=false; // user knows what she is doing
2185              useGreedy = action;
2186              break;
2187            case COMBINE:
2188              defaultSettings=false; // user knows what she is doing
2189              useCombine = action;
2190              break;
2191            case LOCALTREE:
2192              defaultSettings=false; // user knows what she is doing
2193              useLocalTree = action;
2194              break;
2195            case COSTSTRATEGY:
2196              useCosts=action;
2197              break;
2198            case NODESTRATEGY:
2199              nodeStrategy=action;
2200              break;
2201            case PREPROCESS:
2202              preProcess = action;
2203              break;
2204            case USESOLUTION:
2205              useSolution = action;
2206              break;
2207            default:
2208              abort();
2209            }
2210          }
2211        } else {
2212          // action
2213          if (type==EXIT) {
2214#ifdef COIN_HAS_ASL
2215            if(usingAmpl) {
2216              if (info.numberIntegers||info.numberBinary) {
2217                // integer
2218              } else {
2219                // linear
2220              }
2221              writeAmpl(&info);
2222              freeArrays2(&info);
2223              freeArgs(&info);
2224            }
2225#endif
2226            break; // stop all
2227          }
2228          switch (type) {
2229          case DUALSIMPLEX:
2230          case PRIMALSIMPLEX:
2231          case SOLVECONTINUOUS:
2232          case BARRIER:
2233            if (goodModel) {
2234              double objScale = 
2235                parameters[whichParam(OBJSCALE2,numberParameters,parameters)].doubleValue();
2236              if (objScale!=1.0) {
2237                int iColumn;
2238                int numberColumns=lpSolver->numberColumns();
2239                double * dualColumnSolution = 
2240                  lpSolver->dualColumnSolution();
2241                ClpObjective * obj = lpSolver->objectiveAsObject();
2242                assert(dynamic_cast<ClpLinearObjective *> (obj));
2243                double offset;
2244                double * objective = obj->gradient(NULL,NULL,offset,true);
2245                for (iColumn=0;iColumn<numberColumns;iColumn++) {
2246                  dualColumnSolution[iColumn] *= objScale;
2247                  objective[iColumn] *= objScale;;
2248                }
2249                int iRow;
2250                int numberRows=lpSolver->numberRows();
2251                double * dualRowSolution = 
2252                  lpSolver->dualRowSolution();
2253                for (iRow=0;iRow<numberRows;iRow++) 
2254                  dualRowSolution[iRow] *= objScale;
2255                lpSolver->setObjectiveOffset(objScale*lpSolver->objectiveOffset());
2256              }
2257              ClpSolve::SolveType method;
2258              ClpSolve::PresolveType presolveType;
2259              ClpSimplex * model2 = lpSolver;
2260              if (dualize) {
2261                model2 = ((ClpSimplexOther *) model2)->dualOfModel();
2262                printf("Dual of model has %d rows and %d columns\n",
2263                       model2->numberRows(),model2->numberColumns());
2264                model2->setOptimizationDirection(1.0);
2265              }
2266              if (noPrinting)
2267                lpSolver->setLogLevel(0);
2268              ClpSolve solveOptions;
2269              solveOptions.setPresolveActions(presolveOptions);
2270              solveOptions.setSubstitution(substitution);
2271              if (preSolve!=5&&preSolve) {
2272                presolveType=ClpSolve::presolveNumber;
2273                if (preSolve<0) {
2274                  preSolve = - preSolve;
2275                  if (preSolve<=100) {
2276                    presolveType=ClpSolve::presolveNumber;
2277                    printf("Doing %d presolve passes - picking up non-costed slacks\n",
2278                           preSolve);
2279                    solveOptions.setDoSingletonColumn(true);
2280                  } else {
2281                    preSolve -=100;
2282                    presolveType=ClpSolve::presolveNumberCost;
2283                    printf("Doing %d presolve passes - picking up costed slacks\n",
2284                           preSolve);
2285                  }
2286                } 
2287              } else if (preSolve) {
2288                presolveType=ClpSolve::presolveOn;
2289              } else {
2290                presolveType=ClpSolve::presolveOff;
2291              }
2292              solveOptions.setPresolveType(presolveType,preSolve);
2293              if (type==DUALSIMPLEX||type==SOLVECONTINUOUS) {
2294                method=ClpSolve::useDual;
2295              } else if (type==PRIMALSIMPLEX) {
2296                method=ClpSolve::usePrimalorSprint;
2297              } else {
2298                method = ClpSolve::useBarrier;
2299                if (crossover==1) {
2300                  method=ClpSolve::useBarrierNoCross;
2301                } else if (crossover==2) {
2302                  ClpObjective * obj = lpSolver->objectiveAsObject();
2303                  if (obj->type()>1) {
2304                    method=ClpSolve::useBarrierNoCross;
2305                    presolveType=ClpSolve::presolveOff;
2306                    solveOptions.setPresolveType(presolveType,preSolve);
2307                  } 
2308                }
2309              }
2310              solveOptions.setSolveType(method);
2311              if(preSolveFile)
2312                presolveOptions |= 0x40000000;
2313              solveOptions.setSpecialOption(4,presolveOptions);
2314              solveOptions.setSpecialOption(5,printOptions);
2315              if (doVector) {
2316                ClpMatrixBase * matrix = lpSolver->clpMatrix();
2317                if (dynamic_cast< ClpPackedMatrix*>(matrix)) {
2318                  ClpPackedMatrix * clpMatrix = dynamic_cast< ClpPackedMatrix*>(matrix);
2319                  clpMatrix->makeSpecialColumnCopy();
2320                }
2321              }
2322              if (method==ClpSolve::useDual) {
2323                // dual
2324                if (doCrash)
2325                  solveOptions.setSpecialOption(0,1,doCrash); // crash
2326                else if (doIdiot)
2327                  solveOptions.setSpecialOption(0,2,doIdiot); // possible idiot
2328              } else if (method==ClpSolve::usePrimalorSprint) {
2329                // primal
2330                // if slp turn everything off
2331                if (slpValue>0) {
2332                  doCrash=false;
2333                  doSprint=0;
2334                  doIdiot=-1;
2335                  solveOptions.setSpecialOption(1,10,slpValue); // slp
2336                  method=ClpSolve::usePrimal;
2337                }
2338                if (doCrash) {
2339                  solveOptions.setSpecialOption(1,1,doCrash); // crash
2340                } else if (doSprint>0) {
2341                  // sprint overrides idiot
2342                  solveOptions.setSpecialOption(1,3,doSprint); // sprint
2343                } else if (doIdiot>0) {
2344                  solveOptions.setSpecialOption(1,2,doIdiot); // idiot
2345                } else if (slpValue<=0) {
2346                  if (doIdiot==0) {
2347                    if (doSprint==0)
2348                      solveOptions.setSpecialOption(1,4); // all slack
2349                    else
2350                      solveOptions.setSpecialOption(1,9); // all slack or sprint
2351                  } else {
2352                    if (doSprint==0)
2353                      solveOptions.setSpecialOption(1,8); // all slack or idiot
2354                    else
2355                      solveOptions.setSpecialOption(1,7); // initiative
2356                  }
2357                }
2358                if (basisHasValues==-1)
2359                  solveOptions.setSpecialOption(1,11); // switch off values
2360              } else if (method==ClpSolve::useBarrier||method==ClpSolve::useBarrierNoCross) {
2361                int barrierOptions = choleskyType;
2362                if (scaleBarrier)
2363                  barrierOptions |= 8;
2364                if (doKKT)
2365                  barrierOptions |= 16;
2366                if (gamma)
2367                  barrierOptions |= 32*gamma;
2368                if (crossover==3) 
2369                  barrierOptions |= 256; // try presolve in crossover
2370                solveOptions.setSpecialOption(4,barrierOptions);
2371              }
2372#ifdef COIN_HAS_LINK
2373              OsiSolverInterface * coinSolver = model.solver();
2374              OsiSolverLink * linkSolver = dynamic_cast< OsiSolverLink*> (coinSolver);
2375              if (!linkSolver) {
2376                model2->initialSolve(solveOptions);
2377              } else {
2378                // special solver
2379                double * solution = linkSolver->nonlinearSLP(slpValue,1.0e-5);
2380                if (solution) {
2381                  memcpy(model2->primalColumnSolution(),solution,
2382                         CoinMin(model2->numberColumns(),linkSolver->coinModel()->numberColumns())*sizeof(double));
2383                  delete [] solution;
2384                } else {
2385                  printf("No nonlinear solution\n");
2386                }
2387              }
2388#else
2389              model2->initialSolve(solveOptions);
2390#endif
2391              basisHasValues=1;
2392              if (dualize) {
2393                int returnCode=((ClpSimplexOther *) lpSolver)->restoreFromDual(model2);
2394                delete model2;
2395                if (returnCode)
2396                  lpSolver->primal(1);
2397                model2=lpSolver;
2398              }
2399#ifdef COIN_HAS_ASL
2400              if (usingAmpl) {
2401                double value = model2->getObjValue()*model2->getObjSense();
2402                char buf[300];
2403                int pos=0;
2404                int iStat = model2->status();
2405                if (iStat==0) {
2406                  pos += sprintf(buf+pos,"optimal," );
2407                } else if (iStat==1) {
2408                  // infeasible
2409                  pos += sprintf(buf+pos,"infeasible,");
2410                } else if (iStat==2) {
2411                  // unbounded
2412                  pos += sprintf(buf+pos,"unbounded,");
2413                } else if (iStat==3) {
2414                  pos += sprintf(buf+pos,"stopped on iterations or time,");
2415                } else if (iStat==4) {
2416                  iStat = 7;
2417                  pos += sprintf(buf+pos,"stopped on difficulties,");
2418                } else if (iStat==5) {
2419                  iStat = 3;
2420                  pos += sprintf(buf+pos,"stopped on ctrl-c,");
2421                } else {
2422                  pos += sprintf(buf+pos,"status unknown,");
2423                  iStat=6;
2424                }
2425                info.problemStatus=iStat;
2426                info.objValue = value;
2427                pos += sprintf(buf+pos," objective %.*g",ampl_obj_prec(),
2428                               value);
2429                sprintf(buf+pos,"\n%d iterations",
2430                        model2->getIterationCount());
2431                free(info.primalSolution);
2432                int numberColumns=model2->numberColumns();
2433                info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
2434                CoinCopyN(model2->primalColumnSolution(),numberColumns,info.primalSolution);
2435                int numberRows = model2->numberRows();
2436                free(info.dualSolution);
2437                info.dualSolution = (double *) malloc(numberRows*sizeof(double));
2438                CoinCopyN(model2->dualRowSolution(),numberRows,info.dualSolution);
2439                CoinWarmStartBasis * basis = model2->getBasis();
2440                free(info.rowStatus);
2441                info.rowStatus = (int *) malloc(numberRows*sizeof(int));
2442                free(info.columnStatus);
2443                info.columnStatus = (int *) malloc(numberColumns*sizeof(int));
2444                // Put basis in
2445                int i;
2446                // free,basic,ub,lb are 0,1,2,3
2447                for (i=0;i<numberRows;i++) {
2448                  CoinWarmStartBasis::Status status = basis->getArtifStatus(i);
2449                  info.rowStatus[i]=status;
2450                }
2451                for (i=0;i<numberColumns;i++) {
2452                  CoinWarmStartBasis::Status status = basis->getStructStatus(i);
2453                  info.columnStatus[i]=status;
2454                }
2455                // put buffer into info
2456                strcpy(info.buffer,buf);
2457                delete basis;
2458              }
2459#endif
2460            } else {
2461              std::cout<<"** Current model not valid"<<std::endl;
2462            }
2463            break;
2464          case STATISTICS:
2465            if (goodModel) {
2466              // If presolve on look at presolved
2467              bool deleteModel2=false;
2468              ClpSimplex * model2 = lpSolver;
2469              if (preSolve) {
2470                ClpPresolve pinfo;
2471                int presolveOptions2 = presolveOptions&~0x40000000;
2472                if ((presolveOptions2&0xffff)!=0)
2473                  pinfo.setPresolveActions(presolveOptions2);
2474                pinfo.setSubstitution(substitution);
2475                if ((printOptions&1)!=0)
2476                  pinfo.statistics();
2477                double presolveTolerance = 
2478                  parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
2479                model2 = 
2480                  pinfo.presolvedModel(*lpSolver,presolveTolerance,
2481                                       true,preSolve);
2482                if (model2) {
2483                  printf("Statistics for presolved model\n");
2484                  deleteModel2=true;
2485                } else {
2486                  printf("Presolved model looks infeasible - will use unpresolved\n");
2487                  model2 = lpSolver;
2488                }
2489              } else {
2490                printf("Statistics for unpresolved model\n");
2491                model2 =  lpSolver;
2492              }
2493              statistics(lpSolver,model2);
2494              if (deleteModel2)
2495                delete model2;
2496            } else {
2497              std::cout<<"** Current model not valid"<<std::endl;
2498            }
2499            break;
2500          case TIGHTEN:
2501            if (goodModel) {
2502              int numberInfeasibilities = lpSolver->tightenPrimalBounds();
2503              if (numberInfeasibilities)
2504                std::cout<<"** Analysis indicates model infeasible"<<std::endl;
2505            } else {
2506              std::cout<<"** Current model not valid"<<std::endl;
2507            }
2508            break;
2509          case PLUSMINUS:
2510            if (goodModel) {
2511              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
2512              ClpPackedMatrix* clpMatrix =
2513                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
2514              if (clpMatrix) {
2515                ClpPlusMinusOneMatrix * newMatrix = new ClpPlusMinusOneMatrix(*(clpMatrix->matrix()));
2516                if (newMatrix->getIndices()) {
2517                  lpSolver->replaceMatrix(newMatrix);
2518                  delete saveMatrix;
2519                  std::cout<<"Matrix converted to +- one matrix"<<std::endl;
2520                } else {
2521                  std::cout<<"Matrix can not be converted to +- 1 matrix"<<std::endl;
2522                }
2523              } else {
2524                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
2525              }
2526            } else {
2527              std::cout<<"** Current model not valid"<<std::endl;
2528            }
2529            break;
2530          case OUTDUPROWS:
2531            if (goodModel) {
2532              int numberRows = clpSolver->getNumRows();
2533              //int nOut = outDupRow(clpSolver);
2534              CglDuplicateRow dupcuts(clpSolver);
2535              storedCuts = dupcuts.outDuplicates(clpSolver)!=0;
2536              int nOut = numberRows-clpSolver->getNumRows();
2537              if (nOut&&!noPrinting)
2538                printf("%d rows eliminated\n",nOut);
2539            } else {
2540              std::cout<<"** Current model not valid"<<std::endl;
2541            }
2542            break;
2543          case NETWORK:
2544            if (goodModel) {
2545              ClpMatrixBase * saveMatrix = lpSolver->clpMatrix();
2546              ClpPackedMatrix* clpMatrix =
2547                dynamic_cast< ClpPackedMatrix*>(saveMatrix);
2548              if (clpMatrix) {
2549                ClpNetworkMatrix * newMatrix = new ClpNetworkMatrix(*(clpMatrix->matrix()));
2550                if (newMatrix->getIndices()) {
2551                  lpSolver->replaceMatrix(newMatrix);
2552                  delete saveMatrix;
2553                  std::cout<<"Matrix converted to network matrix"<<std::endl;
2554                } else {
2555                  std::cout<<"Matrix can not be converted to network matrix"<<std::endl;
2556                }
2557              } else {
2558                std::cout<<"Matrix not a ClpPackedMatrix"<<std::endl;
2559              }
2560            } else {
2561              std::cout<<"** Current model not valid"<<std::endl;
2562            }
2563            break;
2564          case MIPLIB:
2565            // User can set options - main difference is lack of model and CglPreProcess
2566            goodModel=true;
2567/*
2568  Run branch-and-cut. First set a few options -- node comparison, scaling. If
2569  the solver is Clp, consider running some presolve code (not yet converted
2570  this to generic OSI) with branch-and-cut. If presolve is disabled, or the
2571  solver is not Clp, simply run branch-and-cut. Print elapsed time at the end.
2572*/
2573          case BAB: // branchAndBound
2574          case STRENGTHEN:
2575            if (goodModel) {
2576              bool miplib = type==MIPLIB;
2577              int logLevel = parameters[slog].intValue();
2578              // Reduce printout
2579              if (logLevel<=1)
2580                model.solver()->setHintParam(OsiDoReducePrint,true,OsiHintTry);
2581#ifdef COIN_HAS_LINK
2582              {
2583                OsiSolverInterface * solver = model.solver();
2584                OsiClpSolverInterface * si =
2585                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2586                assert (si != NULL);
2587                // See if quadratic
2588                if (!complicatedInteger) {
2589                  ClpSimplex * lpSolver = si->getModelPtr();
2590                  ClpQuadraticObjective * obj = (dynamic_cast< ClpQuadraticObjective*>(lpSolver->objectiveAsObject()));
2591                  if (obj) {
2592                    preProcess=0;
2593                    parameters[whichParam(TESTOSI,numberParameters,parameters)].setIntValue(0);
2594                    // create coin model
2595                    coinModel = lpSolver->createCoinModel();
2596                    assert (coinModel);
2597                    // load from coin model
2598                    OsiSolverLink solver1;
2599                    OsiSolverInterface * solver2 = solver1.clone();
2600                    model.assignSolver(solver2,true);
2601                    OsiSolverLink * si =
2602                      dynamic_cast<OsiSolverLink *>(model.solver()) ;
2603                    assert (si != NULL);
2604                    si->setDefaultMeshSize(0.001);
2605                    // need some relative granularity
2606                    si->setDefaultBound(100.0);
2607                    si->setDefaultMeshSize(0.01);
2608                    si->setDefaultBound(1000.0);
2609                    si->setIntegerPriority(1000);
2610                    si->setBiLinearPriority(10000);
2611                    si->setSpecialOptions2(2+4+8);
2612                    CoinModel * model2 = (CoinModel *) coinModel;
2613                    si->load(*model2,true,info.logLevel);
2614                    // redo
2615                    solver = model.solver();
2616                    clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
2617                    lpSolver = clpSolver->getModelPtr();
2618                    clpSolver->messageHandler()->setLogLevel(0) ;
2619                    testOsiParameters=0;
2620                    complicatedInteger=2;  // allow cuts
2621                    OsiSolverInterface * coinSolver = model.solver();
2622                    OsiSolverLink * linkSolver = dynamic_cast< OsiSolverLink*> (coinSolver);
2623                    if (linkSolver->quadraticModel()) {
2624                      ClpSimplex * qp = linkSolver->quadraticModel();
2625                      //linkSolver->nonlinearSLP(CoinMax(slpValue,10),1.0e-5);
2626                      qp->nonlinearSLP(CoinMax(slpValue,20),1.0e-5);
2627                      qp->primal(1);
2628                      OsiSolverLinearizedQuadratic solver2(qp);
2629                      const double * solution=NULL;
2630                      // Reduce printout
2631                      solver2.setHintParam(OsiDoReducePrint,true,OsiHintTry);
2632                      CbcModel model2(solver2);
2633                      // Now do requested saves and modifications
2634                      CbcModel * cbcModel = & model2;
2635                      OsiSolverInterface * osiModel = model2.solver();
2636                      OsiClpSolverInterface * osiclpModel = dynamic_cast< OsiClpSolverInterface*> (osiModel);
2637                      ClpSimplex * clpModel = osiclpModel->getModelPtr();
2638                     
2639                      // Set changed values
2640                     
2641                      CglProbing probing;
2642                      probing.setMaxProbe(10);
2643                      probing.setMaxLook(10);
2644                      probing.setMaxElements(200);
2645                      probing.setMaxProbeRoot(50);
2646                      probing.setMaxLookRoot(10);
2647                      probing.setRowCuts(3);
2648                      probing.setUsingObjective(true);
2649                      cbcModel->addCutGenerator(&probing,-1,"Probing",true,false,false,-100,-1,-1);
2650                      cbcModel->cutGenerator(0)->setTiming(true);
2651                     
2652                      CglGomory gomory;
2653                      gomory.setLimitAtRoot(512);
2654                      cbcModel->addCutGenerator(&gomory,-98,"Gomory",true,false,false,-100,-1,-1);
2655                      cbcModel->cutGenerator(1)->setTiming(true);
2656                     
2657                      CglKnapsackCover knapsackCover;
2658                      cbcModel->addCutGenerator(&knapsackCover,-98,"KnapsackCover",true,false,false,-100,-1,-1);
2659                      cbcModel->cutGenerator(2)->setTiming(true);
2660                     
2661                      CglRedSplit redSplit;
2662                      cbcModel->addCutGenerator(&redSplit,-99,"RedSplit",true,false,false,-100,-1,-1);
2663                      cbcModel->cutGenerator(3)->setTiming(true);
2664                     
2665                      CglClique clique;
2666                      clique.setStarCliqueReport(false);
2667                      clique.setRowCliqueReport(false);
2668                      clique.setMinViolation(0.1);
2669                      cbcModel->addCutGenerator(&clique,-98,"Clique",true,false,false,-100,-1,-1);
2670                      cbcModel->cutGenerator(4)->setTiming(true);
2671                     
2672                      CglMixedIntegerRounding2 mixedIntegerRounding2;
2673                      cbcModel->addCutGenerator(&mixedIntegerRounding2,-98,"MixedIntegerRounding2",true,false,false,-100,-1,-1);
2674                      cbcModel->cutGenerator(5)->setTiming(true);
2675                     
2676                      CglFlowCover flowCover;
2677                      cbcModel->addCutGenerator(&flowCover,-98,"FlowCover",true,false,false,-100,-1,-1);
2678                      cbcModel->cutGenerator(6)->setTiming(true);
2679                     
2680                      CglTwomir twomir;
2681                      twomir.setMaxElements(250);
2682                      cbcModel->addCutGenerator(&twomir,-99,"Twomir",true,false,false,-100,-1,-1);
2683                      cbcModel->cutGenerator(7)->setTiming(true);
2684                     
2685                      CbcHeuristicFPump heuristicFPump(*cbcModel);
2686                      heuristicFPump.setWhen(13);
2687                      heuristicFPump.setMaximumPasses(20);
2688                      heuristicFPump.setMaximumRetries(7);
2689                      heuristicFPump.setAbsoluteIncrement(4332.64);
2690                      cbcModel->addHeuristic(&heuristicFPump);
2691                      heuristicFPump.setInitialWeight(1);
2692                     
2693                      CbcRounding rounding(*cbcModel);
2694                      cbcModel->addHeuristic(&rounding);
2695                     
2696                      CbcHeuristicLocal heuristicLocal(*cbcModel);
2697                      heuristicLocal.setSearchType(1);
2698                      cbcModel->addHeuristic(&heuristicLocal);
2699                     
2700                      CbcHeuristicGreedyCover heuristicGreedyCover(*cbcModel);
2701                      cbcModel->addHeuristic(&heuristicGreedyCover);
2702                     
2703                      CbcHeuristicGreedyEquality heuristicGreedyEquality(*cbcModel);
2704                      cbcModel->addHeuristic(&heuristicGreedyEquality);
2705                     
2706                      CbcCompareDefault compare;
2707                      cbcModel->setNodeComparison(compare);
2708                      cbcModel->setNumberBeforeTrust(5);
2709                      cbcModel->setSpecialOptions(2);
2710                      cbcModel->messageHandler()->setLogLevel(1);
2711                      cbcModel->setMaximumCutPassesAtRoot(-100);
2712                      cbcModel->setMaximumCutPasses(1);
2713                      cbcModel->setMinimumDrop(0.05);
2714                      // For branchAndBound this may help
2715                      clpModel->defaultFactorizationFrequency();
2716                      clpModel->setDualBound(1.0001e+08);
2717                      clpModel->setPerturbation(50);
2718                      osiclpModel->setSpecialOptions(193);
2719                      osiclpModel->messageHandler()->setLogLevel(0);
2720                      osiclpModel->setIntParam(OsiMaxNumIterationHotStart,100);
2721                      osiclpModel->setHintParam(OsiDoReducePrint,true,OsiHintTry);
2722                      // You can save some time by switching off message building
2723                      // clpModel->messagesPointer()->setDetailMessages(100,10000,(int *) NULL);
2724                     
2725                      // Solve
2726                     
2727                      cbcModel->initialSolve();
2728                      if (clpModel->tightenPrimalBounds()!=0) {
2729                        std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
2730                        exit(1);
2731                      }
2732                      clpModel->dual();  // clean up
2733                      cbcModel->initialSolve();
2734                      cbcModel->branchAndBound();
2735                      OsiSolverLinearizedQuadratic * solver3 = dynamic_cast<OsiSolverLinearizedQuadratic *> (model2.solver());
2736                      assert (solver3);
2737                      solution = solver3->bestSolution();
2738                      double bestObjectiveValue = solver3->bestObjectiveValue();
2739                      linkSolver->setBestObjectiveValue(bestObjectiveValue);
2740                      linkSolver->setBestSolution(solution,solver3->getNumCols());
2741                      CbcHeuristicDynamic3 dynamic(model);
2742                      model.addHeuristic(&dynamic);
2743                      // if convex
2744                      if ((linkSolver->specialOptions2()&4)!=0) {
2745                        int numberColumns = coinModel->numberColumns();
2746                        // add OA cut
2747                        double offset;
2748                        double * gradient = new double [numberColumns+1];
2749                        memcpy(gradient,qp->objectiveAsObject()->gradient(qp,solution,offset,true,2),
2750                               numberColumns*sizeof(double));
2751                        double rhs = 0.0;
2752                        int * column = new int[numberColumns+1];
2753                        int n=0;
2754                        for (int i=0;i<numberColumns;i++) {
2755                          double value = gradient[i];
2756                          if (fabs(value)>1.0e-12) {
2757                            gradient[n]=value;
2758                            rhs += value*solution[i];
2759                            column[n++]=i;
2760                          }
2761                        }
2762                        gradient[n]=-1.0;
2763                        column[n++]=numberColumns;
2764                        storedAmpl.addCut(-COIN_DBL_MAX,offset+1.0e-7,n,column,gradient);
2765                        delete [] gradient;
2766                        delete [] column;
2767                      }
2768                      // could do three way branching round a) continuous b) best solution
2769                      printf("obj %g\n",bestObjectiveValue);
2770                      linkSolver->initialSolve();
2771                    }
2772                  }
2773                }
2774                si->setSpecialOptions(0x40000000);
2775              }
2776#endif
2777              if (!miplib) {
2778                if (!preSolve) {
2779                  model.solver()->setHintParam(OsiDoPresolveInInitial,false,OsiHintTry);
2780                  model.solver()->setHintParam(OsiDoPresolveInResolve,false,OsiHintTry);
2781                }
2782                model.initialSolve();
2783                OsiSolverInterface * solver = model.solver();
2784                OsiClpSolverInterface * si =
2785                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2786                ClpSimplex * clpSolver = si->getModelPtr();
2787                clpSolver->setSpecialOptions(clpSolver->specialOptions()|0x01000000); // say is Cbc (and in branch and bound)
2788                if (!complicatedInteger&&clpSolver->tightenPrimalBounds()!=0) {
2789                  std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
2790                  exit(1);
2791                }
2792                if (clpSolver->dualBound()==1.0e10) {
2793                  // user did not set - so modify
2794                  // get largest scaled away from bound
2795                  double largest=1.0e-12;
2796                  int numberRows = clpSolver->numberRows();
2797                  const double * rowPrimal = clpSolver->primalRowSolution();
2798                  const double * rowLower = clpSolver->rowLower();
2799                  const double * rowUpper = clpSolver->rowUpper();
2800                  const double * rowScale = clpSolver->rowScale();
2801                  int iRow;
2802                  for (iRow=0;iRow<numberRows;iRow++) {
2803                    double value = rowPrimal[iRow];
2804                    double above = value-rowLower[iRow];
2805                    double below = rowUpper[iRow]-value;
2806                    if (rowScale) {
2807                      double multiplier = rowScale[iRow];
2808                      above *= multiplier;
2809                      below *= multiplier;
2810                    }
2811                    if (above<1.0e12)
2812                      largest = CoinMax(largest,above);
2813                    if (below<1.0e12)
2814                      largest = CoinMax(largest,below);
2815                  }
2816                 
2817                  int numberColumns = clpSolver->numberColumns();
2818                  const double * columnPrimal = clpSolver->primalColumnSolution();
2819                  const double * columnLower = clpSolver->columnLower();
2820                  const double * columnUpper = clpSolver->columnUpper();
2821                  const double * columnScale = clpSolver->columnScale();
2822                  int iColumn;
2823                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2824                    double value = columnPrimal[iColumn];
2825                    double above = value-columnLower[iColumn];
2826                    double below = columnUpper[iColumn]-value;
2827                    if (columnScale) {
2828                      double multiplier = 1.0/columnScale[iColumn];
2829                      above *= multiplier;
2830                      below *= multiplier;
2831                    }
2832                    if (above<1.0e12)
2833                      largest = CoinMax(largest,above);
2834                    if (below<1.0e12)
2835                      largest = CoinMax(largest,below);
2836                  }
2837                  //if (!noPrinting)
2838                  //std::cout<<"Largest (scaled) away from bound "<<largest<<std::endl;
2839                  clpSolver->setDualBound(CoinMax(1.0001e8,CoinMin(1000.0*largest,1.00001e10)));
2840                }
2841                clpSolver->dual();  // clean up
2842              }
2843              // If user made settings then use them
2844              if (!defaultSettings) {
2845                OsiSolverInterface * solver = model.solver();
2846                if (!doScaling)
2847                  solver->setHintParam(OsiDoScale,false,OsiHintTry);
2848                OsiClpSolverInterface * si =
2849                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2850                assert (si != NULL);
2851                // get clp itself
2852                ClpSimplex * modelC = si->getModelPtr();
2853                //if (modelC->tightenPrimalBounds()!=0) {
2854                //std::cout<<"Problem is infeasible!"<<std::endl;
2855                //break;
2856                //}
2857                // bounds based on continuous
2858                if (tightenFactor&&!complicatedInteger) {
2859                  if (modelC->tightenPrimalBounds(tightenFactor)!=0) {
2860                    std::cout<<"Problem is infeasible!"<<std::endl;
2861                    break;
2862                  }
2863                }
2864                if (djFix<1.0e20) {
2865                  // do some fixing
2866                  int numberColumns = modelC->numberColumns();
2867                  int i;
2868                  const char * type = modelC->integerInformation();
2869                  double * lower = modelC->columnLower();
2870                  double * upper = modelC->columnUpper();
2871                  double * solution = modelC->primalColumnSolution();
2872                  double * dj = modelC->dualColumnSolution();
2873                  int numberFixed=0;
2874                  for (i=0;i<numberColumns;i++) {
2875                    if (type[i]) {
2876                      double value = solution[i];
2877                      if (value<lower[i]+1.0e-5&&dj[i]>djFix) {
2878                        solution[i]=lower[i];
2879                        upper[i]=lower[i];
2880                        numberFixed++;
2881                      } else if (value>upper[i]-1.0e-5&&dj[i]<-djFix) {
2882                        solution[i]=upper[i];
2883                        lower[i]=upper[i];
2884                        numberFixed++;
2885                      }
2886                    }
2887                  }
2888                  printf("%d columns fixed\n",numberFixed);
2889                }
2890              }
2891              // See if we want preprocessing
2892              OsiSolverInterface * saveSolver=NULL;
2893              CglPreProcess process;
2894              delete babModel;
2895              babModel = new CbcModel(model);
2896              OsiSolverInterface * solver3 = clpSolver->clone();
2897              babModel->assignSolver(solver3);
2898              OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2899              int numberChanged=0;
2900              if (clpSolver2->messageHandler()->logLevel())
2901                clpSolver2->messageHandler()->setLogLevel(1);
2902              if (logLevel>-1)
2903                clpSolver2->messageHandler()->setLogLevel(logLevel);
2904              lpSolver = clpSolver2->getModelPtr();
2905              if (lpSolver->factorizationFrequency()==200&&!miplib) {
2906                // User did not touch preset
2907                int numberRows = lpSolver->numberRows();
2908                const int cutoff1=10000;
2909                const int cutoff2=100000;
2910                const int base=75;
2911                const int freq0 = 50;
2912                const int freq1=200;
2913                const int freq2=400;
2914                const int maximum=1000;
2915                int frequency;
2916                if (numberRows<cutoff1)
2917                  frequency=base+numberRows/freq0;
2918                else if (numberRows<cutoff2)
2919                  frequency=base+cutoff1/freq0 + (numberRows-cutoff1)/freq1;
2920                else
2921                  frequency=base+cutoff1/freq0 + (cutoff2-cutoff1)/freq1 + (numberRows-cutoff2)/freq2;
2922                lpSolver->setFactorizationFrequency(CoinMin(maximum,frequency));
2923              }
2924              time2 = CoinCpuTime();
2925              totalTime += time2-time1;
2926              time1 = time2;
2927              double timeLeft = babModel->getMaximumSeconds();
2928              int numberOriginalColumns = babModel->solver()->getNumCols();
2929              if (preProcess==7) {
2930                // use strategy instead
2931                preProcess=0;
2932                useStrategy=true;
2933#ifdef COIN_HAS_LINK
2934                // empty out any cuts
2935                if (storedAmpl.sizeRowCuts()) {
2936                  printf("Emptying ampl stored cuts as internal preprocessing\n");
2937                  CglStored temp;
2938                  storedAmpl=temp;
2939                }
2940#endif
2941              }
2942              if (preProcess&&type==BAB) {
2943                // See if sos from mps file
2944                if (numberSOS==0&&clpSolver->numberSOS()&&doSOS) {
2945                  // SOS
2946                  numberSOS = clpSolver->numberSOS();
2947                  const CoinSet * setInfo = clpSolver->setInfo();
2948                  sosStart = new int [numberSOS+1];
2949                  sosType = new char [numberSOS];
2950                  int i;
2951                  int nTotal=0;
2952                  sosStart[0]=0;
2953                  for ( i=0;i<numberSOS;i++) {
2954                    int type = setInfo[i].setType();
2955                    int n=setInfo[i].numberEntries();
2956                    sosType[i]=type;
2957                    nTotal += n;
2958                    sosStart[i+1] = nTotal;
2959                  }
2960                  sosIndices = new int[nTotal];
2961                  sosReference = new double [nTotal];
2962                  for (i=0;i<numberSOS;i++) {
2963                    int n=setInfo[i].numberEntries();
2964                    const int * which = setInfo[i].which();
2965                    const double * weights = setInfo[i].weights();
2966                    int base = sosStart[i];
2967                    for (int j=0;j<n;j++) {
2968                      int k=which[j];
2969                      sosIndices[j+base]=k;
2970                      sosReference[j+base] = weights ? weights[j] : (double) j;
2971                    }
2972                  }
2973                }
2974                saveSolver=babModel->solver()->clone();
2975                /* Do not try and produce equality cliques and
2976                   do up to 10 passes */
2977                OsiSolverInterface * solver2;
2978                {
2979                  // Tell solver we are in Branch and Cut
2980                  saveSolver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo) ;
2981                  // Default set of cut generators
2982                  CglProbing generator1;
2983                  generator1.setUsingObjective(1);
2984                  generator1.setMaxPass(3);
2985                  generator1.setMaxProbeRoot(saveSolver->getNumCols());
2986                  generator1.setMaxElements(100);
2987                  generator1.setMaxLookRoot(50);
2988                  generator1.setRowCuts(3);
2989                  // Add in generators
2990                  process.addCutGenerator(&generator1);
2991                  int translate[]={9999,0,0,-1,2,3,-2};
2992                  process.messageHandler()->setLogLevel(babModel->logLevel());
2993#ifdef COIN_HAS_ASL
2994                  if (info.numberSos&&doSOS&&usingAmpl) {
2995                    // SOS
2996                    numberSOS = info.numberSos;
2997                    sosStart = info.sosStart;
2998                    sosIndices = info.sosIndices;
2999                  }
3000#endif
3001                  if (numberSOS&&doSOS) {
3002                    // SOS
3003                    int numberColumns = saveSolver->getNumCols();
3004                    char * prohibited = new char[numberColumns];
3005                    memset(prohibited,0,numberColumns);
3006                    int n=sosStart[numberSOS];
3007                    for (int i=0;i<n;i++) {
3008                      int iColumn = sosIndices[i];
3009                      prohibited[iColumn]=1;
3010                    }
3011                    process.passInProhibited(prohibited,numberColumns);
3012                    delete [] prohibited;
3013                  }
3014                  int numberPasses = 10;
3015                  if (tunePreProcess>=1000000) {
3016                    numberPasses = (tunePreProcess/1000000)-1;
3017                    tunePreProcess = tunePreProcess % 1000000;
3018                  } else if (tunePreProcess>=1000) {
3019                    numberPasses = (tunePreProcess/1000)-1;
3020                    tunePreProcess = tunePreProcess % 1000;
3021                  }
3022                  if (doSprint>0) {
3023                    // Sprint for primal solves
3024                    ClpSolve::SolveType method = ClpSolve::usePrimalorSprint;
3025                    ClpSolve::PresolveType presolveType = ClpSolve::presolveOff;
3026                    int numberPasses = 5;
3027                    int options[] = {0,3,0,0,0,0};
3028                    int extraInfo[] = {-1,20,-1,-1,-1,-1};
3029                    extraInfo[1]=doSprint;
3030                    int independentOptions[] = {0,0,3};
3031                    ClpSolve clpSolve(method,presolveType,numberPasses,
3032                                      options,extraInfo,independentOptions);
3033                    // say use in OsiClp
3034                    clpSolve.setSpecialOption(6,1);
3035                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (saveSolver);
3036                    osiclp->setSolveOptions(clpSolve);
3037                    osiclp->setHintParam(OsiDoDualInResolve,false);
3038                    // switch off row copy
3039                    osiclp->getModelPtr()->setSpecialOptions(osiclp->getModelPtr()->specialOptions()|256);
3040                    osiclp->getModelPtr()->setInfeasibilityCost(1.0e11);
3041                  }
3042                  solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],numberPasses,
3043                                                         tunePreProcess);
3044                  // Tell solver we are not in Branch and Cut
3045                  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
3046                  if (solver2)
3047                    solver2->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
3048                }
3049#ifdef COIN_HAS_ASL
3050                if (!solver2&&usingAmpl) {
3051                  // infeasible
3052                  info.problemStatus=1;
3053                  info.objValue = 1.0e100;
3054                  sprintf(info.buffer,"infeasible/unbounded by pre-processing");
3055                  info.primalSolution=NULL;
3056                  info.dualSolution=NULL;
3057                  break;
3058                }
3059#endif
3060                if (!noPrinting) {
3061                  if (!solver2) {
3062                    printf("Pre-processing says infeasible or unbounded\n");
3063                    break;
3064                  } else {
3065                    printf("processed model has %d rows, %d columns and %d elements\n",
3066                           solver2->getNumRows(),solver2->getNumCols(),solver2->getNumElements());
3067                  }
3068                }
3069                //solver2->resolve();
3070                if (preProcess==2) {
3071                  OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (solver2);
3072                  ClpSimplex * lpSolver = clpSolver2->getModelPtr();
3073                  lpSolver->writeMps("presolved.mps",0,1,lpSolver->optimizationDirection());
3074                  printf("Preprocessed model (minimization) on presolved.mps\n");
3075                }
3076                // we have to keep solver2 so pass clone
3077                solver2 = solver2->clone();
3078                babModel->assignSolver(solver2);
3079                babModel->initialSolve();
3080                babModel->setMaximumSeconds(timeLeft-(CoinCpuTime()-time1));
3081              }
3082              // now tighten bounds
3083              if (!miplib) {
3084                OsiClpSolverInterface * si =
3085                  dynamic_cast<OsiClpSolverInterface *>(babModel->solver()) ;
3086                assert (si != NULL);
3087                // get clp itself
3088                ClpSimplex * modelC = si->getModelPtr();
3089                if (noPrinting)
3090                  modelC->setLogLevel(0);
3091                if (!complicatedInteger&&modelC->tightenPrimalBounds()!=0) {
3092                  std::cout<<"Problem is infeasible!"<<std::endl;
3093                  break;
3094                }
3095                modelC->dual();
3096              }
3097#if 0
3098              numberDebugValues=599;
3099              debugValues = new double[numberDebugValues];
3100              CoinZeroN(debugValues,numberDebugValues);
3101              debugValues[3]=1.0;
3102              debugValues[6]=25.0;
3103              debugValues[9]=4.0;
3104              debugValues[26]=4.0;
3105              debugValues[27]=6.0;
3106              debugValues[35]=8.0;
3107              debugValues[53]=21.0;
3108              debugValues[56]=4.0;
3109#endif
3110              if (debugValues) {
3111                // for debug
3112                std::string problemName ;
3113                babModel->solver()->getStrParam(OsiProbName,problemName) ;
3114                babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
3115                twomirGen.probname_=strdup(problemName.c_str());
3116                // checking seems odd
3117                //redsplitGen.set_given_optsol(babModel->solver()->getRowCutDebuggerAlways()->optimalSolution(),
3118                //                         babModel->getNumCols());
3119              }
3120              int testOsiOptions = parameters[whichParam(TESTOSI,numberParameters,parameters)].intValue();
3121#ifdef COIN_HAS_LINK
3122              // If linked then see if expansion wanted
3123              {
3124                OsiSolverLink * solver3 = dynamic_cast<OsiSolverLink *> (babModel->solver());
3125                if (solver3) {
3126                  int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000;
3127                  if (options) {
3128                    /*
3129                      1 - force mini branch and bound
3130                      2 - set priorities high on continuous
3131                      4 - try adding OA cuts
3132                      8 - try doing quadratic linearization
3133                      16 - try expanding knapsacks
3134                    */
3135                    if ((options&16)) {
3136                      int numberColumns = saveCoinModel.numberColumns();
3137                      int numberRows = saveCoinModel.numberRows();
3138                      whichColumn = new int[numberColumns];
3139                      knapsackStart=new int[numberRows+1];
3140                      knapsackRow=new int[numberRows];
3141                      numberKnapsack=10000;
3142                      int extra1 = parameters[whichParam(EXTRA1,numberParameters,parameters)].intValue();
3143                      int extra2 = parameters[whichParam(EXTRA2,numberParameters,parameters)].intValue();
3144                      int logLevel = parameters[log].intValue();
3145                      OsiSolverInterface * solver = expandKnapsack(saveCoinModel,whichColumn,knapsackStart,
3146                                                                   knapsackRow,numberKnapsack,
3147                                                                   storedAmpl,logLevel,extra1,extra2);
3148                      if (solver) {
3149                        clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
3150                        assert (clpSolver);
3151                        lpSolver = clpSolver->getModelPtr();
3152                        babModel->assignSolver(solver);
3153                        testOsiOptions=0;
3154                        // Priorities already done
3155                        free(info.priorities);
3156                        info.priorities=NULL;
3157                      } else {
3158                        numberKnapsack=0;
3159                        delete [] whichColumn;
3160                        delete [] knapsackStart;
3161                        delete [] knapsackRow;
3162                        whichColumn=NULL;
3163                        knapsackStart=NULL;
3164                        knapsackRow=NULL;
3165                      }
3166                    }
3167                  }
3168                }
3169              }
3170#endif
3171              if (useCosts&&testOsiOptions<0) {
3172                int numberColumns = babModel->getNumCols();
3173                int * sort = new int[numberColumns];
3174                double * dsort = new double[numberColumns];
3175                int * priority = new int [numberColumns];
3176                const double * objective = babModel->getObjCoefficients();
3177                const double * lower = babModel->getColLower() ;
3178                const double * upper = babModel->getColUpper() ;
3179                const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
3180                const int * columnLength = matrix->getVectorLengths();
3181                int iColumn;
3182                int n=0;
3183                for (iColumn=0;iColumn<numberColumns;iColumn++) {
3184                  if (babModel->isInteger(iColumn)) {
3185                    sort[n]=n;
3186                    if (useCosts==1)
3187                      dsort[n++]=-fabs(objective[iColumn]);
3188                    else if (useCosts==2)
3189                      dsort[n++]=iColumn;
3190                    else if (useCosts==3)
3191                      dsort[n++]=upper[iColumn]-lower[iColumn];
3192                    else if (useCosts==4)
3193                      dsort[n++]=-(upper[iColumn]-lower[iColumn]);
3194                    else if (useCosts==5)
3195                      dsort[n++]=-columnLength[iColumn];
3196                  }
3197                }
3198                CoinSort_2(dsort,dsort+n,sort);
3199                int level=0;
3200                double last = -1.0e100;
3201                for (int i=0;i<n;i++) {
3202                  int iPut=sort[i];
3203                  if (dsort[i]!=last) {
3204                    level++;
3205                    last=dsort[i];
3206                  }
3207                  priority[iPut]=level;
3208                }
3209                babModel->passInPriorities( priority,false);
3210                delete [] priority;
3211                delete [] sort;
3212                delete [] dsort;
3213              }
3214              // FPump done first as it only works if no solution
3215              CbcHeuristicFPump heuristic4(*babModel);
3216              if (useFpump) {
3217                heuristic4.setMaximumPasses(parameters[whichParam(FPUMPITS,numberParameters,parameters)].intValue());
3218                int pumpTune=parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].intValue();
3219                if (pumpTune>0) {
3220                  /*
3221                    >=10000000 for using obj
3222                    >=1000000 use as accumulate switch
3223                    >=1000 use index+1 as number of large loops
3224                    >=100 use 0.05 objvalue as increment
3225                    >=10 use +0.1 objvalue for cutoff (add)
3226                    1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
3227                    4 and static continuous, 5 as 3 but no internal integers
3228                    6 as 3 but all slack basis!
3229                  */
3230                  int logLevel = parameters[log].intValue();
3231                  double value = babModel->solver()->getObjSense()*babModel->solver()->getObjValue();
3232                  int w = pumpTune/10;
3233                  int c = w % 10;
3234                  w /= 10;
3235                  int i = w % 10;
3236                  w /= 10;
3237                  int r = w;
3238                  int accumulate = r/1000;
3239                  r -= 1000*accumulate;
3240                  if (accumulate>=10) {
3241                    int which = accumulate/10;
3242                    accumulate -= 10*which;
3243                    which--;
3244                    // weights and factors
3245                    double weight[]={0.1,0.1,0.5,0.5,1.0,1.0,5.0,5.0};
3246                    double factor[] = {0.1,0.5,0.1,0.5,0.1,0.5,0.1,0.5};
3247                    heuristic4.setInitialWeight(weight[which]);
3248                    heuristic4.setWeightFactor(factor[which]);
3249                  }
3250                  // fake cutoff
3251                  if (logLevel>1)
3252                    printf("Setting ");
3253                  if (c) {
3254                    double cutoff;
3255                    babModel->solver()->getDblParam(OsiDualObjectiveLimit,cutoff);
3256                    cutoff = CoinMin(cutoff,value + 0.1*fabs(value)*c);
3257                    heuristic4.setFakeCutoff(cutoff);
3258                    if (logLevel>1)
3259                      printf("fake cutoff of %g ",cutoff);
3260                  }
3261                  if (i||r) {
3262                    // also set increment
3263                    heuristic4.setAbsoluteIncrement((0.01*i+0.005)*(fabs(value)+1.0e-12));
3264                    heuristic4.setAccumulate(accumulate);
3265                    heuristic4.setMaximumRetries(r+1);
3266                    if (logLevel>1) {
3267                      if (i) 
3268                        printf("increment of %g ",heuristic4.absoluteIncrement());
3269                      if (accumulate) 
3270                        printf("accumulate of %d ",accumulate);
3271                      printf("%d retries ",r+2);
3272                    }
3273                  }
3274                  pumpTune = pumpTune%100;
3275                  if (logLevel>1)
3276                    printf("and setting when to %d\n",pumpTune+10);
3277                  if (pumpTune==6)
3278                    pumpTune =13;
3279                  heuristic4.setWhen(pumpTune+10);
3280                }
3281                babModel->addHeuristic(&heuristic4);
3282              }
3283              if (!miplib) {
3284                CbcRounding heuristic1(*babModel);
3285                if (useRounding)
3286                  babModel->addHeuristic(&heuristic1) ;
3287                CbcHeuristicLocal heuristic2(*babModel);
3288                heuristic2.setSearchType(1);
3289                if (useCombine)
3290                  babModel->addHeuristic(&heuristic2);
3291                CbcHeuristicGreedyCover heuristic3(*babModel);
3292                CbcHeuristicGreedyEquality heuristic3a(*babModel);
3293                if (useGreedy) {
3294                  babModel->addHeuristic(&heuristic3);
3295                  babModel->addHeuristic(&heuristic3a);
3296                }
3297                if (useLocalTree) {
3298                  CbcTreeLocal localTree(babModel,NULL,10,0,0,10000,2000);
3299                  babModel->passInTreeHandler(localTree);
3300                }
3301              }
3302              CbcHeuristicRINS heuristic5(*babModel);
3303              if (useRINS)
3304                babModel->addHeuristic(&heuristic5) ;
3305              if (type==MIPLIB) {
3306                if (babModel->numberStrong()==5&&babModel->numberBeforeTrust()==5) 
3307                  babModel->setNumberBeforeTrust(50);
3308              }
3309              // add cut generators if wanted
3310              int switches[20];
3311              int numberGenerators=0;
3312              int translate[]={-100,-1,-99,-98,1,1,1,1};
3313              if (probingAction) {
3314                if (probingAction==5||probingAction==7)
3315                  probingGen.setRowCuts(-3); // strengthening etc just at root
3316                if (probingAction==6||probingAction==7) {
3317                  // Number of unsatisfied variables to look at
3318                  probingGen.setMaxProbe(1000);
3319                  probingGen.setMaxProbeRoot(1000);
3320                  // How far to follow the consequences
3321                  probingGen.setMaxLook(50);
3322                  probingGen.setMaxLookRoot(50);
3323                }
3324                babModel->addCutGenerator(&probingGen,translate[probingAction],"Probing");
3325                switches[numberGenerators++]=0;
3326              }
3327              if (gomoryAction&&(complicatedInteger!=1||gomoryAction==1)) {
3328                babModel->addCutGenerator(&gomoryGen,translate[gomoryAction],"Gomory");
3329                switches[numberGenerators++]=-1;
3330              }
3331              if (knapsackAction) {
3332                babModel->addCutGenerator(&knapsackGen,translate[knapsackAction],"Knapsack");
3333                switches[numberGenerators++]=0;
3334              }
3335              if (redsplitAction&&!complicatedInteger) {
3336                babModel->addCutGenerator(&redsplitGen,translate[redsplitAction],"Reduce-and-split");
3337                switches[numberGenerators++]=1;
3338              }
3339              if (cliqueAction) {
3340                babModel->addCutGenerator(&cliqueGen,translate[cliqueAction],"Clique");
3341                switches[numberGenerators++]=0;
3342              }
3343              if (mixedAction) {
3344                babModel->addCutGenerator(&mixedGen,translate[mixedAction],"MixedIntegerRounding2");
3345                switches[numberGenerators++]=-1;
3346              }
3347              if (flowAction) {
3348                babModel->addCutGenerator(&flowGen,translate[flowAction],"FlowCover");
3349                switches[numberGenerators++]=1;
3350              }
3351              if (twomirAction&&!complicatedInteger) {
3352                babModel->addCutGenerator(&twomirGen,translate[twomirAction],"TwoMirCuts");
3353                switches[numberGenerators++]=1;
3354              }
3355              if (landpAction) {
3356                babModel->addCutGenerator(&landpGen,translate[landpAction],"LiftAndProject");
3357                switches[numberGenerators++]=1;
3358              }
3359              if (storedCuts) 
3360                babModel->setSpecialOptions(babModel->specialOptions()|64);
3361              // Say we want timings
3362              numberGenerators = babModel->numberCutGenerators();
3363              int iGenerator;
3364              int cutDepth=
3365                parameters[whichParam(CUTDEPTH,numberParameters,parameters)].intValue();
3366              for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
3367                CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
3368                int howOften = generator->howOften();
3369                if (howOften==-98||howOften==-99) 
3370                  generator->setSwitchOffIfLessThan(switches[iGenerator]);
3371                generator->setTiming(true);
3372                if (cutDepth>=0)
3373                  generator->setWhatDepth(cutDepth) ;
3374              }
3375              // Could tune more
3376              if (!miplib) {
3377                babModel->setMinimumDrop(min(5.0e-2,
3378                                             fabs(babModel->getMinimizationObjValue())*1.0e-3+1.0e-4));
3379                if (cutPass==-1234567) {
3380                  if (babModel->getNumCols()<500)
3381                    babModel->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
3382                  else if (babModel->getNumCols()<5000)
3383                    babModel->setMaximumCutPassesAtRoot(100); // use minimum drop
3384                  else
3385                    babModel->setMaximumCutPassesAtRoot(20);
3386                } else {
3387                  babModel->setMaximumCutPassesAtRoot(cutPass);
3388                }
3389                babModel->setMaximumCutPasses(1);
3390              }
3391              // Do more strong branching if small
3392              //if (babModel->getNumCols()<5000)
3393              //babModel->setNumberStrong(20);
3394              // Switch off strong branching if wanted
3395              //if (babModel->getNumCols()>10*babModel->getNumRows())
3396              //babModel->setNumberStrong(0);
3397              if (!noPrinting) {
3398                int iLevel = parameters[log].intValue();
3399                if (iLevel<0) {
3400                  babModel->setPrintingMode(1);
3401                  iLevel = -iLevel;
3402                }
3403                babModel->messageHandler()->setLogLevel(iLevel);
3404                if (babModel->getNumCols()>2000||babModel->getNumRows()>1500||
3405                    babModel->messageHandler()->logLevel()>1)
3406                  babModel->setPrintFrequency(100);
3407              }
3408             
3409              babModel->solver()->setIntParam(OsiMaxNumIterationHotStart,
3410                    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].intValue());
3411              OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3412              // go faster stripes
3413              if (osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
3414                osiclp->setupForRepeatedUse(2,parameters[slog].intValue());
3415              } else {
3416                osiclp->setupForRepeatedUse(0,parameters[slog].intValue());
3417              }
3418              double increment=babModel->getCutoffIncrement();;
3419              int * changed = NULL;
3420              if (!miplib&&increment==normalIncrement)
3421                changed=analyze( osiclp,numberChanged,increment,false);
3422              if (debugValues) {
3423                if (numberDebugValues==babModel->getNumCols()) {
3424                  // for debug
3425                  babModel->solver()->activateRowCutDebugger(debugValues) ;
3426                } else {
3427                  printf("debug file has incorrect number of columns\n");
3428                }
3429              }
3430              babModel->setCutoffIncrement(CoinMax(babModel->getCutoffIncrement(),increment));
3431              // Turn this off if you get problems
3432              // Used to be automatically set
3433              int mipOptions = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()%10000;
3434              if (mipOptions!=(1))
3435                printf("mip options %d\n",mipOptions);
3436              osiclp->setSpecialOptions(mipOptions);
3437              if (gapRatio < 1.0e100) {
3438                double value = babModel->solver()->getObjValue() ;
3439                double value2 = gapRatio*(1.0e-5+fabs(value)) ;
3440                babModel->setAllowableGap(value2) ;
3441                std::cout << "Continuous " << value
3442                          << ", so allowable gap set to "
3443                          << value2 << std::endl ;
3444              }
3445              // probably faster to use a basis to get integer solutions
3446              babModel->setSpecialOptions(babModel->specialOptions()|2);
3447              currentBranchModel = babModel;
3448              OsiSolverInterface * strengthenedModel=NULL;
3449              if (type==BAB||type==MIPLIB) {
3450                int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
3451                if (moreMipOptions>=0) {
3452                  printf("more mip options %d\n",moreMipOptions);
3453                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3454                  if (moreMipOptions==10000) {
3455                    // test memory saving
3456                    moreMipOptions -= 10000;
3457                    ClpSimplex * lpSolver = osiclp->getModelPtr();
3458                    lpSolver->setPersistenceFlag(1);
3459                    // switch off row copy if few rows
3460                    if (lpSolver->numberRows()<150)
3461                      lpSolver->setSpecialOptions(lpSolver->specialOptions()|256);
3462                  }
3463                  if (((moreMipOptions+1)%1000000)!=0)
3464                    babModel->setSearchStrategy(moreMipOptions%1000000);
3465                  // go faster stripes
3466                  if( moreMipOptions >=999999) {
3467                    if (osiclp) {
3468                      int save = osiclp->specialOptions();
3469                      osiclp->setupForRepeatedUse(2,0);
3470                      osiclp->setSpecialOptions(save|osiclp->specialOptions());
3471                    }
3472                  } 
3473                }
3474              }
3475              if (type==BAB) {
3476#ifdef COIN_HAS_ASL
3477                if (usingAmpl) {
3478                  priorities=info.priorities;
3479                  branchDirection=info.branchDirection;
3480                  pseudoDown=info.pseudoDown;
3481                  pseudoUp=info.pseudoUp;
3482                  solutionIn=info.primalSolution;
3483                  prioritiesIn = info.priorities;
3484                  if (info.numberSos&&doSOS) {
3485                    // SOS
3486                    numberSOS = info.numberSos;
3487                    sosStart = info.sosStart;
3488                    sosIndices = info.sosIndices;
3489                    sosType = info.sosType;
3490                    sosReference = info.sosReference;
3491                    sosPriority = info.sosPriority;
3492                  }
3493                }
3494#endif               
3495                const int * originalColumns = preProcess ? process.originalColumns() : NULL;
3496                if (solutionIn&&useSolution) {
3497                  if (preProcess) {
3498                    int numberColumns = babModel->getNumCols();
3499                    // extend arrays in case SOS
3500                    int n = originalColumns[numberColumns-1]+1;
3501                    int nSmaller = CoinMin(n,numberOriginalColumns);
3502                    double * solutionIn2 = new double [n];
3503                    int * prioritiesIn2 = new int[n];
3504                    int i;
3505                    for (i=0;i<nSmaller;i++) {
3506                      solutionIn2[i]=solutionIn[i];
3507                      prioritiesIn2[i]=prioritiesIn[i];
3508                    }
3509                    for (;i<n;i++) {
3510                      solutionIn2[i]=0.0;
3511                      prioritiesIn2[i]=1000000;
3512                    }
3513                    int iLast=-1;
3514                    for (i=0;i<numberColumns;i++) {
3515                      int iColumn = originalColumns[i];
3516                      assert (iColumn>iLast);
3517                      iLast=iColumn;
3518                      solutionIn2[i]=solutionIn2[iColumn];
3519                      if (prioritiesIn)
3520                        prioritiesIn2[i]=prioritiesIn2[iColumn];
3521                    }
3522                    babModel->setHotstartSolution(solutionIn2,prioritiesIn2);
3523                    delete [] solutionIn2;
3524                    delete [] prioritiesIn2;
3525                  } else {
3526                    babModel->setHotstartSolution(solutionIn,prioritiesIn);
3527                  }
3528                }
3529                OsiSolverInterface * testOsiSolver= (testOsiOptions>=0) ? babModel->solver() : NULL;
3530                if (!testOsiSolver) {
3531                  // *************************************************************
3532                  // CbcObjects
3533                  if (preProcess&&process.numberSOS()) {
3534                    int numberSOS = process.numberSOS();
3535                    int numberIntegers = babModel->numberIntegers();
3536                    /* model may not have created objects
3537                       If none then create
3538                    */
3539                    if (!numberIntegers||!babModel->numberObjects()) {
3540                      int type = (pseudoUp) ? 1 : 0;
3541                      babModel->findIntegers(true,type);
3542                      numberIntegers = babModel->numberIntegers();
3543                    }
3544                    OsiObject ** oldObjects = babModel->objects();
3545                    // Do sets and priorities
3546                    OsiObject ** objects = new OsiObject * [numberSOS];
3547                    // set old objects to have low priority
3548                    int numberOldObjects = babModel->numberObjects();
3549                    int numberColumns = babModel->getNumCols();
3550                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
3551                      oldObjects[iObj]->setPriority(numberColumns+1);
3552                      int iColumn = oldObjects[iObj]->columnNumber();
3553                      assert (iColumn>=0);
3554                      if (iColumn>=numberOriginalColumns)
3555                        continue;
3556                      if (originalColumns)
3557                        iColumn = originalColumns[iColumn];
3558                      if (branchDirection) {
3559                        CbcSimpleInteger * obj =
3560                          dynamic_cast <CbcSimpleInteger *>(oldObjects[iObj]) ;
3561                        if (obj) { 
3562                          obj->setPreferredWay(branchDirection[iColumn]);
3563                        } else {
3564                          CbcObject * obj =
3565                            dynamic_cast <CbcObject *>(oldObjects[iObj]) ;
3566                          assert (obj);
3567                          obj->setPreferredWay(branchDirection[iColumn]);
3568                        }
3569                      }
3570                      if (pseudoUp) {
3571                        CbcSimpleIntegerPseudoCost * obj1a =
3572                          dynamic_cast <CbcSimpleIntegerPseudoCost *>(oldObjects[iObj]) ;
3573                        assert (obj1a);
3574                        if (pseudoDown[iColumn]>0.0)
3575                          obj1a->setDownPseudoCost(pseudoDown[iColumn]);
3576                        if (pseudoUp[iColumn]>0.0)
3577                          obj1a->setUpPseudoCost(pseudoUp[iColumn]);
3578                      }
3579                    }
3580                    const int * starts = process.startSOS();
3581                    const int * which = process.whichSOS();
3582                    const int * type = process.typeSOS();
3583                    const double * weight = process.weightSOS();
3584                    int iSOS;
3585                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
3586                      int iStart = starts[iSOS];
3587                      int n=starts[iSOS+1]-iStart;
3588                      objects[iSOS] = new CbcSOS(babModel,n,which+iStart,weight+iStart,
3589                                                 iSOS,type[iSOS]);
3590                      // branch on long sets first
3591                      objects[iSOS]->setPriority(numberColumns-n);
3592                    }
3593                    babModel->addObjects(numberSOS,objects);
3594                    for (iSOS=0;iSOS<numberSOS;iSOS++)
3595                      delete objects[iSOS];
3596                    delete [] objects;
3597                  } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
3598                    // do anyway for priorities etc
3599                    int numberIntegers = babModel->numberIntegers();
3600                    /* model may not have created objects
3601                       If none then create
3602                    */
3603                    if (!numberIntegers||!babModel->numberObjects()) {
3604                      int type = (pseudoUp) ? 1 : 0;
3605                      babModel->findIntegers(true,type);
3606                    }
3607                    if (numberSOS) {
3608                      // Do sets and priorities
3609                      OsiObject ** objects = new OsiObject * [numberSOS];
3610                      int iSOS;
3611                      if (originalColumns) {
3612                        // redo sequence numbers
3613                        int numberColumns = babModel->getNumCols();
3614                        int nOld = originalColumns[numberColumns-1]+1;
3615                        int * back = new int[nOld];
3616                        int i;
3617                        for (i=0;i<nOld;i++)
3618                          back[i]=-1;
3619                        for (i=0;i<numberColumns;i++)
3620                          back[originalColumns[i]]=i;
3621                        // Really need better checks
3622                        int nMissing=0;
3623                        int n=sosStart[numberSOS];
3624                        for (i=0;i<n;i++) {
3625                          int iColumn = sosIndices[i];
3626                          int jColumn = back[iColumn];
3627                          if (jColumn>=0) 
3628                            sosIndices[i] = jColumn;
3629                          else 
3630                            nMissing++;
3631                        }
3632                        delete [] back;
3633                        if (nMissing)
3634                          printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
3635                      }
3636                      for (iSOS =0;iSOS<numberSOS;iSOS++) {
3637                        int iStart = sosStart[iSOS];
3638                        int n=sosStart[iSOS+1]-iStart;
3639                        objects[iSOS] = new CbcSOS(babModel,n,sosIndices+iStart,sosReference+iStart,
3640                                                   iSOS,sosType[iSOS]);
3641                        if (sosPriority)
3642                          objects[iSOS]->setPriority(sosPriority[iSOS]);
3643                        else if (!prioritiesIn)
3644                          objects[iSOS]->setPriority(10);  // rather than 1000
3645                      }
3646                      // delete any existing SOS objects
3647                      int numberObjects=babModel->numberObjects();
3648                      OsiObject ** oldObjects=babModel->objects();
3649                      int nNew=0;
3650                      for (int i=0;i<numberObjects;i++) {
3651                        OsiObject * objThis = oldObjects[i];
3652                        CbcSOS * obj1 =
3653                          dynamic_cast <CbcSOS *>(objThis) ;
3654                        OsiSOS * obj2 =
3655                          dynamic_cast <OsiSOS *>(objThis) ;
3656                        if (!obj1&&!obj2) {
3657                          oldObjects[nNew++]=objThis;
3658                        } else {
3659                          delete objThis;
3660                        }
3661                      }
3662                      babModel->setNumberObjects(nNew);
3663                      babModel->addObjects(numberSOS,objects);
3664                      for (iSOS=0;iSOS<numberSOS;iSOS++)
3665                        delete objects[iSOS];
3666                      delete [] objects;
3667                    }
3668                  }
3669                  OsiObject ** objects = babModel->objects();
3670                  int numberObjects = babModel->numberObjects();
3671                  for (int iObj = 0;iObj<numberObjects;iObj++) {
3672                    // skip sos
3673                    CbcSOS * objSOS =
3674                      dynamic_cast <CbcSOS *>(objects[iObj]) ;
3675                    if (objSOS)
3676                      continue;
3677                    int iColumn = objects[iObj]->columnNumber();
3678                    assert (iColumn>=0);
3679                    if (originalColumns)
3680                      iColumn = originalColumns[iColumn];
3681                    if (branchDirection) {
3682                      CbcSimpleInteger * obj =
3683                        dynamic_cast <CbcSimpleInteger *>(objects[iObj]) ;
3684                      if (obj) { 
3685                        obj->setPreferredWay(branchDirection[iColumn]);
3686                      } else {
3687                        CbcObject * obj =
3688                          dynamic_cast <CbcObject *>(objects[iObj]) ;
3689                        assert (obj);
3690                        obj->setPreferredWay(branchDirection[iColumn]);
3691                      }
3692                    }
3693                    if (priorities) {
3694                      int iPriority = priorities[iColumn];
3695                      if (iPriority>0)
3696                        objects[iObj]->setPriority(iPriority);
3697                    }
3698                    if (pseudoUp&&pseudoUp[iColumn]) {
3699                      CbcSimpleIntegerPseudoCost * obj1a =
3700                        dynamic_cast <CbcSimpleIntegerPseudoCost *>(objects[iObj]) ;
3701                      assert (obj1a);
3702                      if (pseudoDown[iColumn]>0.0)
3703                        obj1a->setDownPseudoCost(pseudoDown[iColumn]);
3704                      if (pseudoUp[iColumn]>0.0)
3705                        obj1a->setUpPseudoCost(pseudoUp[iColumn]);
3706                    }
3707                  }
3708                  // *************************************************************
3709                } else {
3710                  // *************************************************************
3711                  // OsiObjects
3712                  // Find if none
3713                  int numberIntegers = testOsiSolver->getNumIntegers();
3714                  /* model may not have created objects
3715                     If none then create
3716                  */
3717                  if (!numberIntegers||!testOsiSolver->numberObjects()) {
3718                    //int type = (pseudoUp) ? 1 : 0;
3719                    testOsiSolver->findIntegers(false);
3720                    numberIntegers = testOsiSolver->getNumIntegers();
3721                  }
3722                  if (preProcess&&process.numberSOS()) {
3723                    int numberSOS = process.numberSOS();
3724                    OsiObject ** oldObjects = testOsiSolver->objects();
3725                    // Do sets and priorities
3726                    OsiObject ** objects = new OsiObject * [numberSOS];
3727                    // set old objects to have low priority
3728                    int numberOldObjects = testOsiSolver->numberObjects();
3729                    int numberColumns = testOsiSolver->getNumCols();
3730                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
3731                      oldObjects[iObj]->setPriority(numberColumns+1);
3732                      int iColumn = oldObjects[iObj]->columnNumber();
3733                      assert (iColumn>=0);
3734                      if (iColumn>=numberOriginalColumns)
3735                        continue;
3736                      if (originalColumns)
3737                        iColumn = originalColumns[iColumn];
3738                      if (branchDirection) {
3739                        OsiSimpleInteger * obj =
3740                          dynamic_cast <OsiSimpleInteger *>(oldObjects[iObj]) ;
3741                        if (obj) { 
3742                          obj->setPreferredWay(branchDirection[iColumn]);
3743                        } else {
3744                          OsiObject2 * obj =
3745                            dynamic_cast <OsiObject2 *>(oldObjects[iObj]) ;
3746                          if (obj)
3747                            obj->setPreferredWay(branchDirection[iColumn]);
3748                        }
3749                      }
3750                      if (pseudoUp) {
3751                        abort();
3752                      }
3753                    }
3754                    const int * starts = process.startSOS();
3755                    const int * which = process.whichSOS();
3756                    const int * type = process.typeSOS();
3757                    const double * weight = process.weightSOS();
3758                    int iSOS;
3759                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
3760                      int iStart = starts[iSOS];
3761                      int n=starts[iSOS+1]-iStart;
3762                      objects[iSOS] = new OsiSOS(testOsiSolver,n,which+iStart,weight+iStart,
3763                                                 type[iSOS]);
3764                      // branch on long sets first
3765                      objects[iSOS]->setPriority(numberColumns-n);
3766                    }
3767                    testOsiSolver->addObjects(numberSOS,objects);
3768                    for (iSOS=0;iSOS<numberSOS;iSOS++)
3769                      delete objects[iSOS];
3770                    delete [] objects;
3771                  } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
3772                    if (numberSOS) {
3773                      // Do sets and priorities
3774                      OsiObject ** objects = new OsiObject * [numberSOS];
3775                      int iSOS;
3776                      if (originalColumns) {
3777                        // redo sequence numbers
3778                        int numberColumns = testOsiSolver->getNumCols();
3779                        int nOld = originalColumns[numberColumns-1]+1;
3780                        int * back = new int[nOld];
3781                        int i;
3782                        for (i=0;i<nOld;i++)
3783                          back[i]=-1;
3784                        for (i=0;i<numberColumns;i++)
3785                          back[originalColumns[i]]=i;
3786                        // Really need better checks
3787                        int nMissing=0;
3788                        int n=sosStart[numberSOS];
3789                        for (i=0;i<n;i++) {
3790                          int iColumn = sosIndices[i];
3791                          int jColumn = back[iColumn];
3792                          if (jColumn>=0) 
3793                            sosIndices[i] = jColumn;
3794                          else 
3795                            nMissing++;
3796                        }
3797                        delete [] back;
3798                        if (nMissing)
3799                          printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
3800                      }
3801                      for (iSOS =0;iSOS<numberSOS;iSOS++) {
3802                        int iStart = sosStart[iSOS];
3803                        int n=sosStart[iSOS+1]-iStart;
3804                        objects[iSOS] = new OsiSOS(testOsiSolver,n,sosIndices+iStart,sosReference+iStart,
3805                                                   sosType[iSOS]);
3806                        if (sosPriority)
3807                          objects[iSOS]->setPriority(sosPriority[iSOS]);
3808                        else if (!prioritiesIn)
3809                          objects[iSOS]->setPriority(10);  // rather than 1000
3810                      }
3811                      // delete any existing SOS objects
3812                      int numberObjects=testOsiSolver->numberObjects();
3813                      OsiObject ** oldObjects=testOsiSolver->objects();
3814                      int nNew=0;
3815                      for (int i=0;i<numberObjects;i++) {
3816                        OsiObject * objThis = oldObjects[i];
3817                        OsiSOS * obj1 =
3818                          dynamic_cast <OsiSOS *>(objThis) ;
3819                        OsiSOS * obj2 =
3820                          dynamic_cast <OsiSOS *>(objThis) ;
3821                        if (!obj1&&!obj2) {
3822                          oldObjects[nNew++]=objThis;
3823                        } else {
3824                          delete objThis;
3825                        }
3826                      }
3827                      testOsiSolver->setNumberObjects(nNew);
3828                      testOsiSolver->addObjects(numberSOS,objects);
3829                      for (iSOS=0;iSOS<numberSOS;iSOS++)
3830                        delete objects[iSOS];
3831                      delete [] objects;
3832                    }
3833                  }
3834                  OsiObject ** objects = testOsiSolver->objects();
3835                  int numberObjects = testOsiSolver->numberObjects();
3836                  int logLevel = parameters[log].intValue();
3837                  for (int iObj = 0;iObj<numberObjects;iObj++) {
3838                    // skip sos
3839                    OsiSOS * objSOS =
3840                      dynamic_cast <OsiSOS *>(objects[iObj]) ;
3841                    if (objSOS) {
3842                      if (logLevel>2)
3843                        printf("Set %d is SOS - priority %d\n",iObj,objSOS->priority());
3844                      continue;
3845                    }
3846                    int iColumn = objects[iObj]->columnNumber();
3847                    if (iColumn>=0) {
3848                      if (originalColumns)
3849                        iColumn = originalColumns[iColumn];
3850                      if (branchDirection) {
3851                        OsiSimpleInteger * obj =
3852                          dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
3853                        if (obj) { 
3854                          obj->setPreferredWay(branchDirection[iColumn]);
3855                        } else {
3856                          OsiObject2 * obj =
3857                            dynamic_cast <OsiObject2 *>(objects[iObj]) ;
3858                          if (obj)
3859                            obj->setPreferredWay(branchDirection[iColumn]);
3860                        }
3861                      }
3862                    }
3863                    if (priorities) {
3864                      int iPriority = priorities[iColumn];
3865                      if (iPriority>0)
3866                        objects[iObj]->setPriority(iPriority);
3867                    }
3868                    if (logLevel>2)
3869                      printf("Obj %d is int? - priority %d\n",iObj,objects[iObj]->priority());
3870                    if (pseudoUp&&pseudoUp[iColumn]) {
3871                      abort();
3872                    }
3873                  }
3874                  // *************************************************************
3875                }
3876                int statistics = (printOptions>0) ? printOptions: 0;
3877#ifdef COIN_HAS_ASL
3878                if (!usingAmpl) {
3879#endif
3880                  free(priorities);
3881                  priorities=NULL;
3882                  free(branchDirection);
3883                  branchDirection=NULL;
3884                  free(pseudoDown);
3885                  pseudoDown=NULL;
3886                  free(pseudoUp);
3887                  pseudoUp=NULL;
3888                  free(solutionIn);
3889                  solutionIn=NULL;
3890                  free(prioritiesIn);
3891                  prioritiesIn=NULL;
3892                  free(sosStart);
3893                  sosStart=NULL;
3894                  free(sosIndices);
3895                  sosIndices=NULL;
3896                  free(sosType);
3897                  sosType=NULL;
3898                  free(sosReference);
3899                  sosReference=NULL;
3900                  free(cut);
3901                  cut=NULL;
3902                  free(sosPriority);
3903                  sosPriority=NULL;
3904#ifdef COIN_HAS_ASL
3905                }
3906#endif               
3907                if (nodeStrategy) {
3908                  // change default
3909                  if (nodeStrategy>1) {
3910                    // up or down
3911                    int way = ((nodeStrategy%1)==1) ? -1 : +1;
3912                    babModel->setPreferredWay(way);
3913#if 0
3914                    OsiObject ** objects = babModel->objects();
3915                    int numberObjects = babModel->numberObjects();
3916                    for (int iObj = 0;iObj<numberObjects;iObj++) {
3917                      CbcObject * obj =
3918                        dynamic_cast <CbcObject *>(objects[iObj]) ;
3919                      assert (obj);
3920                      obj->setPreferredWay(way);
3921                    }
3922#endif
3923                  }
3924                  if (nodeStrategy==1||nodeStrategy>3) {
3925                    // depth
3926                    CbcCompareDefault compare;
3927                    compare.setWeight(-3.0);
3928                    babModel->setNodeComparison(compare);
3929                  }
3930                }
3931                if (cppValue>=0) {
3932                  int prepro = useStrategy ? -1 : preProcess;
3933                  // generate code
3934                  FILE * fp = fopen("user_driver.cpp","w");
3935                  if (fp) {
3936                    // generate enough to do BAB
3937                    babModel->generateCpp(fp,1);
3938                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3939                    // Make general so do factorization
3940                    int factor = osiclp->getModelPtr()->factorizationFrequency();
3941                    osiclp->getModelPtr()->setFactorizationFrequency(200);
3942                    osiclp->generateCpp(fp);
3943                    osiclp->getModelPtr()->setFactorizationFrequency(factor);
3944                    //solveOptions.generateCpp(fp);
3945                    fclose(fp);
3946                    // now call generate code
3947                    generateCode(babModel,"user_driver.cpp",cppValue,prepro);
3948                  } else {
3949                    std::cout<<"Unable to open file user_driver.cpp"<<std::endl;
3950                  }
3951                }
3952                if (!babModel->numberStrong())
3953                  babModel->setNumberBeforeTrust(0);
3954                if (useStrategy) {
3955                  CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
3956                  strategy.setupPreProcessing(1);
3957                  babModel->setStrategy(strategy);
3958                }
3959                if (testOsiOptions>=0) {
3960                  printf("Testing OsiObject options %d\n",testOsiOptions);
3961                  if (!numberSOS) {
3962                    babModel->solver()->findIntegersAndSOS(false);
3963#ifdef COIN_HAS_LINK
3964                    // If linked then pass in model
3965                    OsiSolverLink * solver3 = dynamic_cast<OsiSolverLink *> (babModel->solver());
3966                    if (solver3) {
3967                      int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000;
3968                      CglStored stored;
3969                      if (options) {
3970                        printf("nlp options %d\n",options);
3971                        /*
3972                          1 - force mini branch and bound
3973                          2 - set priorities high on continuous
3974                          4 - try adding OA cuts
3975                          8 - try doing quadratic linearization
3976                          16 - try expanding knapsacks
3977                          32 - OA cuts strictly concave
3978                        */
3979                        if ((options&2)) {
3980                          solver3->setBiLinearPriorities(10,tightenFactor > 0.0 ? tightenFactor : 1.0);
3981                        } else if (tightenFactor>0.0) {
3982                          // set grid size for all continuous bi-linear
3983                          solver3->setMeshSizes(tightenFactor);
3984                        }
3985                        if ((options&4)) {
3986                          solver3->setSpecialOptions2(solver3->specialOptions2()|(8+4));
3987                          // say convex
3988                          solver3->sayConvex((options&32)==0);
3989                        }
3990                        if ((options&1)!=0&&slpValue>0)
3991                          solver3->setFixedPriority(slpValue);
3992                        double cutoff=COIN_DBL_MAX;
3993                        if ((options&8))
3994                          cutoff=solver3->linearizedBAB(&stored);
3995                        if (cutoff<babModel->getCutoff())
3996                          babModel->setCutoff(cutoff);
3997                      }
3998                      solver3->setCbcModel(babModel);
3999                      babModel->addCutGenerator(&stored,1,"Stored");
4000                      CglTemporary temp;
4001                      babModel->addCutGenerator(&temp,1,"OnceOnly");
4002                      //choose.setNumberBeforeTrusted(2000);
4003                      //choose.setNumberStrong(20);
4004                    }
4005#endif
4006                  } else {
4007                    // move across
4008                    babModel->deleteObjects(false);
4009                    //babModel->addObjects(babModel->solver()->numberObjects(),babModel->solver()->objects());
4010                  }
4011                  CbcBranchDefaultDecision decision;
4012                  if (babModel->numberStrong()) {
4013                    OsiChooseStrong choose(babModel->solver());
4014                    choose.setNumberBeforeTrusted(babModel->numberBeforeTrust());
4015                    choose.setNumberStrong(babModel->numberStrong());
4016                    choose.setShadowPriceMode(testOsiOptions);
4017                    decision.setChooseMethod(choose);
4018                  } else {
4019                    OsiChooseVariable choose(babModel->solver());
4020                    decision.setChooseMethod(choose);
4021                  }
4022                  babModel->setBranchingMethod(decision);
4023                  if (useCosts&&testOsiOptions>=0) {
4024                    int numberColumns = babModel->getNumCols();
4025                    int * sort = new int[numberColumns];
4026                    double * dsort = new double[numberColumns];
4027                    int * priority = new int [numberColumns];
4028                    const double * objective = babModel->getObjCoefficients();
4029                    const double * lower = babModel->getColLower() ;
4030                    const double * upper = babModel->getColUpper() ;
4031                    const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
4032                    const int * columnLength = matrix->getVectorLengths();
4033                    int iColumn;
4034                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4035                      sort[iColumn]=iColumn;
4036                      if (useCosts==1)
4037                        dsort[iColumn]=-fabs(objective[iColumn]);
4038                      else if (useCosts==2)
4039                        dsort[iColumn]=iColumn;
4040                      else if (useCosts==3)
4041                        dsort[iColumn]=upper[iColumn]-lower[iColumn];
4042                      else if (useCosts==4)
4043                        dsort[iColumn]=-(upper[iColumn]-lower[iColumn]);
4044                      else if (useCosts==5)
4045                        dsort[iColumn]=-columnLength[iColumn];
4046                    }
4047                    CoinSort_2(dsort,dsort+numberColumns,sort);
4048                    int level=0;
4049                    double last = -1.0e100;
4050                    for (int i=0;i<numberColumns;i++) {
4051                      int iPut=sort[i];
4052                      if (dsort[i]!=last) {
4053                        level++;
4054                        last=dsort[i];
4055                      }
4056                      priority[iPut]=level;
4057                    }
4058                    OsiObject ** objects = babModel->objects();
4059                    int numberObjects = babModel->numberObjects();
4060                    for (int iObj = 0;iObj<numberObjects;iObj++) {
4061                      OsiObject * obj = objects[iObj] ;
4062                      int iColumn = obj->columnNumber();
4063                      if (iColumn>=0)
4064                        obj->setPriority(priority[iColumn]);
4065                    }
4066                    delete [] priority;
4067                    delete [] sort;
4068                    delete [] dsort;
4069                  }
4070                }
4071                checkSOS(babModel, babModel->solver());
4072                if (doSprint>0) {
4073                  // Sprint for primal solves
4074                  ClpSolve::SolveType method = ClpSolve::usePrimalorSprint;
4075                  ClpSolve::PresolveType presolveType = ClpSolve::presolveOff;
4076                  int numberPasses = 5;
4077                  int options[] = {0,3,0,0,0,0};
4078                  int extraInfo[] = {-1,20,-1,-1,-1,-1};
4079                  extraInfo[1]=doSprint;
4080                  int independentOptions[] = {0,0,3};
4081                  ClpSolve clpSolve(method,presolveType,numberPasses,
4082                                    options,extraInfo,independentOptions);
4083                  // say use in OsiClp
4084                  clpSolve.setSpecialOption(6,1);
4085                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4086                  osiclp->setSolveOptions(clpSolve);
4087                  osiclp->setHintParam(OsiDoDualInResolve,false);
4088                  // switch off row copy
4089                  osiclp->getModelPtr()->setSpecialOptions(osiclp->getModelPtr()->specialOptions()|256);
4090                  osiclp->getModelPtr()->setInfeasibilityCost(1.0e11);
4091                }
4092#ifdef COIN_HAS_LINK
4093                if (storedAmpl.sizeRowCuts()) {
4094                  if (preProcess) {
4095                    const int * originalColumns = process.originalColumns();
4096                    int numberColumns = babModel->getNumCols();
4097                    int * newColumn = new int[numberOriginalColumns];
4098                    int i;
4099                    for (i=0;i<numberOriginalColumns;i++) 
4100                      newColumn[i]=-1;
4101                    for (i=0;i<numberColumns;i++) {
4102                      int iColumn = originalColumns[i];
4103                      newColumn[iColumn]=i;
4104                    }
4105                    int * buildColumn = new int[numberColumns];
4106                    // Build up valid cuts
4107                    int nBad=0;
4108                    int nCuts = storedAmpl.sizeRowCuts();
4109                    CglStored newCuts;
4110                    for (i=0;i<nCuts;i++) {
4111                      const OsiRowCut * cut = storedAmpl.rowCutPointer(i);
4112                      double lb = cut->lb();
4113                      double ub = cut->ub();
4114                      int n=cut->row().getNumElements();
4115                      const int * column = cut->row().getIndices();
4116                      const double * element = cut->row().getElements();
4117                      bool bad=false;
4118                      for (int i=0;i<n;i++) {
4119                        int iColumn = column[i];
4120                        iColumn = newColumn[iColumn];
4121                        if (iColumn>=0) {
4122                          buildColumn[i]=iColumn;
4123                        } else {
4124                          bad=true;
4125                          break;
4126                        }
4127                      }
4128                      if (!bad) {
4129                        newCuts.addCut(lb,ub,n,buildColumn,element);
4130                      } else {
4131                        nBad++;
4132                      }
4133                    }
4134                    storedAmpl=newCuts;
4135                    if (nBad)
4136                      printf("%d cuts dropped\n",nBad);
4137                    delete [] newColumn;
4138                    delete [] buildColumn;
4139                  }
4140                }
4141#endif
4142#ifdef CLP_MALLOC_STATISTICS
4143                malloc_stats();
4144                malloc_stats2();
4145#endif
4146                if (outputFormat==5) {
4147                  osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4148                  lpSolver = osiclp->getModelPtr();
4149                  lpSolver->setPersistenceFlag(1);
4150                }
4151#ifdef COIN_HAS_ASL
4152                // add in lotsizing
4153                if (usingAmpl&&info.special) {
4154                  int numberColumns = babModel->getNumCols();
4155                  int i;
4156                  int n=0;
4157                  if (preProcess) {
4158                    const int * originalColumns = process.originalColumns();
4159                    for (i=0;i<numberColumns;i++) {
4160                      int iColumn = originalColumns[i];
4161                      assert (iColumn>=i);
4162                      int iType = info.special[iColumn];
4163                      if (iType) {
4164                        assert (iType==1);
4165                        n++;
4166                      }
4167                      info.special[i]=iType;
4168                    }
4169                  }
4170                  if (n) {
4171                    int numberIntegers=0;
4172                    int numberOldObjects=0;
4173                    OsiObject ** oldObjects=NULL;
4174                    const double * lower = babModel->solver()->getColLower();
4175                    const double * upper = babModel->solver()->getColUpper();
4176                    if (testOsiOptions<0) {
4177                      // *************************************************************
4178                      // CbcObjects
4179                      numberIntegers = babModel->numberIntegers();
4180                      /* model may not have created objects
4181                         If none then create
4182                      */
4183                      if (!numberIntegers||!babModel->numberObjects()) {
4184                        int type = (pseudoUp) ? 1 : 0;
4185                        babModel->findIntegers(true,type);
4186                        numberIntegers = babModel->numberIntegers();
4187                      }
4188                      oldObjects = babModel->objects();
4189                      numberOldObjects = babModel->numberObjects();
4190                    } else {
4191                      numberIntegers = testOsiSolver->getNumIntegers();
4192                      if (!numberIntegers||!testOsiSolver->numberObjects()) {
4193                        /* model may not have created objects
4194                           If none then create
4195                        */
4196                        testOsiSolver->findIntegers(false);
4197                        numberIntegers = testOsiSolver->getNumIntegers();
4198                      }
4199                      oldObjects = testOsiSolver->objects();
4200                      numberOldObjects = testOsiSolver->numberObjects();
4201                    }
4202                    OsiObject ** objects = new OsiObject * [n];
4203                    n=0;
4204                    // set new objects to have one lower priority
4205                    double ranges[] = {-COIN_DBL_MAX,-1.0,1.0,COIN_DBL_MAX};
4206                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
4207                      int iColumn = oldObjects[iObj]->columnNumber();
4208                      if (iColumn>=0&&info.special[iColumn]) {
4209                        if (lower[iColumn]<=-1.0&&upper[iColumn]>=0.0) {
4210                          ranges[0]=lower[iColumn];
4211                          ranges[3]=upper[iColumn];
4212                          int priority = oldObjects[iObj]->priority();
4213                          if (testOsiOptions<0) {
4214                            objects[n] = new CbcLotsize(babModel,iColumn,2,ranges,true);
4215                          } else {
4216                            objects[n] = new OsiLotsize(testOsiSolver,iColumn,2,ranges,true);
4217                          }
4218                          objects[n++]->setPriority (priority-1);
4219                        }
4220                      }
4221                    }
4222                    if (testOsiOptions<0) {
4223                      babModel->addObjects(n,objects);
4224                    } else {
4225                      testOsiSolver->addObjects(n,objects);
4226                    }
4227                    for (i=0;i<n;i++)
4228                      delete objects[i];
4229                    delete [] objects;
4230                  }
4231                }
4232#endif
4233                if (storedAmpl.sizeRowCuts()) {
4234                  //babModel->addCutGenerator(&storedAmpl,1,"AmplStored");
4235                  int numberRowCuts = storedAmpl.sizeRowCuts();
4236                  for (int i=0;i<numberRowCuts;i++) {
4237                    const OsiRowCut * rowCutPointer = storedAmpl.rowCutPointer(i);
4238                    babModel->makeGlobalCut(rowCutPointer);
4239                  }
4240                }
4241                babModel->branchAndBound(statistics);
4242#ifdef CLP_MALLOC_STATISTICS
4243                malloc_stats();
4244                malloc_stats2();
4245#endif
4246                checkSOS(babModel, babModel->solver());
4247              } else if (type==MIPLIB) {
4248                CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
4249                // Set up pre-processing
4250                int translate2[]={9999,1,1,3,2,4,5};
4251                if (preProcess)
4252                  strategy.setupPreProcessing(translate2[preProcess]);
4253                babModel->setStrategy(strategy);
4254                if (outputFormat==5) {
4255                  osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4256                  lpSolver = osiclp->getModelPtr();
4257                  lpSolver->setPersistenceFlag(1);
4258                }
4259                if (testOsiOptions>=0) {
4260                  printf("Testing OsiObject options %d\n",testOsiOptions);
4261                  CbcBranchDefaultDecision decision;
4262                  OsiChooseStrong choose(babModel->solver());
4263                  choose.setNumberBeforeTrusted(babModel->numberBeforeTrust());
4264                  choose.setNumberStrong(babModel->numberStrong());
4265                  choose.setShadowPriceMode(testOsiOptions);
4266                  //babModel->deleteObjects(false);
4267                  decision.setChooseMethod(choose);
4268                  babModel->setBranchingMethod(decision);
4269                }
4270                *miplibSolver = babModel;
4271                return 777;
4272              } else {
4273                strengthenedModel = babModel->strengthenedModel();
4274              }
4275              currentBranchModel = NULL;
4276              osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4277              if (debugFile=="createAfterPre"&&babModel->bestSolution()) {
4278                lpSolver = osiclp->getModelPtr();
4279                //move best solution (should be there -- but ..)
4280                int n = lpSolver->getNumCols();
4281                memcpy(lpSolver->primalColumnSolution(),babModel->bestSolution(),n*sizeof(double));
4282                saveSolution(osiclp->getModelPtr(),"debug.file");
4283              }
4284              if (!noPrinting) {
4285                // Print more statistics
4286                std::cout<<"Cuts at root node changed objective from "<<babModel->getContinuousObjective()
4287                         <<" to "<<babModel->rootObjectiveAfterCuts()<<std::endl;
4288               
4289                numberGenerators = babModel->numberCutGenerators();
4290                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
4291                  CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
4292                  std::cout<<generator->cutGeneratorName()<<" was tried "
4293                           <<generator->numberTimesEntered()<<" times and created "
4294                           <<generator->numberCutsInTotal()<<" cuts of which "
4295                           <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
4296                  if (generator->timing())
4297                    std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
4298                  else
4299                    std::cout<<std::endl;
4300                }
4301              }
4302              time2 = CoinCpuTime();
4303              totalTime += time2-time1;
4304              // For best solution
4305              double * bestSolution = NULL;
4306              if (babModel->getMinimizationObjValue()<1.0e50&&type==BAB) {
4307                // post process
4308                if (preProcess) {
4309                  int n = saveSolver->getNumCols();
4310                  bestSolution = new double [n];
4311                  process.postProcess(*babModel->solver());
4312                  // Solution now back in saveSolver
4313                  babModel->assignSolver(saveSolver);
4314                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
4315                } else {
4316                  int n = babModel->solver()->getNumCols();
4317                  bestSolution = new double [n];
4318                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
4319                }
4320                // and put back in very original solver
4321                {
4322                  ClpSimplex * original = originalSolver.getModelPtr();
4323                  double * lower = original->columnLower();
4324                  double * upper = original->columnUpper();
4325                  double * solution = original->primalColumnSolution();
4326                  int n = original->numberColumns();
4327                  assert (!n||n==babModel->solver()->getNumCols());
4328                  for (int i=0;i<n;i++) {
4329                    solution[i]=bestSolution[i];
4330                    if (solver1.isInteger(i)) {
4331                      lower[i]=solution[i];
4332                      upper[i]=solution[i];
4333                    }
4334                  }
4335                }
4336                checkSOS(babModel, babModel->solver());
4337              }
4338              if (type==STRENGTHEN&&strengthenedModel)
4339                clpSolver = dynamic_cast< OsiClpSolverInterface*> (strengthenedModel);
4340#ifdef COIN_HAS_ASL
4341              else if (usingAmpl) 
4342                clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4343#endif
4344              lpSolver = clpSolver->getModelPtr();
4345              if (numberChanged) {
4346                for (int i=0;i<numberChanged;i++) {
4347                  int iColumn=changed[i];
4348                  clpSolver->setContinuous(iColumn);
4349                }
4350                delete [] changed;
4351              }
4352              if (type==BAB) {
4353                //move best solution (should be there -- but ..)
4354                int n = lpSolver->getNumCols();
4355                if (bestSolution) {
4356                  memcpy(lpSolver->primalColumnSolution(),bestSolution,n*sizeof(double));
4357                  // now see what that does to row solution
4358                  int numberRows=lpSolver->numberRows();
4359                  double * rowSolution = lpSolver->primalRowSolution();
4360                  memset (rowSolution,0,numberRows*sizeof(double));
4361                  lpSolver->clpMatrix()->times(1.0,bestSolution,rowSolution);
4362                }
4363                if (debugFile=="create"&&bestSolution) {
4364                  saveSolution(lpSolver,"debug.file");
4365                }
4366                delete [] bestSolution;
4367                std::string statusName[]={"Finished","Stopped on ","Difficulties",
4368                                          "","","User ctrl-c"};
4369                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
4370                int iStat = babModel->status();
4371                int iStat2 = babModel->secondaryStatus();
4372                if (!noPrinting)
4373                  std::cout<<"Result - "<<statusName[iStat]<<minor[iStat2]
4374                           <<" objective "<<babModel->getObjValue()<<
4375                    " after "<<babModel->getNodeCount()<<" nodes and "
4376                           <<babModel->getIterationCount()<<
4377                    " iterations - took "<<time2-time1<<" seconds"<<std::endl;
4378#ifdef COIN_HAS_ASL
4379                if (usingAmpl) {
4380                  clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4381                  lpSolver = clpSolver->getModelPtr();
4382                  double value = babModel->getObjValue()*lpSolver->getObjSense();
4383                  char buf[300];
4384                  int pos=0;
4385                  if (iStat==0) {
4386                    if (babModel->getObjValue()<1.0e40) {
4387                      pos += sprintf(buf+pos,"optimal," );
4388                    } else {
4389                      // infeasible
4390                      iStat=1;
4391                      pos += sprintf(buf+pos,"infeasible,");
4392                    }
4393                  } else if (iStat==1) {
4394                    if (iStat2!=6)
4395                      iStat=3;
4396                    else
4397                      iStat=4;
4398                    pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
4399                  } else if (iStat==2) {
4400                    iStat = 7;
4401                    pos += sprintf(buf+pos,"stopped on difficulties,");
4402                  } else if (iStat==5) {
4403                    iStat = 3;
4404                    pos += sprintf(buf+pos,"stopped on ctrl-c,");
4405                  } else {
4406                    pos += sprintf(buf+pos,"status unknown,");
4407                    iStat=6;
4408                  }
4409                  info.problemStatus=iStat;
4410                  info.objValue = value;
4411                  if (babModel->getObjValue()<1.0e40) {
4412                    int precision = ampl_obj_prec();
4413                    if (precision>0)
4414                      pos += sprintf(buf+pos," objective %.*g",precision,
4415                                     value);
4416                    else
4417                      pos += sprintf(buf+pos," objective %g",value);
4418                  }
4419                  sprintf(buf+pos,"\n%d nodes, %d iterations, %g seconds",
4420                          babModel->getNodeCount(),
4421                          babModel->getIterationCount(),
4422                          totalTime);
4423                  if (bestSolution) {
4424                    free(info.primalSolution);
4425                    if (!numberKnapsack) {
4426                      info.primalSolution = (double *) malloc(n*sizeof(double));
4427                      CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
4428                      int numberRows = lpSolver->numberRows();
4429                      free(info.dualSolution);
4430                      info.dualSolution = (double *) malloc(numberRows*sizeof(double));
4431                      CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
4432                    } else {
4433                      // expanded knapsack
4434                      info.dualSolution=NULL;
4435                      int numberColumns = saveCoinModel.numberColumns();
4436                      info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
4437                      // Fills in original solution (coinModel length)
4438                      afterKnapsack(saveCoinModel,  whichColumn,  knapsackStart, 
4439                                    knapsackRow,  numberKnapsack,
4440                                    lpSolver->primalColumnSolution(), info.primalSolution,1);
4441                    }
4442                  } else {
4443                    info.primalSolution=NULL;
4444                    info.dualSolution=NULL;
4445                  }
4446                  // put buffer into info
4447                  strcpy(info.buffer,buf);
4448                }
4449#endif
4450              } else {
4451                std::cout<<"Model strengthened - now has "<<clpSolver->getNumRows()
4452                         <<" rows"<<std::endl;
4453              }
4454              time1 = time2;
4455              delete babModel;
4456              babModel=NULL;
4457            } else {
4458              std::cout << "** Current model not valid" << std::endl ; 
4459            }
4460            break ;
4461          case IMPORT:
4462            {
4463#ifdef COIN_HAS_ASL
4464              if (!usingAmpl) {
4465#endif
4466                free(priorities);
4467                priorities=NULL;
4468                free(branchDirection);
4469                branchDirection=NULL;
4470                free(pseudoDown);
4471                pseudoDown=NULL;
4472                free(pseudoUp);
4473                pseudoUp=NULL;
4474                free(solutionIn);
4475                solutionIn=NULL;
4476                free(prioritiesIn);
4477                prioritiesIn=NULL;
4478                free(sosStart);
4479                sosStart=NULL;
4480                free(sosIndices);
4481                sosIndices=NULL;
4482                free(sosType);
4483                sosType=NULL;
4484                free(sosReference);
4485                sosReference=NULL;
4486                free(cut);
4487                cut=NULL;
4488                free(sosPriority);
4489                sosPriority=NULL;
4490#ifdef COIN_HAS_ASL
4491              }
4492#endif               
4493              delete babModel;
4494              babModel=NULL;
4495              // get next field
4496              field = CoinReadGetString(argc,argv);
4497              if (field=="$") {
4498                field = parameters[iParam].stringValue();
4499              } else if (field=="EOL") {
4500                parameters[iParam].printString();
4501                break;
4502              } else {
4503                parameters[iParam].setStringValue(field);
4504              }
4505              std::string fileName;
4506              bool canOpen=false;
4507              // See if gmpl file
4508              int gmpl=0;
4509              std::string gmplData;
4510              if (field=="-") {
4511                // stdin
4512                canOpen=true;
4513                fileName = "-";
4514              } else {
4515                // See if .lp
4516                {
4517                  const char * c_name = field.c_str();
4518                  int length = strlen(c_name);
4519                  if (length>3&&!strncmp(c_name+length-3,".lp",3))
4520                    gmpl=-1; // .lp
4521                }
4522                bool absolutePath;
4523                if (dirsep=='/') {
4524                  // non Windows (or cygwin)
4525                  absolutePath=(field[0]=='/');
4526                } else {
4527                  //Windows (non cycgwin)
4528                  absolutePath=(field[0]=='\\');
4529                  // but allow for :
4530                  if (strchr(field.c_str(),':'))
4531                    absolutePath=true;
4532                }
4533                if (absolutePath) {
4534                  fileName = field;
4535                } else if (field[0]=='~') {
4536                  char * environVar = getenv("HOME");
4537                  if (environVar) {
4538                    std::string home(environVar);
4539                    field=field.erase(0,1);
4540                    fileName = home+field;
4541                  } else {
4542                    fileName=field;
4543                  }
4544                } else {
4545                  fileName = directory+field;
4546                  // See if gmpl (model & data) - or even lp file
4547                  int length = field.size();
4548                  int percent = field.find('%');
4549                  if (percent<length&&percent>0) {
4550                    gmpl=1;
4551                    fileName = directory+field.substr(0,percent);
4552                    gmplData = directory+field.substr(percent+1);
4553                    if (percent<length-1)
4554                      gmpl=2; // two files
4555                    printf("GMPL model file %s and data file %s\n",
4556                           fileName.c_str(),gmplData.c_str());
4557                  }
4558                }
4559                std::string name=fileName;
4560                if (fileCoinReadable(name)) {
4561                  // can open - lets go for it
4562                  canOpen=true;
4563                  if (gmpl==2) {
4564                    FILE *fp;
4565                    fp=fopen(gmplData.c_str(),"r");
4566                    if (fp) {
4567                      fclose(fp);
4568                    } else {
4569                      canOpen=false;
4570                      std::cout<<"Unable to open file "<<gmplData<<std::endl;
4571                    }
4572                  }
4573                } else {
4574                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4575                }
4576              }
4577              if (canOpen) {
4578                int status;
4579                ClpSimplex * lpSolver = clpSolver->getModelPtr();
4580                if (!gmpl)
4581                  status =clpSolver->readMps(fileName.c_str(),
4582                                                 keepImportNames!=0,
4583                                                 allowImportErrors!=0);
4584                else if (gmpl>0)
4585                  status= lpSolver->readGMPL(fileName.c_str(),
4586                                                  (gmpl==2) ? gmplData.c_str() : NULL,
4587                                                  keepImportNames!=0);
4588                else
4589                  status= lpSolver->readLp(fileName.c_str(),1.0e-12);
4590                if (!status||(status>0&&allowImportErrors)) {
4591                  goodModel=true;
4592                  // sets to all slack (not necessary?)
4593                  lpSolver->createStatus();
4594                  // make sure integer
4595                  int numberColumns = lpSolver->numberColumns();
4596                  for (int i=0;i<numberColumns;i++) {
4597                    if (lpSolver->isInteger(i))
4598                      clpSolver->setInteger(i);
4599                  }
4600                  time2 = CoinCpuTime();
4601                  totalTime += time2-time1;
4602                  time1=time2;
4603                  // Go to canned file if just input file
4604                  if (CbcOrClpRead_mode==2&&argc==2) {
4605                    // only if ends .mps
4606                    char * find = (char *)strstr(fileName.c_str(),".mps");
4607                    if (find&&find[4]=='\0') {
4608                      find[1]='p'; find[2]='a';find[3]='r';
4609                      FILE *fp=fopen(fileName.c_str(),"r");
4610                      if (fp) {
4611                        CbcOrClpReadCommand=fp; // Read from that file
4612                        CbcOrClpRead_mode=-1;
4613                      }
4614                    }
4615                  }
4616                } else {
4617                  // errors
4618                  std::cout<<"There were "<<status<<
4619                    " errors on input"<<std::endl;
4620                }
4621              }
4622            }
4623            break;
4624          case MODELIN:
4625#ifdef COIN_HAS_LINK
4626            {
4627              // get next field
4628              field = CoinReadGetString(argc,argv);
4629              if (field=="$") {
4630                field = parameters[iParam].stringValue();
4631              } else if (field=="EOL") {
4632                parameters[iParam].printString();
4633                break;
4634              } else {
4635                parameters[iParam].setStringValue(field);
4636              }
4637              std::string fileName;
4638              bool canOpen=false;
4639              if (field=="-") {
4640                // stdin
4641                canOpen=true;
4642                fileName = "-";
4643              } else {
4644                bool absolutePath;
4645                if (dirsep=='/') {
4646                  // non Windows (or cygwin)
4647                  absolutePath=(field[0]=='/');
4648                } else {
4649                  //Windows (non cycgwin)
4650                  absolutePath=(field[0]=='\\');
4651                  // but allow for :
4652                  if (strchr(field.c_str(),':'))
4653                    absolutePath=true;
4654                }
4655                if (absolutePath) {
4656                  fileName = field;
4657                } else if (field[0]=='~') {
4658                  char * environVar = getenv("HOME");
4659                  if (environVar) {
4660                    std::string home(environVar);
4661                    field=field.erase(0,1);
4662                    fileName = home+field;
4663                  } else {
4664                    fileName=field;
4665                  }
4666                } else {
4667                  fileName = directory+field;
4668                }
4669                FILE *fp=fopen(fileName.c_str(),"r");
4670                if (fp) {
4671                  // can open - lets go for it
4672                  fclose(fp);
4673                  canOpen=true;
4674                } else {
4675                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4676                }
4677              }
4678              if (canOpen) {
4679                CoinModel coinModel(fileName.c_str(),2);
4680                // load from coin model
4681                OsiSolverLink solver1;
4682                OsiSolverInterface * solver2 = solver1.clone();
4683                model.assignSolver(solver2,true);
4684                OsiSolverLink * si =
4685                  dynamic_cast<OsiSolverLink *>(model.solver()) ;
4686                assert (si != NULL);
4687                si->setDefaultMeshSize(0.001);
4688                // need some relative granularity
4689                si->setDefaultBound(100.0);
4690                si->setDefaultMeshSize(0.01);
4691                si->setDefaultBound(100.0);
4692                si->setIntegerPriority(1000);
4693                si->setBiLinearPriority(10000);
4694                CoinModel * model2 = (CoinModel *) &coinModel;
4695                si->load(*model2);
4696                // redo
4697                solver = model.solver();
4698                clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
4699                lpSolver = clpSolver->getModelPtr();
4700                clpSolver->messageHandler()->setLogLevel(0) ;
4701                testOsiParameters=0;
4702                complicatedInteger=2;
4703              }
4704            }
4705#endif
4706            break;
4707          case EXPORT:
4708            if (goodModel) {
4709              // get next field
4710              field = CoinReadGetString(argc,argv);
4711              if (field=="$") {
4712                field = parameters[iParam].stringValue();
4713              } else if (field=="EOL") {
4714                parameters[iParam].printString();
4715                break;
4716              } else {
4717                parameters[iParam].setStringValue(field);
4718              }
4719              std::string fileName;
4720              bool canOpen=false;
4721              if (field[0]=='/'||field[0]=='\\') {
4722                fileName = field;
4723              } else if (field[0]=='~') {
4724                char * environVar = getenv("HOME");
4725                if (environVar) {
4726                  std::string home(environVar);
4727                  field=field.erase(0,1);
4728                  fileName = home+field;
4729                } else {
4730                  fileName=field;
4731                }
4732              } else {
4733                fileName = directory+field;
4734              }
4735              FILE *fp=fopen(fileName.c_str(),"w");
4736              if (fp) {
4737                // can open - lets go for it
4738                fclose(fp);
4739                canOpen=true;
4740              } else {
4741                std::cout<<"Unable to open file "<<fileName<<std::endl;
4742              }
4743              if (canOpen) {
4744                // If presolve on then save presolved
4745                bool deleteModel2=false;
4746                ClpSimplex * model2 = lpSolver;
4747                if (dualize) {
4748                  model2 = ((ClpSimplexOther *) model2)->dualOfModel();
4749                  printf("Dual of model has %d rows and %d columns\n",
4750                         model2->numberRows(),model2->numberColumns());
4751                  model2->setOptimizationDirection(1.0);
4752                }
4753#ifdef COIN_HAS_ASL
4754                if (info.numberSos&&doSOS&&usingAmpl) {
4755                  // SOS
4756                  numberSOS = info.numberSos;
4757                  sosStart = info.sosStart;
4758                  sosIndices = info.sosIndices;
4759                  sosReference = info.sosReference;
4760                  preSolve=false;
4761                  clpSolver->setSOSData(numberSOS,info.sosType,sosStart,sosIndices,sosReference);
4762                }
4763#endif
4764                if (preSolve) {
4765                  ClpPresolve pinfo;
4766                  int presolveOptions2 = presolveOptions&~0x40000000;
4767                  if ((presolveOptions2&0xffff)!=0)
4768                    pinfo.setPresolveActions(presolveOptions2);
4769                  if ((printOptions&1)!=0)
4770                    pinfo.statistics();
4771                  double presolveTolerance = 
4772                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
4773                  model2 = 
4774                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
4775                                         true,preSolve);
4776                  if (model2) {
4777                    printf("Saving presolved model on %s\n",
4778                           fileName.c_str());
4779                    deleteModel2=true;
4780                  } else {
4781                    printf("Presolved model looks infeasible - saving original on %s\n",
4782                           fileName.c_str());
4783                    deleteModel2=false;
4784                    model2 = lpSolver;
4785
4786                  }
4787                  model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
4788                  if (deleteModel2)
4789                    delete model2;
4790                } else {
4791                  printf("Saving model on %s\n",
4792                           fileName.c_str());
4793                  if (numberSOS) {
4794                    // Convert names
4795                    int iRow;
4796                    int numberRows=model2->numberRows();
4797                    int iColumn;
4798                    int numberColumns=model2->numberColumns();
4799                   
4800                    char ** rowNames = NULL;
4801                    char ** columnNames = NULL;
4802                    if (model2->lengthNames()) {
4803                      rowNames = new char * [numberRows];
4804                      for (iRow=0;iRow<numberRows;iRow++) {
4805                        rowNames[iRow] = 
4806                          strdup(model2->rowName(iRow).c_str());
4807                      }
4808                     
4809                      columnNames = new char * [numberColumns];
4810                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
4811                        columnNames[iColumn] = 
4812                          strdup(model2->columnName(iColumn).c_str());
4813                      }
4814                    }
4815                    clpSolver->writeMpsNative(fileName.c_str(),(const char **) rowNames,(const char **) columnNames,
4816                                              (outputFormat-1)/2,1+((outputFormat-1)&1));
4817                    if (rowNames) {
4818                      for (iRow=0;iRow<numberRows;iRow++) {
4819                        free(rowNames[iRow]);
4820                      }
4821                      delete [] rowNames;
4822                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
4823                        free(columnNames[iColumn]);
4824                      }
4825                      delete [] columnNames;
4826                    }
4827                  } else {
4828                    model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
4829                  }
4830                }
4831                time2 = CoinCpuTime();
4832                totalTime += time2-time1;
4833                time1=time2;
4834              }
4835            } else {
4836              std::cout<<"** Current model not valid"<<std::endl;
4837            }
4838            break;
4839          case BASISIN:
4840            if (goodModel) {
4841              // get next field
4842              field = CoinReadGetString(argc,argv);
4843              if (field=="$") {
4844                field = parameters[iParam].stringValue();
4845              } else if (field=="EOL") {
4846                parameters[iParam].printString();
4847                break;
4848              } else {
4849                parameters[iParam].setStringValue(field);
4850              }
4851              std::string fileName;
4852              bool canOpen=false;
4853              if (field=="-") {
4854                // stdin
4855                canOpen=true;
4856                fileName = "-";
4857              } else {
4858                if (field[0]=='/'||field[0]=='\\') {
4859                  fileName = field;
4860                } else if (field[0]=='~') {
4861                  char * environVar = getenv("HOME");
4862                  if (environVar) {
4863                    std::string home(environVar);
4864                    field=field.erase(0,1);
4865                    fileName = home+field;
4866                  } else {
4867                    fileName=field;
4868                  }
4869                } else {
4870                  fileName = directory+field;
4871                }
4872                FILE *fp=fopen(fileName.c_str(),"r");
4873                if (fp) {
4874                  // can open - lets go for it
4875                  fclose(fp);
4876                  canOpen=true;
4877                } else {
4878                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4879                }
4880              }
4881              if (canOpen) {
4882                int values = lpSolver->readBasis(fileName.c_str());
4883                if (values==0)
4884                  basisHasValues=-1;
4885                else
4886                  basisHasValues=1;
4887              }
4888            } else {
4889              std::cout<<"** Current model not valid"<<std::endl;
4890            }
4891            break;
4892          case PRIORITYIN:
4893            if (goodModel) {
4894              // get next field
4895              field = CoinReadGetString(argc,argv);
4896              if (field=="$") {
4897                field = parameters[iParam].stringValue();
4898              } else if (field=="EOL") {
4899                parameters[iParam].printString();
4900                break;
4901              } else {
4902                parameters[iParam].setStringValue(field);
4903              }
4904              std::string fileName;
4905              if (field[0]=='/'||field[0]=='\\') {
4906                fileName = field;
4907              } else if (field[0]=='~') {
4908                char * environVar = getenv("HOME");
4909                if (environVar) {
4910                  std::string home(environVar);
4911                  field=field.erase(0,1);
4912                  fileName = home+field;
4913                } else {
4914                  fileName=field;
4915                }
4916              } else {
4917                fileName = directory+field;
4918              }
4919              FILE *fp=fopen(fileName.c_str(),"r");
4920              if (fp) {
4921                // can open - lets go for it
4922                std::string headings[]={"name","number","direction","priority","up","down",
4923                                        "solution","priin"};
4924                int got[]={-1,-1,-1,-1,-1,-1,-1,-1};
4925                int order[8];
4926                assert(sizeof(got)==sizeof(order));
4927                int nAcross=0;
4928                char line[1000];
4929                int numberColumns = lpSolver->numberColumns();
4930                if (!fgets(line,1000,fp)) {
4931                  std::cout<<"Odd file "<<fileName<<std::endl;
4932                } else {
4933                  char * pos = line;
4934                  char * put = line;
4935                  while (*pos>=' '&&*pos!='\n') {
4936                    if (*pos!=' '&&*pos!='\t') {
4937                      *put=tolower(*pos);
4938                      put++;
4939                    }
4940                    pos++;
4941                  }
4942                  *put='\0';
4943                  pos=line;
4944                  int i;
4945                  bool good=true;
4946                  while (pos) {
4947                    char * comma = strchr(pos,',');
4948                    if (comma)
4949                      *comma='\0';
4950                    for (i=0;i<(int) (sizeof(got)/sizeof(int));i++) {
4951                      if (headings[i]==pos) {
4952                        if (got[i]<0) {
4953                          order[nAcross]=i;
4954                          got[i]=nAcross++;
4955                        } else {
4956                          // duplicate
4957                          good=false;
4958                        }
4959                        break;
4960                      }
4961                    }
4962                    if (i==(int) (sizeof(got)/sizeof(int)))
4963                      good=false;
4964                    if (comma) {
4965                      *comma=',';
4966                      pos=comma+1;
4967                    } else {
4968                      break;
4969                    }
4970                  }
4971                  if (got[0]<0&&got[1]<0)
4972                    good=false;
4973                  if (got[0]>=0&&got[1]>=0)
4974                    good=false;
4975                  if (got[0]>=0&&!lpSolver->lengthNames())
4976                    good=false;
4977                  if (good) {
4978                    char ** columnNames = new char * [numberColumns];
4979                    pseudoDown= (double *) malloc(numberColumns*sizeof(double));
4980                    pseudoUp = (double *) malloc(numberColumns*sizeof(double));
4981                    branchDirection = (int *) malloc(numberColumns*sizeof(int));
4982                    priorities= (int *) malloc(numberColumns*sizeof(int));
4983                    free(solutionIn);
4984                    solutionIn=NULL;
4985                    free(prioritiesIn);
4986                    prioritiesIn=NULL;
4987                    int iColumn;
4988                    if (got[6]>=0) {
4989                      solutionIn = (double *) malloc(numberColumns*sizeof(double));
4990                      CoinZeroN(solutionIn,numberColumns);
4991                    }
4992                    if (got[7]>=0) {
4993                      prioritiesIn = (int *) malloc(numberColumns*sizeof(int));
4994                      for (iColumn=0;iColumn<numberColumns;iColumn++) 
4995                        prioritiesIn[iColumn]=10000;
4996                    }
4997                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4998                      columnNames[iColumn] = 
4999                        strdup(lpSolver->columnName(iColumn).c_str());
5000                      pseudoDown[iColumn]=0.0;
5001                      pseudoUp[iColumn]=0.0;
5002                      branchDirection[iColumn]=0;
5003                      priorities[iColumn]=0;
5004                    }
5005                    int nBadPseudo=0;
5006