source: branches/devel/Cbc/src/CoinSolve.cpp @ 579

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

add time and try two level bilinear

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