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

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

to compile without ampl

  • Property svn:eol-style set to native
  • Property svn:keywords set to Author Date Id Revision
File size: 222.8 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#ifdef COIN_HAS_LINK
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#endif
2646              if (!miplib) {
2647                if (!preSolve) {
2648                  model.solver()->setHintParam(OsiDoPresolveInInitial,false,OsiHintTry);
2649                  model.solver()->setHintParam(OsiDoPresolveInResolve,false,OsiHintTry);
2650                }
2651                model.initialSolve();
2652                OsiSolverInterface * solver = model.solver();
2653                OsiClpSolverInterface * si =
2654                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2655                ClpSimplex * clpSolver = si->getModelPtr();
2656                if (!complicatedInteger&&clpSolver->tightenPrimalBounds()!=0) {
2657                  std::cout<<"Problem is infeasible - tightenPrimalBounds!"<<std::endl;
2658                  exit(1);
2659                }
2660                if (clpSolver->dualBound()==1.0e10) {
2661                  // user did not set - so modify
2662                  // get largest scaled away from bound
2663                  double largest=1.0e-12;
2664                  int numberRows = clpSolver->numberRows();
2665                  const double * rowPrimal = clpSolver->primalRowSolution();
2666                  const double * rowLower = clpSolver->rowLower();
2667                  const double * rowUpper = clpSolver->rowUpper();
2668                  const double * rowScale = clpSolver->rowScale();
2669                  int iRow;
2670                  for (iRow=0;iRow<numberRows;iRow++) {
2671                    double value = rowPrimal[iRow];
2672                    double above = value-rowLower[iRow];
2673                    double below = rowUpper[iRow]-value;
2674                    if (rowScale) {
2675                      double multiplier = rowScale[iRow];
2676                      above *= multiplier;
2677                      below *= multiplier;
2678                    }
2679                    if (above<1.0e12)
2680                      largest = CoinMax(largest,above);
2681                    if (below<1.0e12)
2682                      largest = CoinMax(largest,below);
2683                  }
2684                 
2685                  int numberColumns = clpSolver->numberColumns();
2686                  const double * columnPrimal = clpSolver->primalColumnSolution();
2687                  const double * columnLower = clpSolver->columnLower();
2688                  const double * columnUpper = clpSolver->columnUpper();
2689                  const double * columnScale = clpSolver->columnScale();
2690                  int iColumn;
2691                  for (iColumn=0;iColumn<numberColumns;iColumn++) {
2692                    double value = columnPrimal[iColumn];
2693                    double above = value-columnLower[iColumn];
2694                    double below = columnUpper[iColumn]-value;
2695                    if (columnScale) {
2696                      double multiplier = 1.0/columnScale[iColumn];
2697                      above *= multiplier;
2698                      below *= multiplier;
2699                    }
2700                    if (above<1.0e12)
2701                      largest = CoinMax(largest,above);
2702                    if (below<1.0e12)
2703                      largest = CoinMax(largest,below);
2704                  }
2705                  //if (!noPrinting)
2706                  //std::cout<<"Largest (scaled) away from bound "<<largest<<std::endl;
2707                  clpSolver->setDualBound(CoinMax(1.0001e8,CoinMin(1000.0*largest,1.00001e10)));
2708                }
2709                clpSolver->dual();  // clean up
2710              }
2711              // If user made settings then use them
2712              if (!defaultSettings) {
2713                OsiSolverInterface * solver = model.solver();
2714                if (!doScaling)
2715                  solver->setHintParam(OsiDoScale,false,OsiHintTry);
2716                OsiClpSolverInterface * si =
2717                  dynamic_cast<OsiClpSolverInterface *>(solver) ;
2718                assert (si != NULL);
2719                // get clp itself
2720                ClpSimplex * modelC = si->getModelPtr();
2721                //if (modelC->tightenPrimalBounds()!=0) {
2722                //std::cout<<"Problem is infeasible!"<<std::endl;
2723                //break;
2724                //}
2725                // bounds based on continuous
2726                if (tightenFactor&&!complicatedInteger) {
2727                  if (modelC->tightenPrimalBounds(tightenFactor)!=0) {
2728                    std::cout<<"Problem is infeasible!"<<std::endl;
2729                    break;
2730                  }
2731                }
2732                if (djFix<1.0e20) {
2733                  // do some fixing
2734                  int numberColumns = modelC->numberColumns();
2735                  int i;
2736                  const char * type = modelC->integerInformation();
2737                  double * lower = modelC->columnLower();
2738                  double * upper = modelC->columnUpper();
2739                  double * solution = modelC->primalColumnSolution();
2740                  double * dj = modelC->dualColumnSolution();
2741                  int numberFixed=0;
2742                  for (i=0;i<numberColumns;i++) {
2743                    if (type[i]) {
2744                      double value = solution[i];
2745                      if (value<lower[i]+1.0e-5&&dj[i]>djFix) {
2746                        solution[i]=lower[i];
2747                        upper[i]=lower[i];
2748                        numberFixed++;
2749                      } else if (value>upper[i]-1.0e-5&&dj[i]<-djFix) {
2750                        solution[i]=upper[i];
2751                        lower[i]=upper[i];
2752                        numberFixed++;
2753                      }
2754                    }
2755                  }
2756                  printf("%d columns fixed\n",numberFixed);
2757                }
2758              }
2759              // See if we want preprocessing
2760              OsiSolverInterface * saveSolver=NULL;
2761              CglPreProcess process;
2762              delete babModel;
2763              babModel = new CbcModel(model);
2764              OsiSolverInterface * solver3 = clpSolver->clone();
2765              babModel->assignSolver(solver3);
2766              OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
2767              int numberChanged=0;
2768              if (clpSolver2->messageHandler()->logLevel())
2769                clpSolver2->messageHandler()->setLogLevel(1);
2770              if (logLevel>-1)
2771                clpSolver2->messageHandler()->setLogLevel(logLevel);
2772              lpSolver = clpSolver2->getModelPtr();
2773              if (lpSolver->factorizationFrequency()==200&&!miplib) {
2774                // User did not touch preset
2775                int numberRows = lpSolver->numberRows();
2776                const int cutoff1=10000;
2777                const int cutoff2=100000;
2778                const int base=75;
2779                const int freq0 = 50;
2780                const int freq1=200;
2781                const int freq2=400;
2782                const int maximum=1000;
2783                int frequency;
2784                if (numberRows<cutoff1)
2785                  frequency=base+numberRows/freq0;
2786                else if (numberRows<cutoff2)
2787                  frequency=base+cutoff1/freq0 + (numberRows-cutoff1)/freq1;
2788                else
2789                  frequency=base+cutoff1/freq0 + (cutoff2-cutoff1)/freq1 + (numberRows-cutoff2)/freq2;
2790                lpSolver->setFactorizationFrequency(CoinMin(maximum,frequency));
2791              }
2792              time2 = CoinCpuTime();
2793              totalTime += time2-time1;
2794              time1 = time2;
2795              double timeLeft = babModel->getMaximumSeconds();
2796              int numberOriginalColumns = babModel->solver()->getNumCols();
2797              if (preProcess==7) {
2798                // use strategy instead
2799                preProcess=0;
2800                useStrategy=true;
2801#ifdef COIN_HAS_LINK
2802                // empty out any cuts
2803                if (storedAmpl.sizeRowCuts()) {
2804                  printf("Emptying ampl stored cuts as internal preprocessing\n");
2805                  CglStored temp;
2806                  storedAmpl=temp;
2807                }
2808#endif
2809              }
2810              if (preProcess&&type==BAB) {
2811                // See if sos from mps file
2812                if (numberSOS==0&&clpSolver->numberSOS()&&doSOS) {
2813                  // SOS
2814                  numberSOS = clpSolver->numberSOS();
2815                  const CoinSet * setInfo = clpSolver->setInfo();
2816                  sosStart = new int [numberSOS+1];
2817                  sosType = new char [numberSOS];
2818                  int i;
2819                  int nTotal=0;
2820                  sosStart[0]=0;
2821                  for ( i=0;i<numberSOS;i++) {
2822                    int type = setInfo[i].setType();
2823                    int n=setInfo[i].numberEntries();
2824                    sosType[i]=type;
2825                    nTotal += n;
2826                    sosStart[i+1] = nTotal;
2827                  }
2828                  sosIndices = new int[nTotal];
2829                  sosReference = new double [nTotal];
2830                  for (i=0;i<numberSOS;i++) {
2831                    int n=setInfo[i].numberEntries();
2832                    const int * which = setInfo[i].which();
2833                    const double * weights = setInfo[i].weights();
2834                    int base = sosStart[i];
2835                    for (int j=0;j<n;j++) {
2836                      int k=which[j];
2837                      sosIndices[j+base]=k;
2838                      sosReference[j+base] = weights ? weights[j] : (double) j;
2839                    }
2840                  }
2841                }
2842                saveSolver=babModel->solver()->clone();
2843                /* Do not try and produce equality cliques and
2844                   do up to 10 passes */
2845                OsiSolverInterface * solver2;
2846                {
2847                  // Tell solver we are in Branch and Cut
2848                  saveSolver->setHintParam(OsiDoInBranchAndCut,true,OsiHintDo) ;
2849                  // Default set of cut generators
2850                  CglProbing generator1;
2851                  generator1.setUsingObjective(1);
2852                  generator1.setMaxPass(3);
2853                  generator1.setMaxProbeRoot(saveSolver->getNumCols());
2854                  generator1.setMaxElements(100);
2855                  generator1.setMaxLookRoot(50);
2856                  generator1.setRowCuts(3);
2857                  // Add in generators
2858                  process.addCutGenerator(&generator1);
2859                  int translate[]={9999,0,0,-1,2,3,-2};
2860                  process.messageHandler()->setLogLevel(babModel->logLevel());
2861#ifdef COIN_HAS_ASL
2862                  if (info.numberSos&&doSOS&&usingAmpl) {
2863                    // SOS
2864                    numberSOS = info.numberSos;
2865                    sosStart = info.sosStart;
2866                    sosIndices = info.sosIndices;
2867                  }
2868#endif
2869                  if (numberSOS&&doSOS) {
2870                    // SOS
2871                    int numberColumns = saveSolver->getNumCols();
2872                    char * prohibited = new char[numberColumns];
2873                    memset(prohibited,0,numberColumns);
2874                    int n=sosStart[numberSOS];
2875                    for (int i=0;i<n;i++) {
2876                      int iColumn = sosIndices[i];
2877                      prohibited[iColumn]=1;
2878                    }
2879                    process.passInProhibited(prohibited,numberColumns);
2880                    delete [] prohibited;
2881                  }
2882                  int numberPasses = 10;
2883                  if (tunePreProcess>=1000) {
2884                    numberPasses = (tunePreProcess/1000)-1;
2885                    tunePreProcess = tunePreProcess % 1000;
2886                  }
2887                  if (doSprint>0) {
2888                    // Sprint for primal solves
2889                    ClpSolve::SolveType method = ClpSolve::usePrimalorSprint;
2890                    ClpSolve::PresolveType presolveType = ClpSolve::presolveOff;
2891                    int numberPasses = 5;
2892                    int options[] = {0,3,0,0,0,0};
2893                    int extraInfo[] = {-1,20,-1,-1,-1,-1};
2894                    extraInfo[1]=doSprint;
2895                    int independentOptions[] = {0,0,3};
2896                    ClpSolve clpSolve(method,presolveType,numberPasses,
2897                                      options,extraInfo,independentOptions);
2898                    // say use in OsiClp
2899                    clpSolve.setSpecialOption(6,1);
2900                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (saveSolver);
2901                    osiclp->setSolveOptions(clpSolve);
2902                    osiclp->setHintParam(OsiDoDualInResolve,false);
2903                    // switch off row copy
2904                    osiclp->getModelPtr()->setSpecialOptions(osiclp->getModelPtr()->specialOptions()|256);
2905                    osiclp->getModelPtr()->setInfeasibilityCost(1.0e11);
2906                  }
2907                  solver2 = process.preProcessNonDefault(*saveSolver,translate[preProcess],numberPasses,
2908                                                         tunePreProcess);
2909                  // Tell solver we are not in Branch and Cut
2910                  saveSolver->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
2911                  if (solver2)
2912                    solver2->setHintParam(OsiDoInBranchAndCut,false,OsiHintDo) ;
2913                }
2914#ifdef COIN_HAS_ASL
2915                if (!solver2&&usingAmpl) {
2916                  // infeasible
2917                  info.problemStatus=1;
2918                  info.objValue = 1.0e100;
2919                  sprintf(info.buffer,"infeasible/unbounded by pre-processing");
2920                  info.primalSolution=NULL;
2921                  info.dualSolution=NULL;
2922                  break;
2923                }
2924#endif
2925                if (!noPrinting) {
2926                  if (!solver2) {
2927                    printf("Pre-processing says infeasible or unbounded\n");
2928                    break;
2929                  } else {
2930                    printf("processed model has %d rows, %d columns and %d elements\n",
2931                           solver2->getNumRows(),solver2->getNumCols(),solver2->getNumElements());
2932                  }
2933                }
2934                //solver2->resolve();
2935                if (preProcess==2) {
2936                  OsiClpSolverInterface * clpSolver2 = dynamic_cast< OsiClpSolverInterface*> (solver2);
2937                  ClpSimplex * lpSolver = clpSolver2->getModelPtr();
2938                  lpSolver->writeMps("presolved.mps",0,1,lpSolver->optimizationDirection());
2939                  printf("Preprocessed model (minimization) on presolved.mps\n");
2940                }
2941                // we have to keep solver2 so pass clone
2942                solver2 = solver2->clone();
2943                babModel->assignSolver(solver2);
2944                babModel->initialSolve();
2945                babModel->setMaximumSeconds(timeLeft-(CoinCpuTime()-time1));
2946              }
2947              // now tighten bounds
2948              if (!miplib) {
2949                OsiClpSolverInterface * si =
2950                  dynamic_cast<OsiClpSolverInterface *>(babModel->solver()) ;
2951                assert (si != NULL);
2952                // get clp itself
2953                ClpSimplex * modelC = si->getModelPtr();
2954                if (noPrinting)
2955                  modelC->setLogLevel(0);
2956                if (!complicatedInteger&&modelC->tightenPrimalBounds()!=0) {
2957                  std::cout<<"Problem is infeasible!"<<std::endl;
2958                  break;
2959                }
2960                modelC->dual();
2961              }
2962#if 0
2963              numberDebugValues=599;
2964              debugValues = new double[numberDebugValues];
2965              CoinZeroN(debugValues,numberDebugValues);
2966              debugValues[3]=1.0;
2967              debugValues[6]=25.0;
2968              debugValues[9]=4.0;
2969              debugValues[26]=4.0;
2970              debugValues[27]=6.0;
2971              debugValues[35]=8.0;
2972              debugValues[53]=21.0;
2973              debugValues[56]=4.0;
2974#endif
2975              if (debugValues) {
2976                // for debug
2977                std::string problemName ;
2978                babModel->solver()->getStrParam(OsiProbName,problemName) ;
2979                babModel->solver()->activateRowCutDebugger(problemName.c_str()) ;
2980                twomirGen.probname_=strdup(problemName.c_str());
2981                // checking seems odd
2982                //redsplitGen.set_given_optsol(babModel->solver()->getRowCutDebuggerAlways()->optimalSolution(),
2983                //                         babModel->getNumCols());
2984              }
2985              int testOsiOptions = parameters[whichParam(TESTOSI,numberParameters,parameters)].intValue();
2986#ifdef COIN_HAS_LINK
2987              // If linked then see if expansion wanted
2988              {
2989                OsiSolverLink * solver3 = dynamic_cast<OsiSolverLink *> (babModel->solver());
2990                if (solver3) {
2991                  int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000;
2992                  if (options) {
2993                    /*
2994                      1 - force mini branch and bound
2995                      2 - set priorities high on continuous
2996                      4 - try adding OA cuts
2997                      8 - try doing quadratic linearization
2998                      16 - try expanding knapsacks
2999                    */
3000                    if ((options&16)) {
3001                      int numberColumns = saveCoinModel.numberColumns();
3002                      int numberRows = saveCoinModel.numberRows();
3003                      whichColumn = new int[numberColumns];
3004                      knapsackStart=new int[numberRows+1];
3005                      knapsackRow=new int[numberRows];
3006                      numberKnapsack=10000;
3007                      int extra1 = parameters[whichParam(EXTRA1,numberParameters,parameters)].intValue();
3008                      int extra2 = parameters[whichParam(EXTRA2,numberParameters,parameters)].intValue();
3009                      int logLevel = parameters[log].intValue();
3010                      OsiSolverInterface * solver = expandKnapsack(saveCoinModel,whichColumn,knapsackStart,
3011                                                                   knapsackRow,numberKnapsack,
3012                                                                   storedAmpl,logLevel,extra1,extra2);
3013                      if (solver) {
3014                        clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
3015                        assert (clpSolver);
3016                        lpSolver = clpSolver->getModelPtr();
3017                        babModel->assignSolver(solver);
3018                        testOsiOptions=0;
3019                        // Priorities already done
3020                        free(info.priorities);
3021                        info.priorities=NULL;
3022                      } else {
3023                        numberKnapsack=0;
3024                        delete [] whichColumn;
3025                        delete [] knapsackStart;
3026                        delete [] knapsackRow;
3027                        whichColumn=NULL;
3028                        knapsackStart=NULL;
3029                        knapsackRow=NULL;
3030                      }
3031                    }
3032                  }
3033                }
3034              }
3035#endif
3036              if (useCosts&&testOsiOptions<0) {
3037                int numberColumns = babModel->getNumCols();
3038                int * sort = new int[numberColumns];
3039                double * dsort = new double[numberColumns];
3040                int * priority = new int [numberColumns];
3041                const double * objective = babModel->getObjCoefficients();
3042                const double * lower = babModel->getColLower() ;
3043                const double * upper = babModel->getColUpper() ;
3044                const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
3045                const int * columnLength = matrix->getVectorLengths();
3046                int iColumn;
3047                int n=0;
3048                for (iColumn=0;iColumn<numberColumns;iColumn++) {
3049                  if (babModel->isInteger(iColumn)) {
3050                    sort[n]=n;
3051                    if (useCosts==1)
3052                      dsort[n++]=-fabs(objective[iColumn]);
3053                    else if (useCosts==2)
3054                      dsort[n++]=iColumn;
3055                    else if (useCosts==3)
3056                      dsort[n++]=upper[iColumn]-lower[iColumn];
3057                    else if (useCosts==4)
3058                      dsort[n++]=-(upper[iColumn]-lower[iColumn]);
3059                    else if (useCosts==5)
3060                      dsort[n++]=-columnLength[iColumn];
3061                  }
3062                }
3063                CoinSort_2(dsort,dsort+n,sort);
3064                int level=0;
3065                double last = -1.0e100;
3066                for (int i=0;i<n;i++) {
3067                  int iPut=sort[i];
3068                  if (dsort[i]!=last) {
3069                    level++;
3070                    last=dsort[i];
3071                  }
3072                  priority[iPut]=level;
3073                }
3074                babModel->passInPriorities( priority,false);
3075                delete [] priority;
3076                delete [] sort;
3077                delete [] dsort;
3078              }
3079              // FPump done first as it only works if no solution
3080              CbcHeuristicFPump heuristic4(*babModel);
3081              if (useFpump) {
3082                heuristic4.setMaximumPasses(parameters[whichParam(FPUMPITS,numberParameters,parameters)].intValue());
3083                int pumpTune=parameters[whichParam(FPUMPTUNE,numberParameters,parameters)].intValue();
3084                if (pumpTune>0) {
3085                  /*
3086                    >=10000000 for using obj
3087                    >=1000000 use as accumulate switch
3088                    >=1000 use index+1 as number of large loops
3089                    >=100 use 0.05 objvalue as increment
3090                    >=10 use +0.1 objvalue for cutoff (add)
3091                    1 == fix ints at bounds, 2 fix all integral ints, 3 and continuous at bounds
3092                    4 and static continuous, 5 as 3 but no internal integers
3093                    6 as 3 but all slack basis!
3094                  */
3095                  int logLevel = parameters[log].intValue();
3096                  double value = babModel->solver()->getObjSense()*babModel->solver()->getObjValue();
3097                  int w = pumpTune/10;
3098                  int c = w % 10;
3099                  w /= 10;
3100                  int i = w % 10;
3101                  w /= 10;
3102                  int r = w;
3103                  int accumulate = r/1000;
3104                  r -= 1000*accumulate;
3105                  if (accumulate>=10) {
3106                    int which = accumulate/10;
3107                    accumulate -= 10*which;
3108                    which--;
3109                    // weights and factors
3110                    double weight[]={0.1,0.1,0.5,0.5,1.0,1.0,5.0,5.0};
3111                    double factor[] = {0.1,0.5,0.1,0.5,0.1,0.5,0.1,0.5};
3112                    heuristic4.setInitialWeight(weight[which]);
3113                    heuristic4.setWeightFactor(factor[which]);
3114                  }
3115                  // fake cutoff
3116                  if (logLevel>1)
3117                    printf("Setting ");
3118                  if (c) {
3119                    double cutoff;
3120                    babModel->solver()->getDblParam(OsiDualObjectiveLimit,cutoff);
3121                    cutoff = CoinMin(cutoff,value + 0.1*fabs(value)*c);
3122                    heuristic4.setFakeCutoff(cutoff);
3123                    if (logLevel>1)
3124                      printf("fake cutoff of %g ",cutoff);
3125                  }
3126                  if (i||r) {
3127                    // also set increment
3128                    heuristic4.setAbsoluteIncrement((0.01*i+0.005)*(fabs(value)+1.0e-12));
3129                    heuristic4.setAccumulate(accumulate);
3130                    heuristic4.setMaximumRetries(r+1);
3131                    if (logLevel>1) {
3132                      if (i) 
3133                        printf("increment of %g ",heuristic4.absoluteIncrement());
3134                      if (accumulate) 
3135                        printf("accumulate of %d ",accumulate);
3136                      printf("%d retries ",r+2);
3137                    }
3138                  }
3139                  pumpTune = pumpTune%100;
3140                  if (logLevel>1)
3141                    printf("and setting when to %d\n",pumpTune+10);
3142                  if (pumpTune==6)
3143                    pumpTune =13;
3144                  heuristic4.setWhen(pumpTune+10);
3145                }
3146                babModel->addHeuristic(&heuristic4);
3147              }
3148              if (!miplib) {
3149                CbcRounding heuristic1(*babModel);
3150                if (useRounding)
3151                  babModel->addHeuristic(&heuristic1) ;
3152                CbcHeuristicLocal heuristic2(*babModel);
3153                heuristic2.setSearchType(1);
3154                if (useCombine)
3155                  babModel->addHeuristic(&heuristic2);
3156                CbcHeuristicGreedyCover heuristic3(*babModel);
3157                CbcHeuristicGreedyEquality heuristic3a(*babModel);
3158                if (useGreedy) {
3159                  babModel->addHeuristic(&heuristic3);
3160                  babModel->addHeuristic(&heuristic3a);
3161                }
3162                if (useLocalTree) {
3163                  CbcTreeLocal localTree(babModel,NULL,10,0,0,10000,2000);
3164                  babModel->passInTreeHandler(localTree);
3165                }
3166              }
3167              CbcHeuristicRINS heuristic5(*babModel);
3168              if (useRINS)
3169                babModel->addHeuristic(&heuristic5) ;
3170              if (type==MIPLIB) {
3171                if (babModel->numberStrong()==5&&babModel->numberBeforeTrust()==5) 
3172                  babModel->setNumberBeforeTrust(50);
3173              }
3174              // add cut generators if wanted
3175              int switches[20];
3176              int numberGenerators=0;
3177              int translate[]={-100,-1,-99,-98,1,1,1,1};
3178              if (probingAction) {
3179                if (probingAction==5||probingAction==7)
3180                  probingGen.setRowCuts(-3); // strengthening etc just at root
3181                if (probingAction==6||probingAction==7) {
3182                  // Number of unsatisfied variables to look at
3183                  probingGen.setMaxProbe(1000);
3184                  probingGen.setMaxProbeRoot(1000);
3185                  // How far to follow the consequences
3186                  probingGen.setMaxLook(50);
3187                  probingGen.setMaxLookRoot(50);
3188                }
3189                babModel->addCutGenerator(&probingGen,translate[probingAction],"Probing");
3190                switches[numberGenerators++]=0;
3191              }
3192              if (gomoryAction&&(complicatedInteger!=1||gomoryAction==1)) {
3193                babModel->addCutGenerator(&gomoryGen,translate[gomoryAction],"Gomory");
3194                switches[numberGenerators++]=-1;
3195              }
3196              if (knapsackAction) {
3197                babModel->addCutGenerator(&knapsackGen,translate[knapsackAction],"Knapsack");
3198                switches[numberGenerators++]=0;
3199              }
3200              if (redsplitAction&&!complicatedInteger) {
3201                babModel->addCutGenerator(&redsplitGen,translate[redsplitAction],"Reduce-and-split");
3202                switches[numberGenerators++]=1;
3203              }
3204              if (cliqueAction) {
3205                babModel->addCutGenerator(&cliqueGen,translate[cliqueAction],"Clique");
3206                switches[numberGenerators++]=0;
3207              }
3208              if (mixedAction) {
3209                babModel->addCutGenerator(&mixedGen,translate[mixedAction],"MixedIntegerRounding2");
3210                switches[numberGenerators++]=-1;
3211              }
3212              if (flowAction) {
3213                babModel->addCutGenerator(&flowGen,translate[flowAction],"FlowCover");
3214                switches[numberGenerators++]=1;
3215              }
3216              if (twomirAction&&!complicatedInteger) {
3217                babModel->addCutGenerator(&twomirGen,translate[twomirAction],"TwoMirCuts");
3218                switches[numberGenerators++]=1;
3219              }
3220              if (landpAction) {
3221                babModel->addCutGenerator(&landpGen,translate[landpAction],"LiftAndProject");
3222                switches[numberGenerators++]=1;
3223              }
3224              if (storedCuts) 
3225                babModel->setSpecialOptions(babModel->specialOptions()|64);
3226              // Say we want timings
3227              numberGenerators = babModel->numberCutGenerators();
3228              int iGenerator;
3229              int cutDepth=
3230                parameters[whichParam(CUTDEPTH,numberParameters,parameters)].intValue();
3231              for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
3232                CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
3233                int howOften = generator->howOften();
3234                if (howOften==-98||howOften==-99) 
3235                  generator->setSwitchOffIfLessThan(switches[iGenerator]);
3236                generator->setTiming(true);
3237                if (cutDepth>=0)
3238                  generator->setWhatDepth(cutDepth) ;
3239              }
3240              // Could tune more
3241              if (!miplib) {
3242                babModel->setMinimumDrop(min(5.0e-2,
3243                                             fabs(babModel->getMinimizationObjValue())*1.0e-3+1.0e-4));
3244                if (cutPass==-1234567) {
3245                  if (babModel->getNumCols()<500)
3246                    babModel->setMaximumCutPassesAtRoot(-100); // always do 100 if possible
3247                  else if (babModel->getNumCols()<5000)
3248                    babModel->setMaximumCutPassesAtRoot(100); // use minimum drop
3249                  else
3250                    babModel->setMaximumCutPassesAtRoot(20);
3251                } else {
3252                  babModel->setMaximumCutPassesAtRoot(cutPass);
3253                }
3254                babModel->setMaximumCutPasses(1);
3255              }
3256              // Do more strong branching if small
3257              //if (babModel->getNumCols()<5000)
3258              //babModel->setNumberStrong(20);
3259              // Switch off strong branching if wanted
3260              //if (babModel->getNumCols()>10*babModel->getNumRows())
3261              //babModel->setNumberStrong(0);
3262              if (!noPrinting) {
3263                int iLevel = parameters[log].intValue();
3264                if (iLevel<0) {
3265                  babModel->setPrintingMode(1);
3266                  iLevel = -iLevel;
3267                }
3268                babModel->messageHandler()->setLogLevel(iLevel);
3269                if (babModel->getNumCols()>2000||babModel->getNumRows()>1500||
3270                    babModel->messageHandler()->logLevel()>1)
3271                  babModel->setPrintFrequency(100);
3272              }
3273             
3274              babModel->solver()->setIntParam(OsiMaxNumIterationHotStart,
3275                    parameters[whichParam(MAXHOTITS,numberParameters,parameters)].intValue());
3276              OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3277              // go faster stripes
3278              if (osiclp->getNumRows()<300&&osiclp->getNumCols()<500) {
3279                osiclp->setupForRepeatedUse(2,parameters[slog].intValue());
3280              } else {
3281                osiclp->setupForRepeatedUse(0,parameters[slog].intValue());
3282              }
3283              double increment=babModel->getCutoffIncrement();;
3284              int * changed = NULL;
3285              if (!miplib&&increment==normalIncrement)
3286                changed=analyze( osiclp,numberChanged,increment,false);
3287              if (debugValues) {
3288                if (numberDebugValues==babModel->getNumCols()) {
3289                  // for debug
3290                  babModel->solver()->activateRowCutDebugger(debugValues) ;
3291                } else {
3292                  printf("debug file has incorrect number of columns\n");
3293                }
3294              }
3295              babModel->setCutoffIncrement(CoinMax(babModel->getCutoffIncrement(),increment));
3296              // Turn this off if you get problems
3297              // Used to be automatically set
3298              int mipOptions = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()%10000;
3299              if (mipOptions!=(1))
3300                printf("mip options %d\n",mipOptions);
3301              osiclp->setSpecialOptions(mipOptions);
3302              if (gapRatio < 1.0e100) {
3303                double value = babModel->solver()->getObjValue() ;
3304                double value2 = gapRatio*(1.0e-5+fabs(value)) ;
3305                babModel->setAllowableGap(value2) ;
3306                std::cout << "Continuous " << value
3307                          << ", so allowable gap set to "
3308                          << value2 << std::endl ;
3309              }
3310              // probably faster to use a basis to get integer solutions
3311              babModel->setSpecialOptions(babModel->specialOptions()|2);
3312              currentBranchModel = babModel;
3313              OsiSolverInterface * strengthenedModel=NULL;
3314              if (type==BAB||type==MIPLIB) {
3315                int moreMipOptions = parameters[whichParam(MOREMIPOPTIONS,numberParameters,parameters)].intValue();
3316                if (moreMipOptions>=0) {
3317                  printf("more mip options %d\n",moreMipOptions);
3318                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3319                  if (moreMipOptions==10000) {
3320                    // test memory saving
3321                    moreMipOptions -= 10000;
3322                    ClpSimplex * lpSolver = osiclp->getModelPtr();
3323                    lpSolver->setPersistenceFlag(1);
3324                    // switch off row copy if few rows
3325                    if (lpSolver->numberRows()<150)
3326                      lpSolver->setSpecialOptions(lpSolver->specialOptions()|256);
3327                  }
3328                  if (((moreMipOptions+1)%1000000)!=0)
3329                    babModel->setSearchStrategy(moreMipOptions%1000000);
3330                  // go faster stripes
3331                  if( moreMipOptions >=999999) {
3332                    if (osiclp) {
3333                      int save = osiclp->specialOptions();
3334                      osiclp->setupForRepeatedUse(2,0);
3335                      osiclp->setSpecialOptions(save|osiclp->specialOptions());
3336                    }
3337                  } 
3338                }
3339              }
3340              if (type==BAB) {
3341#ifdef COIN_HAS_ASL
3342                if (usingAmpl) {
3343                  priorities=info.priorities;
3344                  branchDirection=info.branchDirection;
3345                  pseudoDown=info.pseudoDown;
3346                  pseudoUp=info.pseudoUp;
3347                  solutionIn=info.primalSolution;
3348                  prioritiesIn = info.priorities;
3349                  if (info.numberSos&&doSOS) {
3350                    // SOS
3351                    numberSOS = info.numberSos;
3352                    sosStart = info.sosStart;
3353                    sosIndices = info.sosIndices;
3354                    sosType = info.sosType;
3355                    sosReference = info.sosReference;
3356                    sosPriority = info.sosPriority;
3357                  }
3358                }
3359#endif               
3360                const int * originalColumns = preProcess ? process.originalColumns() : NULL;
3361                if (solutionIn&&useSolution) {
3362                  if (preProcess) {
3363                    int numberColumns = babModel->getNumCols();
3364                    // extend arrays in case SOS
3365                    int n = originalColumns[numberColumns-1]+1;
3366                    int nSmaller = CoinMin(n,numberOriginalColumns);
3367                    double * solutionIn2 = new double [n];
3368                    int * prioritiesIn2 = new int[n];
3369                    int i;
3370                    for (i=0;i<nSmaller;i++) {
3371                      solutionIn2[i]=solutionIn[i];
3372                      prioritiesIn2[i]=prioritiesIn[i];
3373                    }
3374                    for (;i<n;i++) {
3375                      solutionIn2[i]=0.0;
3376                      prioritiesIn2[i]=1000000;
3377                    }
3378                    int iLast=-1;
3379                    for (i=0;i<numberColumns;i++) {
3380                      int iColumn = originalColumns[i];
3381                      assert (iColumn>iLast);
3382                      iLast=iColumn;
3383                      solutionIn2[i]=solutionIn2[iColumn];
3384                      if (prioritiesIn)
3385                        prioritiesIn2[i]=prioritiesIn2[iColumn];
3386                    }
3387                    babModel->setHotstartSolution(solutionIn2,prioritiesIn2);
3388                    delete [] solutionIn2;
3389                    delete [] prioritiesIn2;
3390                  } else {
3391                    babModel->setHotstartSolution(solutionIn,prioritiesIn);
3392                  }
3393                }
3394                OsiSolverInterface * testOsiSolver= (testOsiOptions>=0) ? babModel->solver() : NULL;
3395                if (!testOsiSolver) {
3396                  // *************************************************************
3397                  // CbcObjects
3398                  if (preProcess&&process.numberSOS()) {
3399                    int numberSOS = process.numberSOS();
3400                    int numberIntegers = babModel->numberIntegers();
3401                    /* model may not have created objects
3402                       If none then create
3403                    */
3404                    if (!numberIntegers||!babModel->numberObjects()) {
3405                      int type = (pseudoUp) ? 1 : 0;
3406                      babModel->findIntegers(true,type);
3407                      numberIntegers = babModel->numberIntegers();
3408                    }
3409                    OsiObject ** oldObjects = babModel->objects();
3410                    // Do sets and priorities
3411                    OsiObject ** objects = new OsiObject * [numberSOS];
3412                    // set old objects to have low priority
3413                    int numberOldObjects = babModel->numberObjects();
3414                    int numberColumns = babModel->getNumCols();
3415                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
3416                      oldObjects[iObj]->setPriority(numberColumns+1);
3417                      int iColumn = oldObjects[iObj]->columnNumber();
3418                      assert (iColumn>=0);
3419                      if (iColumn>=numberOriginalColumns)
3420                        continue;
3421                      if (originalColumns)
3422                        iColumn = originalColumns[iColumn];
3423                      if (branchDirection) {
3424                        CbcSimpleInteger * obj =
3425                          dynamic_cast <CbcSimpleInteger *>(oldObjects[iObj]) ;
3426                        if (obj) { 
3427                          obj->setPreferredWay(branchDirection[iColumn]);
3428                        } else {
3429                          CbcObject * obj =
3430                            dynamic_cast <CbcObject *>(oldObjects[iObj]) ;
3431                          assert (obj);
3432                          obj->setPreferredWay(branchDirection[iColumn]);
3433                        }
3434                      }
3435                      if (pseudoUp) {
3436                        CbcSimpleIntegerPseudoCost * obj1a =
3437                          dynamic_cast <CbcSimpleIntegerPseudoCost *>(oldObjects[iObj]) ;
3438                        assert (obj1a);
3439                        if (pseudoDown[iColumn]>0.0)
3440                          obj1a->setDownPseudoCost(pseudoDown[iColumn]);
3441                        if (pseudoUp[iColumn]>0.0)
3442                          obj1a->setUpPseudoCost(pseudoUp[iColumn]);
3443                      }
3444                    }
3445                    const int * starts = process.startSOS();
3446                    const int * which = process.whichSOS();
3447                    const int * type = process.typeSOS();
3448                    const double * weight = process.weightSOS();
3449                    int iSOS;
3450                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
3451                      int iStart = starts[iSOS];
3452                      int n=starts[iSOS+1]-iStart;
3453                      objects[iSOS] = new CbcSOS(babModel,n,which+iStart,weight+iStart,
3454                                                 iSOS,type[iSOS]);
3455                      // branch on long sets first
3456                      objects[iSOS]->setPriority(numberColumns-n);
3457                    }
3458                    babModel->addObjects(numberSOS,objects);
3459                    for (iSOS=0;iSOS<numberSOS;iSOS++)
3460                      delete objects[iSOS];
3461                    delete [] objects;
3462                  } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
3463                    // do anyway for priorities etc
3464                    int numberIntegers = babModel->numberIntegers();
3465                    /* model may not have created objects
3466                       If none then create
3467                    */
3468                    if (!numberIntegers||!babModel->numberObjects()) {
3469                      int type = (pseudoUp) ? 1 : 0;
3470                      babModel->findIntegers(true,type);
3471                    }
3472                    if (numberSOS) {
3473                      // Do sets and priorities
3474                      OsiObject ** objects = new OsiObject * [numberSOS];
3475                      int iSOS;
3476                      if (originalColumns) {
3477                        // redo sequence numbers
3478                        int numberColumns = babModel->getNumCols();
3479                        int nOld = originalColumns[numberColumns-1]+1;
3480                        int * back = new int[nOld];
3481                        int i;
3482                        for (i=0;i<nOld;i++)
3483                          back[i]=-1;
3484                        for (i=0;i<numberColumns;i++)
3485                          back[originalColumns[i]]=i;
3486                        // Really need better checks
3487                        int nMissing=0;
3488                        int n=sosStart[numberSOS];
3489                        for (i=0;i<n;i++) {
3490                          int iColumn = sosIndices[i];
3491                          int jColumn = back[iColumn];
3492                          if (jColumn>=0) 
3493                            sosIndices[i] = jColumn;
3494                          else 
3495                            nMissing++;
3496                        }
3497                        delete [] back;
3498                        if (nMissing)
3499                          printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
3500                      }
3501                      for (iSOS =0;iSOS<numberSOS;iSOS++) {
3502                        int iStart = sosStart[iSOS];
3503                        int n=sosStart[iSOS+1]-iStart;
3504                        objects[iSOS] = new CbcSOS(babModel,n,sosIndices+iStart,sosReference+iStart,
3505                                                   iSOS,sosType[iSOS]);
3506                        if (sosPriority)
3507                          objects[iSOS]->setPriority(sosPriority[iSOS]);
3508                        else if (!prioritiesIn)
3509                          objects[iSOS]->setPriority(10);  // rather than 1000
3510                      }
3511                      // delete any existing SOS objects
3512                      int numberObjects=babModel->numberObjects();
3513                      OsiObject ** oldObjects=babModel->objects();
3514                      int nNew=0;
3515                      for (int i=0;i<numberObjects;i++) {
3516                        OsiObject * objThis = oldObjects[i];
3517                        CbcSOS * obj1 =
3518                          dynamic_cast <CbcSOS *>(objThis) ;
3519                        OsiSOS * obj2 =
3520                          dynamic_cast <OsiSOS *>(objThis) ;
3521                        if (!obj1&&!obj2) {
3522                          oldObjects[nNew++]=objThis;
3523                        } else {
3524                          delete objThis;
3525                        }
3526                      }
3527                      babModel->setNumberObjects(nNew);
3528                      babModel->addObjects(numberSOS,objects);
3529                      for (iSOS=0;iSOS<numberSOS;iSOS++)
3530                        delete objects[iSOS];
3531                      delete [] objects;
3532                    }
3533                  }
3534                  OsiObject ** objects = babModel->objects();
3535                  int numberObjects = babModel->numberObjects();
3536                  for (int iObj = 0;iObj<numberObjects;iObj++) {
3537                    // skip sos
3538                    CbcSOS * objSOS =
3539                      dynamic_cast <CbcSOS *>(objects[iObj]) ;
3540                    if (objSOS)
3541                      continue;
3542                    int iColumn = objects[iObj]->columnNumber();
3543                    assert (iColumn>=0);
3544                    if (originalColumns)
3545                      iColumn = originalColumns[iColumn];
3546                    if (branchDirection) {
3547                      CbcSimpleInteger * obj =
3548                        dynamic_cast <CbcSimpleInteger *>(objects[iObj]) ;
3549                      if (obj) { 
3550                        obj->setPreferredWay(branchDirection[iColumn]);
3551                      } else {
3552                        CbcObject * obj =
3553                          dynamic_cast <CbcObject *>(objects[iObj]) ;
3554                        assert (obj);
3555                        obj->setPreferredWay(branchDirection[iColumn]);
3556                      }
3557                    }
3558                    if (priorities) {
3559                      int iPriority = priorities[iColumn];
3560                      if (iPriority>0)
3561                        objects[iObj]->setPriority(iPriority);
3562                    }
3563                    if (pseudoUp&&pseudoUp[iColumn]) {
3564                      CbcSimpleIntegerPseudoCost * obj1a =
3565                        dynamic_cast <CbcSimpleIntegerPseudoCost *>(objects[iObj]) ;
3566                      assert (obj1a);
3567                      if (pseudoDown[iColumn]>0.0)
3568                        obj1a->setDownPseudoCost(pseudoDown[iColumn]);
3569                      if (pseudoUp[iColumn]>0.0)
3570                        obj1a->setUpPseudoCost(pseudoUp[iColumn]);
3571                    }
3572                  }
3573                  // *************************************************************
3574                } else {
3575                  // *************************************************************
3576                  // OsiObjects
3577                  // Find if none
3578                  int numberIntegers = testOsiSolver->getNumIntegers();
3579                  /* model may not have created objects
3580                     If none then create
3581                  */
3582                  if (!numberIntegers||!testOsiSolver->numberObjects()) {
3583                    //int type = (pseudoUp) ? 1 : 0;
3584                    testOsiSolver->findIntegers(false);
3585                    numberIntegers = testOsiSolver->getNumIntegers();
3586                  }
3587                  if (preProcess&&process.numberSOS()) {
3588                    int numberSOS = process.numberSOS();
3589                    OsiObject ** oldObjects = testOsiSolver->objects();
3590                    // Do sets and priorities
3591                    OsiObject ** objects = new OsiObject * [numberSOS];
3592                    // set old objects to have low priority
3593                    int numberOldObjects = testOsiSolver->numberObjects();
3594                    int numberColumns = testOsiSolver->getNumCols();
3595                    for (int iObj = 0;iObj<numberOldObjects;iObj++) {
3596                      oldObjects[iObj]->setPriority(numberColumns+1);
3597                      int iColumn = oldObjects[iObj]->columnNumber();
3598                      assert (iColumn>=0);
3599                      if (iColumn>=numberOriginalColumns)
3600                        continue;
3601                      if (originalColumns)
3602                        iColumn = originalColumns[iColumn];
3603                      if (branchDirection) {
3604                        OsiSimpleInteger * obj =
3605                          dynamic_cast <OsiSimpleInteger *>(oldObjects[iObj]) ;
3606                        if (obj) { 
3607                          obj->setPreferredWay(branchDirection[iColumn]);
3608                        } else {
3609                          OsiObject2 * obj =
3610                            dynamic_cast <OsiObject2 *>(oldObjects[iObj]) ;
3611                          if (obj)
3612                            obj->setPreferredWay(branchDirection[iColumn]);
3613                        }
3614                      }
3615                      if (pseudoUp) {
3616                        abort();
3617                      }
3618                    }
3619                    const int * starts = process.startSOS();
3620                    const int * which = process.whichSOS();
3621                    const int * type = process.typeSOS();
3622                    const double * weight = process.weightSOS();
3623                    int iSOS;
3624                    for (iSOS =0;iSOS<numberSOS;iSOS++) {
3625                      int iStart = starts[iSOS];
3626                      int n=starts[iSOS+1]-iStart;
3627                      objects[iSOS] = new OsiSOS(testOsiSolver,n,which+iStart,weight+iStart,
3628                                                 type[iSOS]);
3629                      // branch on long sets first
3630                      objects[iSOS]->setPriority(numberColumns-n);
3631                    }
3632                    testOsiSolver->addObjects(numberSOS,objects);
3633                    for (iSOS=0;iSOS<numberSOS;iSOS++)
3634                      delete objects[iSOS];
3635                    delete [] objects;
3636                  } else if (priorities||branchDirection||pseudoDown||pseudoUp||numberSOS) {
3637                    if (numberSOS) {
3638                      // Do sets and priorities
3639                      OsiObject ** objects = new OsiObject * [numberSOS];
3640                      int iSOS;
3641                      if (originalColumns) {
3642                        // redo sequence numbers
3643                        int numberColumns = testOsiSolver->getNumCols();
3644                        int nOld = originalColumns[numberColumns-1]+1;
3645                        int * back = new int[nOld];
3646                        int i;
3647                        for (i=0;i<nOld;i++)
3648                          back[i]=-1;
3649                        for (i=0;i<numberColumns;i++)
3650                          back[originalColumns[i]]=i;
3651                        // Really need better checks
3652                        int nMissing=0;
3653                        int n=sosStart[numberSOS];
3654                        for (i=0;i<n;i++) {
3655                          int iColumn = sosIndices[i];
3656                          int jColumn = back[iColumn];
3657                          if (jColumn>=0) 
3658                            sosIndices[i] = jColumn;
3659                          else 
3660                            nMissing++;
3661                        }
3662                        delete [] back;
3663                        if (nMissing)
3664                          printf("%d SOS variables vanished due to pre processing? - check validity?\n",nMissing);
3665                      }
3666                      for (iSOS =0;iSOS<numberSOS;iSOS++) {
3667                        int iStart = sosStart[iSOS];
3668                        int n=sosStart[iSOS+1]-iStart;
3669                        objects[iSOS] = new OsiSOS(testOsiSolver,n,sosIndices+iStart,sosReference+iStart,
3670                                                   sosType[iSOS]);
3671                        if (sosPriority)
3672                          objects[iSOS]->setPriority(sosPriority[iSOS]);
3673                        else if (!prioritiesIn)
3674                          objects[iSOS]->setPriority(10);  // rather than 1000
3675                      }
3676                      // delete any existing SOS objects
3677                      int numberObjects=testOsiSolver->numberObjects();
3678                      OsiObject ** oldObjects=testOsiSolver->objects();
3679                      int nNew=0;
3680                      for (int i=0;i<numberObjects;i++) {
3681                        OsiObject * objThis = oldObjects[i];
3682                        OsiSOS * obj1 =
3683                          dynamic_cast <OsiSOS *>(objThis) ;
3684                        OsiSOS * obj2 =
3685                          dynamic_cast <OsiSOS *>(objThis) ;
3686                        if (!obj1&&!obj2) {
3687                          oldObjects[nNew++]=objThis;
3688                        } else {
3689                          delete objThis;
3690                        }
3691                      }
3692                      testOsiSolver->setNumberObjects(nNew);
3693                      testOsiSolver->addObjects(numberSOS,objects);
3694                      for (iSOS=0;iSOS<numberSOS;iSOS++)
3695                        delete objects[iSOS];
3696                      delete [] objects;
3697                    }
3698                  }
3699                  OsiObject ** objects = testOsiSolver->objects();
3700                  int numberObjects = testOsiSolver->numberObjects();
3701                  int logLevel = parameters[log].intValue();
3702                  for (int iObj = 0;iObj<numberObjects;iObj++) {
3703                    // skip sos
3704                    OsiSOS * objSOS =
3705                      dynamic_cast <OsiSOS *>(objects[iObj]) ;
3706                    if (objSOS) {
3707                      if (logLevel>2)
3708                        printf("Set %d is SOS - priority %d\n",iObj,objSOS->priority());
3709                      continue;
3710                    }
3711                    int iColumn = objects[iObj]->columnNumber();
3712                    if (iColumn>=0) {
3713                      if (originalColumns)
3714                        iColumn = originalColumns[iColumn];
3715                      if (branchDirection) {
3716                        OsiSimpleInteger * obj =
3717                          dynamic_cast <OsiSimpleInteger *>(objects[iObj]) ;
3718                        if (obj) { 
3719                          obj->setPreferredWay(branchDirection[iColumn]);
3720                        } else {
3721                          OsiObject2 * obj =
3722                            dynamic_cast <OsiObject2 *>(objects[iObj]) ;
3723                          if (obj)
3724                            obj->setPreferredWay(branchDirection[iColumn]);
3725                        }
3726                      }
3727                    }
3728                    if (priorities) {
3729                      int iPriority = priorities[iColumn];
3730                      if (iPriority>0)
3731                        objects[iObj]->setPriority(iPriority);
3732                    }
3733                    if (logLevel>2)
3734                      printf("Obj %d is int? - priority %d\n",iObj,objects[iObj]->priority());
3735                    if (pseudoUp&&pseudoUp[iColumn]) {
3736                      abort();
3737                    }
3738                  }
3739                  // *************************************************************
3740                }
3741                int statistics = (printOptions>0) ? printOptions: 0;
3742#ifdef COIN_HAS_ASL
3743                if (!usingAmpl) {
3744#endif
3745                  free(priorities);
3746                  priorities=NULL;
3747                  free(branchDirection);
3748                  branchDirection=NULL;
3749                  free(pseudoDown);
3750                  pseudoDown=NULL;
3751                  free(pseudoUp);
3752                  pseudoUp=NULL;
3753                  free(solutionIn);
3754                  solutionIn=NULL;
3755                  free(prioritiesIn);
3756                  prioritiesIn=NULL;
3757                  free(sosStart);
3758                  sosStart=NULL;
3759                  free(sosIndices);
3760                  sosIndices=NULL;
3761                  free(sosType);
3762                  sosType=NULL;
3763                  free(sosReference);
3764                  sosReference=NULL;
3765                  free(cut);
3766                  cut=NULL;
3767                  free(sosPriority);
3768                  sosPriority=NULL;
3769#ifdef COIN_HAS_ASL
3770                }
3771#endif               
3772                if (nodeStrategy) {
3773                  // change default
3774                  if (nodeStrategy>1) {
3775                    // up or down
3776                    int way = ((nodeStrategy%1)==1) ? -1 : +1;
3777                    babModel->setPreferredWay(way);
3778#if 0
3779                    OsiObject ** objects = babModel->objects();
3780                    int numberObjects = babModel->numberObjects();
3781                    for (int iObj = 0;iObj<numberObjects;iObj++) {
3782                      CbcObject * obj =
3783                        dynamic_cast <CbcObject *>(objects[iObj]) ;
3784                      assert (obj);
3785                      obj->setPreferredWay(way);
3786                    }
3787#endif
3788                  }
3789                  if (nodeStrategy==1||nodeStrategy>3) {
3790                    // depth
3791                    CbcCompareDefault compare;
3792                    compare.setWeight(-3.0);
3793                    babModel->setNodeComparison(compare);
3794                  }
3795                }
3796                if (cppValue>=0) {
3797                  int prepro = useStrategy ? -1 : preProcess;
3798                  // generate code
3799                  FILE * fp = fopen("user_driver.cpp","w");
3800                  if (fp) {
3801                    // generate enough to do BAB
3802                    babModel->generateCpp(fp,1);
3803                    OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3804                    // Make general so do factorization
3805                    int factor = osiclp->getModelPtr()->factorizationFrequency();
3806                    osiclp->getModelPtr()->setFactorizationFrequency(200);
3807                    osiclp->generateCpp(fp);
3808                    osiclp->getModelPtr()->setFactorizationFrequency(factor);
3809                    //solveOptions.generateCpp(fp);
3810                    fclose(fp);
3811                    // now call generate code
3812                    generateCode(babModel,"user_driver.cpp",cppValue,prepro);
3813                  } else {
3814                    std::cout<<"Unable to open file user_driver.cpp"<<std::endl;
3815                  }
3816                }
3817                if (!babModel->numberStrong())
3818                  babModel->setNumberBeforeTrust(0);
3819                if (useStrategy) {
3820                  CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
3821                  strategy.setupPreProcessing(1);
3822                  babModel->setStrategy(strategy);
3823                }
3824                if (testOsiOptions>=0) {
3825                  printf("Testing OsiObject options %d\n",testOsiOptions);
3826                  if (!numberSOS) {
3827                    babModel->solver()->findIntegersAndSOS(false);
3828#ifdef COIN_HAS_LINK
3829                    // If linked then pass in model
3830                    OsiSolverLink * solver3 = dynamic_cast<OsiSolverLink *> (babModel->solver());
3831                    if (solver3) {
3832                      int options = parameters[whichParam(MIPOPTIONS,numberParameters,parameters)].intValue()/10000;
3833                      CglStored stored;
3834                      if (options) {
3835                        printf("nlp options %d\n",options);
3836                        /*
3837                          1 - force mini branch and bound
3838                          2 - set priorities high on continuous
3839                          4 - try adding OA cuts
3840                          8 - try doing quadratic linearization
3841                          16 - try expanding knapsacks
3842                          32 - OA cuts strictly concave
3843                        */
3844                        if ((options&2)) {
3845                          solver3->setBiLinearPriorities(10,tightenFactor > 0.0 ? tightenFactor : 1.0);
3846                        } else if (tightenFactor>0.0) {
3847                          // set grid size for all continuous bi-linear
3848                          solver3->setMeshSizes(tightenFactor);
3849                        }
3850                        if ((options&4)) {
3851                          solver3->setSpecialOptions2(solver3->specialOptions2()|(8+4));
3852                          // say convex
3853                          solver3->sayConvex((options&32)==0);
3854                        }
3855                        if ((options&1)!=0&&slpValue>0)
3856                          solver3->setFixedPriority(slpValue);
3857                        double cutoff=COIN_DBL_MAX;
3858                        if ((options&8))
3859                          cutoff=solver3->linearizedBAB(&stored);
3860                        if (cutoff<babModel->getCutoff())
3861                          babModel->setCutoff(cutoff);
3862                      }
3863                      solver3->setCbcModel(babModel);
3864                      babModel->addCutGenerator(&stored,1,"Stored");
3865                      CglTemporary temp;
3866                      babModel->addCutGenerator(&temp,1,"OnceOnly");
3867                      //choose.setNumberBeforeTrusted(2000);
3868                      //choose.setNumberStrong(20);
3869                    }
3870#endif
3871                  } else {
3872                    // move across
3873                    babModel->deleteObjects(false);
3874                    //babModel->addObjects(babModel->solver()->numberObjects(),babModel->solver()->objects());
3875                  }
3876                  CbcBranchDefaultDecision decision;
3877                  if (babModel->numberStrong()) {
3878                    OsiChooseStrong choose(babModel->solver());
3879                    choose.setNumberBeforeTrusted(babModel->numberBeforeTrust());
3880                    choose.setNumberStrong(babModel->numberStrong());
3881                    choose.setShadowPriceMode(testOsiOptions);
3882                    decision.setChooseMethod(choose);
3883                  } else {
3884                    OsiChooseVariable choose(babModel->solver());
3885                    decision.setChooseMethod(choose);
3886                  }
3887                  babModel->setBranchingMethod(decision);
3888                  if (useCosts&&testOsiOptions>=0) {
3889                    int numberColumns = babModel->getNumCols();
3890                    int * sort = new int[numberColumns];
3891                    double * dsort = new double[numberColumns];
3892                    int * priority = new int [numberColumns];
3893                    const double * objective = babModel->getObjCoefficients();
3894                    const double * lower = babModel->getColLower() ;
3895                    const double * upper = babModel->getColUpper() ;
3896                    const CoinPackedMatrix * matrix = babModel->solver()->getMatrixByCol();
3897                    const int * columnLength = matrix->getVectorLengths();
3898                    int iColumn;
3899                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
3900                      sort[iColumn]=iColumn;
3901                      if (useCosts==1)
3902                        dsort[iColumn]=-fabs(objective[iColumn]);
3903                      else if (useCosts==2)
3904                        dsort[iColumn]=iColumn;
3905                      else if (useCosts==3)
3906                        dsort[iColumn]=upper[iColumn]-lower[iColumn];
3907                      else if (useCosts==4)
3908                        dsort[iColumn]=-(upper[iColumn]-lower[iColumn]);
3909                      else if (useCosts==5)
3910                        dsort[iColumn]=-columnLength[iColumn];
3911                    }
3912                    CoinSort_2(dsort,dsort+numberColumns,sort);
3913                    int level=0;
3914                    double last = -1.0e100;
3915                    for (int i=0;i<numberColumns;i++) {
3916                      int iPut=sort[i];
3917                      if (dsort[i]!=last) {
3918                        level++;
3919                        last=dsort[i];
3920                      }
3921                      priority[iPut]=level;
3922                    }
3923                    OsiObject ** objects = babModel->objects();
3924                    int numberObjects = babModel->numberObjects();
3925                    for (int iObj = 0;iObj<numberObjects;iObj++) {
3926                      OsiObject * obj = objects[iObj] ;
3927                      int iColumn = obj->columnNumber();
3928                      if (iColumn>=0)
3929                        obj->setPriority(priority[iColumn]);
3930                    }
3931                    delete [] priority;
3932                    delete [] sort;
3933                    delete [] dsort;
3934                  }
3935                }
3936                checkSOS(babModel, babModel->solver());
3937                if (doSprint>0) {
3938                  // Sprint for primal solves
3939                  ClpSolve::SolveType method = ClpSolve::usePrimalorSprint;
3940                  ClpSolve::PresolveType presolveType = ClpSolve::presolveOff;
3941                  int numberPasses = 5;
3942                  int options[] = {0,3,0,0,0,0};
3943                  int extraInfo[] = {-1,20,-1,-1,-1,-1};
3944                  extraInfo[1]=doSprint;
3945                  int independentOptions[] = {0,0,3};
3946                  ClpSolve clpSolve(method,presolveType,numberPasses,
3947                                    options,extraInfo,independentOptions);
3948                  // say use in OsiClp
3949                  clpSolve.setSpecialOption(6,1);
3950                  OsiClpSolverInterface * osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
3951                  osiclp->setSolveOptions(clpSolve);
3952                  osiclp->setHintParam(OsiDoDualInResolve,false);
3953                  // switch off row copy
3954                  osiclp->getModelPtr()->setSpecialOptions(osiclp->getModelPtr()->specialOptions()|256);
3955                  osiclp->getModelPtr()->setInfeasibilityCost(1.0e11);
3956                }
3957#ifdef COIN_HAS_LINK
3958                if (storedAmpl.sizeRowCuts()) {
3959                  if (preProcess) {
3960                    const int * originalColumns = process.originalColumns();
3961                    int numberColumns = babModel->getNumCols();
3962                    int * newColumn = new int[numberOriginalColumns];
3963                    int i;
3964                    for (i=0;i<numberOriginalColumns;i++) 
3965                      newColumn[i]=-1;
3966                    for (i=0;i<numberColumns;i++) {
3967                      int iColumn = originalColumns[i];
3968                      newColumn[iColumn]=i;
3969                    }
3970                    int * buildColumn = new int[numberColumns];
3971                    // Build up valid cuts
3972                    int nBad=0;
3973                    int nCuts = storedAmpl.sizeRowCuts();
3974                    CglStored newCuts;
3975                    for (i=0;i<nCuts;i++) {
3976                      const OsiRowCut * cut = storedAmpl.rowCutPointer(i);
3977                      double lb = cut->lb();
3978                      double ub = cut->ub();
3979                      int n=cut->row().getNumElements();
3980                      const int * column = cut->row().getIndices();
3981                      const double * element = cut->row().getElements();
3982                      bool bad=false;
3983                      for (int i=0;i<n;i++) {
3984                        int iColumn = column[i];
3985                        iColumn = newColumn[iColumn];
3986                        if (iColumn>=0) {
3987                          buildColumn[i]=iColumn;
3988                        } else {
3989                          bad=true;
3990                          break;
3991                        }
3992                      }
3993                      if (!bad) {
3994                        newCuts.addCut(lb,ub,n,buildColumn,element);
3995                      } else {
3996                        nBad++;
3997                      }
3998                    }
3999                    storedAmpl=newCuts;
4000                    if (nBad)
4001                      printf("%d cuts dropped\n",nBad);
4002                    delete [] newColumn;
4003                    delete [] buildColumn;
4004                  }
4005                  if (storedAmpl.sizeRowCuts()) {
4006                    //babModel->addCutGenerator(&storedAmpl,1,"AmplStored");
4007                    int numberRowCuts = storedAmpl.sizeRowCuts();
4008                    for (int i=0;i<numberRowCuts;i++) {
4009                      const OsiRowCut * rowCutPointer = storedAmpl.rowCutPointer(i);
4010                      babModel->makeGlobalCut(rowCutPointer);
4011                    }
4012                  }
4013                }
4014#endif
4015#ifdef CLP_MALLOC_STATISTICS
4016                malloc_stats();
4017                malloc_stats2();
4018#endif
4019                if (outputFormat==5) {
4020                  osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4021                  lpSolver = osiclp->getModelPtr();
4022                  lpSolver->setPersistenceFlag(1);
4023                }
4024                babModel->branchAndBound(statistics);
4025#ifdef CLP_MALLOC_STATISTICS
4026                malloc_stats();
4027                malloc_stats2();
4028#endif
4029                checkSOS(babModel, babModel->solver());
4030              } else if (type==MIPLIB) {
4031                CbcStrategyDefault strategy(true,babModel->numberStrong(),babModel->numberBeforeTrust());
4032                // Set up pre-processing
4033                int translate2[]={9999,1,1,3,2,4,5};
4034                if (preProcess)
4035                  strategy.setupPreProcessing(translate2[preProcess]);
4036                babModel->setStrategy(strategy);
4037                if (outputFormat==5) {
4038                  osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4039                  lpSolver = osiclp->getModelPtr();
4040                  lpSolver->setPersistenceFlag(1);
4041                }
4042                if (testOsiOptions>=0) {
4043                  printf("Testing OsiObject options %d\n",testOsiOptions);
4044                  CbcBranchDefaultDecision decision;
4045                  OsiChooseStrong choose(babModel->solver());
4046                  choose.setNumberBeforeTrusted(babModel->numberBeforeTrust());
4047                  choose.setNumberStrong(babModel->numberStrong());
4048                  choose.setShadowPriceMode(testOsiOptions);
4049                  //babModel->deleteObjects(false);
4050                  decision.setChooseMethod(choose);
4051                  babModel->setBranchingMethod(decision);
4052                }
4053                CbcClpUnitTest(*babModel);
4054                goodModel=false;
4055                break;
4056              } else {
4057                strengthenedModel = babModel->strengthenedModel();
4058              }
4059              currentBranchModel = NULL;
4060              osiclp = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4061              if (debugFile=="createAfterPre"&&babModel->bestSolution()) {
4062                lpSolver = osiclp->getModelPtr();
4063                //move best solution (should be there -- but ..)
4064                int n = lpSolver->getNumCols();
4065                memcpy(lpSolver->primalColumnSolution(),babModel->bestSolution(),n*sizeof(double));
4066                saveSolution(osiclp->getModelPtr(),"debug.file");
4067              }
4068              if (!noPrinting) {
4069                // Print more statistics
4070                std::cout<<"Cuts at root node changed objective from "<<babModel->getContinuousObjective()
4071                         <<" to "<<babModel->rootObjectiveAfterCuts()<<std::endl;
4072               
4073                numberGenerators = babModel->numberCutGenerators();
4074                for (iGenerator=0;iGenerator<numberGenerators;iGenerator++) {
4075                  CbcCutGenerator * generator = babModel->cutGenerator(iGenerator);
4076                  std::cout<<generator->cutGeneratorName()<<" was tried "
4077                           <<generator->numberTimesEntered()<<" times and created "
4078                           <<generator->numberCutsInTotal()<<" cuts of which "
4079                           <<generator->numberCutsActive()<<" were active after adding rounds of cuts";
4080                  if (generator->timing())
4081                    std::cout<<" ( "<<generator->timeInCutGenerator()<<" seconds)"<<std::endl;
4082                  else
4083                    std::cout<<std::endl;
4084                }
4085              }
4086              time2 = CoinCpuTime();
4087              totalTime += time2-time1;
4088              // For best solution
4089              double * bestSolution = NULL;
4090              if (babModel->getMinimizationObjValue()<1.0e50&&type==BAB) {
4091                // post process
4092                if (preProcess) {
4093                  int n = saveSolver->getNumCols();
4094                  bestSolution = new double [n];
4095                  process.postProcess(*babModel->solver());
4096                  // Solution now back in saveSolver
4097                  babModel->assignSolver(saveSolver);
4098                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
4099                } else {
4100                  int n = babModel->solver()->getNumCols();
4101                  bestSolution = new double [n];
4102                  memcpy(bestSolution,babModel->solver()->getColSolution(),n*sizeof(double));
4103                }
4104                checkSOS(babModel, babModel->solver());
4105              }
4106              if (type==STRENGTHEN&&strengthenedModel)
4107                clpSolver = dynamic_cast< OsiClpSolverInterface*> (strengthenedModel);
4108#ifdef COIN_HAS_ASL
4109              else if (usingAmpl) 
4110                clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4111#endif
4112              lpSolver = clpSolver->getModelPtr();
4113              if (numberChanged) {
4114                for (int i=0;i<numberChanged;i++) {
4115                  int iColumn=changed[i];
4116                  clpSolver->setContinuous(iColumn);
4117                }
4118                delete [] changed;
4119              }
4120              if (type==BAB) {
4121                //move best solution (should be there -- but ..)
4122                int n = lpSolver->getNumCols();
4123                if (bestSolution) {
4124                  memcpy(lpSolver->primalColumnSolution(),bestSolution,n*sizeof(double));
4125                  // now see what that does to row solution
4126                  int numberRows=lpSolver->numberRows();
4127                  double * rowSolution = lpSolver->primalRowSolution();
4128                  memset (rowSolution,0,numberRows*sizeof(double));
4129                  lpSolver->clpMatrix()->times(1.0,bestSolution,rowSolution);
4130                }
4131                if (debugFile=="create"&&bestSolution) {
4132                  saveSolution(lpSolver,"debug.file");
4133                }
4134                delete [] bestSolution;
4135                std::string statusName[]={"Finished","Stopped on ","Difficulties",
4136                                          "","","User ctrl-c"};
4137                std::string minor[]={"","","gap","nodes","time","","solutions","user ctrl-c"};
4138                int iStat = babModel->status();
4139                int iStat2 = babModel->secondaryStatus();
4140                if (!noPrinting)
4141                  std::cout<<"Result - "<<statusName[iStat]<<minor[iStat2]
4142                           <<" objective "<<babModel->getObjValue()<<
4143                    " after "<<babModel->getNodeCount()<<" nodes and "
4144                           <<babModel->getIterationCount()<<
4145                    " iterations - took "<<time2-time1<<" seconds"<<std::endl;
4146#ifdef COIN_HAS_ASL
4147                if (usingAmpl) {
4148                  clpSolver = dynamic_cast< OsiClpSolverInterface*> (babModel->solver());
4149                  lpSolver = clpSolver->getModelPtr();
4150                  double value = babModel->getObjValue()*lpSolver->getObjSense();
4151                  char buf[300];
4152                  int pos=0;
4153                  if (iStat==0) {
4154                    if (babModel->getObjValue()<1.0e40) {
4155                      pos += sprintf(buf+pos,"optimal," );
4156                    } else {
4157                      // infeasible
4158                      iStat=1;
4159                      pos += sprintf(buf+pos,"infeasible,");
4160                    }
4161                  } else if (iStat==1) {
4162                    if (iStat2!=6)
4163                      iStat=3;
4164                    else
4165                      iStat=4;
4166                    pos += sprintf(buf+pos,"stopped on %s,",minor[iStat2].c_str());
4167                  } else if (iStat==2) {
4168                    iStat = 7;
4169                    pos += sprintf(buf+pos,"stopped on difficulties,");
4170                  } else if (iStat==5) {
4171                    iStat = 3;
4172                    pos += sprintf(buf+pos,"stopped on ctrl-c,");
4173                  } else {
4174                    pos += sprintf(buf+pos,"status unknown,");
4175                    iStat=6;
4176                  }
4177                  info.problemStatus=iStat;
4178                  info.objValue = value;
4179                  if (babModel->getObjValue()<1.0e40) {
4180                    int precision = ampl_obj_prec();
4181                    if (precision>0)
4182                      pos += sprintf(buf+pos," objective %.*g",precision,
4183                                     value);
4184                    else
4185                      pos += sprintf(buf+pos," objective %g",value);
4186                  }
4187                  sprintf(buf+pos,"\n%d nodes, %d iterations, %g seconds",
4188                          babModel->getNodeCount(),
4189                          babModel->getIterationCount(),
4190                          totalTime);
4191                  if (bestSolution) {
4192                    free(info.primalSolution);
4193                    if (!numberKnapsack) {
4194                      info.primalSolution = (double *) malloc(n*sizeof(double));
4195                      CoinCopyN(lpSolver->primalColumnSolution(),n,info.primalSolution);
4196                      int numberRows = lpSolver->numberRows();
4197                      free(info.dualSolution);
4198                      info.dualSolution = (double *) malloc(numberRows*sizeof(double));
4199                      CoinCopyN(lpSolver->dualRowSolution(),numberRows,info.dualSolution);
4200                    } else {
4201                      // expanded knapsack
4202                      info.dualSolution=NULL;
4203                      int numberColumns = saveCoinModel.numberColumns();
4204                      info.primalSolution = (double *) malloc(numberColumns*sizeof(double));
4205                      // Fills in original solution (coinModel length)
4206                      afterKnapsack(saveCoinModel,  whichColumn,  knapsackStart, 
4207                                    knapsackRow,  numberKnapsack,
4208                                    lpSolver->primalColumnSolution(), info.primalSolution,1);
4209                    }
4210                  } else {
4211                    info.primalSolution=NULL;
4212                    info.dualSolution=NULL;
4213                  }
4214                  // put buffer into info
4215                  strcpy(info.buffer,buf);
4216                }
4217#endif
4218              } else {
4219                std::cout<<"Model strengthened - now has "<<clpSolver->getNumRows()
4220                         <<" rows"<<std::endl;
4221              }
4222              time1 = time2;
4223              delete babModel;
4224              babModel=NULL;
4225            } else {
4226              std::cout << "** Current model not valid" << std::endl ; 
4227            }
4228            break ;
4229          case IMPORT:
4230            {
4231#ifdef COIN_HAS_ASL
4232              if (!usingAmpl) {
4233#endif
4234                free(priorities);
4235                priorities=NULL;
4236                free(branchDirection);
4237                branchDirection=NULL;
4238                free(pseudoDown);
4239                pseudoDown=NULL;
4240                free(pseudoUp);
4241                pseudoUp=NULL;
4242                free(solutionIn);
4243                solutionIn=NULL;
4244                free(prioritiesIn);
4245                prioritiesIn=NULL;
4246                free(sosStart);
4247                sosStart=NULL;
4248                free(sosIndices);
4249                sosIndices=NULL;
4250                free(sosType);
4251                sosType=NULL;
4252                free(sosReference);
4253                sosReference=NULL;
4254                free(cut);
4255                cut=NULL;
4256                free(sosPriority);
4257                sosPriority=NULL;
4258#ifdef COIN_HAS_ASL
4259              }
4260#endif               
4261              delete babModel;
4262              babModel=NULL;
4263              // get next field
4264              field = CoinReadGetString(argc,argv);
4265              if (field=="$") {
4266                field = parameters[iParam].stringValue();
4267              } else if (field=="EOL") {
4268                parameters[iParam].printString();
4269                break;
4270              } else {
4271                parameters[iParam].setStringValue(field);
4272              }
4273              std::string fileName;
4274              bool canOpen=false;
4275              if (field=="-") {
4276                // stdin
4277                canOpen=true;
4278                fileName = "-";
4279              } else {
4280                bool absolutePath;
4281                if (dirsep=='/') {
4282                  // non Windows (or cygwin)
4283                  absolutePath=(field[0]=='/');
4284                } else {
4285                  //Windows (non cycgwin)
4286                  absolutePath=(field[0]=='\\');
4287                  // but allow for :
4288                  if (strchr(field.c_str(),':'))
4289                    absolutePath=true;
4290                }
4291                if (absolutePath) {
4292                  fileName = field;
4293                } else if (field[0]=='~') {
4294                  char * environVar = getenv("HOME");
4295                  if (environVar) {
4296                    std::string home(environVar);
4297                    field=field.erase(0,1);
4298                    fileName = home+field;
4299                  } else {
4300                    fileName=field;
4301                  }
4302                } else {
4303                  fileName = directory+field;
4304                }
4305                FILE *fp=fopen(fileName.c_str(),"r");
4306                if (fp) {
4307                  // can open - lets go for it
4308                  fclose(fp);
4309                  canOpen=true;
4310                } else {
4311                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4312                }
4313              }
4314              if (canOpen) {
4315                int status =clpSolver->readMps(fileName.c_str(),
4316                                                   keepImportNames!=0,
4317                                                   allowImportErrors!=0);
4318                if (!status||(status>0&&allowImportErrors)) {
4319                  if (keepImportNames) {
4320                    lengthName = lpSolver->lengthNames();
4321                    rowNames = *(lpSolver->rowNames());
4322                    columnNames = *(lpSolver->columnNames());
4323                  } else {
4324                    lengthName=0;
4325                  }
4326                  goodModel=true;
4327                  //Set integers in clpsolver
4328                  const char * info = lpSolver->integerInformation();
4329                  if (info) {
4330                    int numberColumns = lpSolver->numberColumns();
4331                    int i;
4332                    for (i=0;i<numberColumns;i++) {
4333                      if (info[i]) 
4334                        clpSolver->setInteger(i);
4335                    }
4336                  }
4337                  // sets to all slack (not necessary?)
4338                  lpSolver->createStatus();
4339                  time2 = CoinCpuTime();
4340                  totalTime += time2-time1;
4341                  time1=time2;
4342                  // Go to canned file if just input file
4343                  if (CbcOrClpRead_mode==2&&argc==2) {
4344                    // only if ends .mps
4345                    std::string::size_type loc = fileName.find(".mps") ;
4346                    if (loc != std::string::npos &&
4347                        fileName.length() == loc+3)
4348                    { fileName.replace(loc+1,3,"par") ;
4349                      FILE *fp=fopen(fileName.c_str(),"r");
4350                      if (fp) {
4351                        CbcOrClpReadCommand=fp; // Read from that file
4352                        CbcOrClpRead_mode=-1;
4353                      }
4354                    }
4355                  }
4356                } else {
4357                  // errors
4358                  std::cout<<"There were "<<status<<
4359                    " errors on input"<<std::endl;
4360                }
4361              }
4362            }
4363            break;
4364          case MODELIN:
4365#ifdef COIN_HAS_LINK
4366            {
4367              // get next field
4368              field = CoinReadGetString(argc,argv);
4369              if (field=="$") {
4370                field = parameters[iParam].stringValue();
4371              } else if (field=="EOL") {
4372                parameters[iParam].printString();
4373                break;
4374              } else {
4375                parameters[iParam].setStringValue(field);
4376              }
4377              std::string fileName;
4378              bool canOpen=false;
4379              if (field=="-") {
4380                // stdin
4381                canOpen=true;
4382                fileName = "-";
4383              } else {
4384                bool absolutePath;
4385                if (dirsep=='/') {
4386                  // non Windows (or cygwin)
4387                  absolutePath=(field[0]=='/');
4388                } else {
4389                  //Windows (non cycgwin)
4390                  absolutePath=(field[0]=='\\');
4391                  // but allow for :
4392                  if (strchr(field.c_str(),':'))
4393                    absolutePath=true;
4394                }
4395                if (absolutePath) {
4396                  fileName = field;
4397                } else if (field[0]=='~') {
4398                  char * environVar = getenv("HOME");
4399                  if (environVar) {
4400                    std::string home(environVar);
4401                    field=field.erase(0,1);
4402                    fileName = home+field;
4403                  } else {
4404                    fileName=field;
4405                  }
4406                } else {
4407                  fileName = directory+field;
4408                }
4409                FILE *fp=fopen(fileName.c_str(),"r");
4410                if (fp) {
4411                  // can open - lets go for it
4412                  fclose(fp);
4413                  canOpen=true;
4414                } else {
4415                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4416                }
4417              }
4418              if (canOpen) {
4419                CoinModel coinModel(fileName.c_str(),2);
4420                // load from coin model
4421                OsiSolverLink solver1;
4422                OsiSolverInterface * solver2 = solver1.clone();
4423                model.assignSolver(solver2,true);
4424                OsiSolverLink * si =
4425                  dynamic_cast<OsiSolverLink *>(model.solver()) ;
4426                assert (si != NULL);
4427                si->setDefaultMeshSize(0.001);
4428                // need some relative granularity
4429                si->setDefaultBound(100.0);
4430                si->setDefaultMeshSize(0.01);
4431                si->setDefaultBound(100.0);
4432                si->setIntegerPriority(1000);
4433                si->setBiLinearPriority(10000);
4434                CoinModel * model2 = (CoinModel *) &coinModel;
4435                si->load(*model2);
4436                // redo
4437                solver = model.solver();
4438                clpSolver = dynamic_cast< OsiClpSolverInterface*> (solver);
4439                lpSolver = clpSolver->getModelPtr();
4440                clpSolver->messageHandler()->setLogLevel(0) ;
4441                testOsiParameters=0;
4442                complicatedInteger=2;
4443              }
4444            }
4445#endif
4446            break;
4447          case EXPORT:
4448            if (goodModel) {
4449              // get next field
4450              field = CoinReadGetString(argc,argv);
4451              if (field=="$") {
4452                field = parameters[iParam].stringValue();
4453              } else if (field=="EOL") {
4454                parameters[iParam].printString();
4455                break;
4456              } else {
4457                parameters[iParam].setStringValue(field);
4458              }
4459              std::string fileName;
4460              bool canOpen=false;
4461              if (field[0]=='/'||field[0]=='\\') {
4462                fileName = field;
4463              } else if (field[0]=='~') {
4464                char * environVar = getenv("HOME");
4465                if (environVar) {
4466                  std::string home(environVar);
4467                  field=field.erase(0,1);
4468                  fileName = home+field;
4469                } else {
4470                  fileName=field;
4471                }
4472              } else {
4473                fileName = directory+field;
4474              }
4475              FILE *fp=fopen(fileName.c_str(),"w");
4476              if (fp) {
4477                // can open - lets go for it
4478                fclose(fp);
4479                canOpen=true;
4480              } else {
4481                std::cout<<"Unable to open file "<<fileName<<std::endl;
4482              }
4483              if (canOpen) {
4484                // If presolve on then save presolved
4485                bool deleteModel2=false;
4486                ClpSimplex * model2 = lpSolver;
4487                if (dualize) {
4488                  model2 = ((ClpSimplexOther *) model2)->dualOfModel();
4489                  printf("Dual of model has %d rows and %d columns\n",
4490                         model2->numberRows(),model2->numberColumns());
4491                  model2->setOptimizationDirection(1.0);
4492                }
4493#ifdef COIN_HAS_ASL
4494                if (info.numberSos&&doSOS&&usingAmpl) {
4495                  // SOS
4496                  numberSOS = info.numberSos;
4497                  sosStart = info.sosStart;
4498                  sosIndices = info.sosIndices;
4499                  sosReference = info.sosReference;
4500                  preSolve=false;
4501                  clpSolver->setSOSData(numberSOS,info.sosType,sosStart,sosIndices,sosReference);
4502                }
4503#endif
4504                if (preSolve) {
4505                  ClpPresolve pinfo;
4506                  int presolveOptions2 = presolveOptions&~0x40000000;
4507                  if ((presolveOptions2&0xffff)!=0)
4508                    pinfo.setPresolveActions(presolveOptions2);
4509                  if ((printOptions&1)!=0)
4510                    pinfo.statistics();
4511                  double presolveTolerance = 
4512                    parameters[whichParam(PRESOLVETOLERANCE,numberParameters,parameters)].doubleValue();
4513                  model2 = 
4514                    pinfo.presolvedModel(*lpSolver,presolveTolerance,
4515                                         true,preSolve);
4516                  if (model2) {
4517                    printf("Saving presolved model on %s\n",
4518                           fileName.c_str());
4519                    deleteModel2=true;
4520                  } else {
4521                    printf("Presolved model looks infeasible - saving original on %s\n",
4522                           fileName.c_str());
4523                    deleteModel2=false;
4524                    model2 = lpSolver;
4525
4526                  }
4527                  model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
4528                  if (deleteModel2)
4529                    delete model2;
4530                } else {
4531                  printf("Saving model on %s\n",
4532                           fileName.c_str());
4533                  if (numberSOS) {
4534                    // Convert names
4535                    int iRow;
4536                    int numberRows=model2->numberRows();
4537                    int iColumn;
4538                    int numberColumns=model2->numberColumns();
4539                   
4540                    char ** rowNames = NULL;
4541                    char ** columnNames = NULL;
4542                    if (model2->lengthNames()) {
4543                      rowNames = new char * [numberRows];
4544                      for (iRow=0;iRow<numberRows;iRow++) {
4545                        rowNames[iRow] = 
4546                          strdup(model2->rowName(iRow).c_str());
4547                      }
4548                     
4549                      columnNames = new char * [numberColumns];
4550                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
4551                        columnNames[iColumn] = 
4552                          strdup(model2->columnName(iColumn).c_str());
4553                      }
4554                    }
4555                    clpSolver->writeMpsNative(fileName.c_str(),(const char **) rowNames,(const char **) columnNames,
4556                                              (outputFormat-1)/2,1+((outputFormat-1)&1));
4557                    if (rowNames) {
4558                      for (iRow=0;iRow<numberRows;iRow++) {
4559                        free(rowNames[iRow]);
4560                      }
4561                      delete [] rowNames;
4562                      for (iColumn=0;iColumn<numberColumns;iColumn++) {
4563                        free(columnNames[iColumn]);
4564                      }
4565                      delete [] columnNames;
4566                    }
4567                  } else {
4568                    model2->writeMps(fileName.c_str(),(outputFormat-1)/2,1+((outputFormat-1)&1));
4569                  }
4570                }
4571                time2 = CoinCpuTime();
4572                totalTime += time2-time1;
4573                time1=time2;
4574              }
4575            } else {
4576              std::cout<<"** Current model not valid"<<std::endl;
4577            }
4578            break;
4579          case BASISIN:
4580            if (goodModel) {
4581              // get next field
4582              field = CoinReadGetString(argc,argv);
4583              if (field=="$") {
4584                field = parameters[iParam].stringValue();
4585              } else if (field=="EOL") {
4586                parameters[iParam].printString();
4587                break;
4588              } else {
4589                parameters[iParam].setStringValue(field);
4590              }
4591              std::string fileName;
4592              bool canOpen=false;
4593              if (field=="-") {
4594                // stdin
4595                canOpen=true;
4596                fileName = "-";
4597              } else {
4598                if (field[0]=='/'||field[0]=='\\') {
4599                  fileName = field;
4600                } else if (field[0]=='~') {
4601                  char * environVar = getenv("HOME");
4602                  if (environVar) {
4603                    std::string home(environVar);
4604                    field=field.erase(0,1);
4605                    fileName = home+field;
4606                  } else {
4607                    fileName=field;
4608                  }
4609                } else {
4610                  fileName = directory+field;
4611                }
4612                FILE *fp=fopen(fileName.c_str(),"r");
4613                if (fp) {
4614                  // can open - lets go for it
4615                  fclose(fp);
4616                  canOpen=true;
4617                } else {
4618                  std::cout<<"Unable to open file "<<fileName<<std::endl;
4619                }
4620              }
4621              if (canOpen) {
4622                int values = lpSolver->readBasis(fileName.c_str());
4623                if (values==0)
4624                  basisHasValues=-1;
4625                else
4626                  basisHasValues=1;
4627              }
4628            } else {
4629              std::cout<<"** Current model not valid"<<std::endl;
4630            }
4631            break;
4632          case PRIORITYIN:
4633            if (goodModel) {
4634              // get next field
4635              field = CoinReadGetString(argc,argv);
4636              if (field=="$") {
4637                field = parameters[iParam].stringValue();
4638              } else if (field=="EOL") {
4639                parameters[iParam].printString();
4640                break;
4641              } else {
4642                parameters[iParam].setStringValue(field);
4643              }
4644              std::string fileName;
4645              if (field[0]=='/'||field[0]=='\\') {
4646                fileName = field;
4647              } else if (field[0]=='~') {
4648                char * environVar = getenv("HOME");
4649                if (environVar) {
4650                  std::string home(environVar);
4651                  field=field.erase(0,1);
4652                  fileName = home+field;
4653                } else {
4654                  fileName=field;
4655                }
4656              } else {
4657                fileName = directory+field;
4658              }
4659              FILE *fp=fopen(fileName.c_str(),"r");
4660              if (fp) {
4661                // can open - lets go for it
4662                std::string headings[]={"name","number","direction","priority","up","down",
4663                                        "solution","priin"};
4664                int got[]={-1,-1,-1,-1,-1,-1,-1,-1};
4665                int order[8];
4666                assert(sizeof(got)==sizeof(order));
4667                int nAcross=0;
4668                char line[1000];
4669                int numberColumns = lpSolver->numberColumns();
4670                if (!fgets(line,1000,fp)) {
4671                  std::cout<<"Odd file "<<fileName<<std::endl;
4672                } else {
4673                  char * pos = line;
4674                  char * put = line;
4675                  while (*pos>=' '&&*pos!='\n') {
4676                    if (*pos!=' '&&*pos!='\t') {
4677                      *put=tolower(*pos);
4678                      put++;
4679                    }
4680                    pos++;
4681                  }
4682                  *put='\0';
4683                  pos=line;
4684                  int i;
4685                  bool good=true;
4686                  while (pos) {
4687                    char * comma = strchr(pos,',');
4688                    if (comma)
4689                      *comma='\0';
4690                    for (i=0;i<(int) (sizeof(got)/sizeof(int));i++) {
4691                      if (headings[i]==pos) {
4692                        if (got[i]<0) {
4693                          order[nAcross]=i;
4694                          got[i]=nAcross++;
4695                        } else {
4696                          // duplicate
4697                          good=false;
4698                        }
4699                        break;
4700                      }
4701                    }
4702                    if (i==(int) (sizeof(got)/sizeof(int)))
4703                      good=false;
4704                    if (comma) {
4705                      *comma=',';
4706                      pos=comma+1;
4707                    } else {
4708                      break;
4709                    }
4710                  }
4711                  if (got[0]<0&&got[1]<0)
4712                    good=false;
4713                  if (got[0]>=0&&got[1]>=0)
4714                    good=false;
4715                  if (got[0]>=0&&!lpSolver->lengthNames())
4716                    good=false;
4717                  if (good) {
4718                    char ** columnNames = new char * [numberColumns];
4719                    pseudoDown= (double *) malloc(numberColumns*sizeof(double));
4720                    pseudoUp = (double *) malloc(numberColumns*sizeof(double));
4721                    branchDirection = (int *) malloc(numberColumns*sizeof(int));
4722                    priorities= (int *) malloc(numberColumns*sizeof(int));
4723                    free(solutionIn);
4724                    solutionIn=NULL;
4725                    free(prioritiesIn);
4726                    prioritiesIn=NULL;
4727                    int iColumn;
4728                    if (got[6]>=0) {
4729                      solutionIn = (double *) malloc(numberColumns*sizeof(double));
4730                      CoinZeroN(solutionIn,numberColumns);
4731                    }
4732                    if (got[7]>=0) {
4733                      prioritiesIn = (int *) malloc(numberColumns*sizeof(int));
4734                      for (iColumn=0;iColumn<numberColumns;iColumn++) 
4735                        prioritiesIn[iColumn]=10000;
4736                    }
4737                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4738                      columnNames[iColumn] = 
4739                        strdup(lpSolver->columnName(iColumn).c_str());
4740                      pseudoDown[iColumn]=0.0;
4741                      pseudoUp[iColumn]=0.0;
4742                      branchDirection[iColumn]=0;
4743                      priorities[iColumn]=0;
4744                    }
4745                    int nBadPseudo=0;
4746                    int nBadDir=0;
4747                    int nBadPri=0;
4748                    int nBadName=0;
4749                    int nBadLine=0;
4750                    int nLine=0;
4751                    while (fgets(line,1000,fp)) {
4752                      nLine++;
4753                      iColumn = -1;
4754                      double up =0.0;
4755                      double down=0.0;
4756                      int pri=0;
4757                      int dir=0;
4758                      double solValue=COIN_DBL_MAX;
4759                      int priValue=1000000;
4760                      char * pos = line;
4761                      char * put = line;
4762                      while (*pos>=' '&&*pos!='\n') {
4763                        if (*pos!=' '&&*pos!='\t') {
4764                          *put=tolower(*pos);
4765                          put++;
4766                        }
4767                        pos++;
4768                      }
4769                      *put='\0';
4770                      pos=line;
4771                      for (int i=0;i<nAcross;i++) {
4772                        char * comma = strchr(pos,',');
4773                        if (comma) {
4774                          *comma='\0';
4775                        } else if (i<nAcross-1) {
4776                          nBadLine++;
4777                          break;
4778                        }
4779                        switch (order[i]) {
4780                          // name
4781                        case 0:
4782                          for (iColumn=0;iColumn<numberColumns;iColumn++) {
4783                            if (!strcmp(columnNames[iColumn],pos))
4784                              break;
4785                          }
4786                          if (iColumn==numberColumns)
4787                            iColumn=-1;
4788                          break;
4789                          // number
4790                        case 1:
4791                          iColumn = atoi(pos);
4792                          if (iColumn<0||iColumn>=numberColumns)
4793                            iColumn=-1;
4794                          break;
4795                          // direction
4796                        case 2:
4797                          if (*pos=='D')
4798                            dir=-1;
4799                          else if (*pos=='U')
4800                            dir=1;
4801                          else if (*pos=='N')
4802                            dir=0;
4803                          else if (*pos=='1'&&*(pos+1)=='\0')
4804                            dir=1;
4805                          else if (*pos=='0'&&*(pos+1)=='\0')
4806                            dir=0;
4807                          else if (*pos=='1'&&*(pos+1)=='1'&&*(pos+2)=='\0')
4808                            dir=-1;
4809                          else
4810                            dir=-2; // bad
4811                          break;
4812                          // priority
4813                        case 3:
4814                          pri=atoi(pos);
4815                          break;
4816                          // up
4817                        case 4:
4818                          up = atof(pos);
4819                          break;
4820                          // down
4821                        case 5:
4822                          down = atof(pos);
4823                          break;
4824                          // sol value
4825                        case 6:
4826                          solValue = atof(pos);
4827                          break;
4828                          // priority in value
4829                        case 7:
4830                          priValue = atoi(pos);
4831                          break;
4832                        }
4833                        if (comma) {
4834                          *comma=',';
4835                          pos=comma+1;
4836                        }
4837                      }
4838                      if (iColumn>=0) {
4839                        if (down<0.0) {
4840                          nBadPseudo++;
4841                          down=0.0;
4842                        }
4843                        if (up<0.0) {
4844                          nBadPseudo++;
4845                          up=0.0;
4846                        }
4847                        if (!up)
4848                          up=down;
4849                        if (!down)
4850                          down=up;
4851                        if (dir<-1||dir>1) {
4852                          nBadDir++;
4853                          dir=0;
4854                        }
4855                        if (pri<0) {
4856                          nBadPri++;
4857                          pri=0;
4858                        }
4859                        pseudoDown[iColumn]=down;
4860                        pseudoUp[iColumn]=up;
4861                        branchDirection[iColumn]=dir;
4862                        priorities[iColumn]=pri;
4863                        if (solValue!=COIN_DBL_MAX) {
4864                          assert (solutionIn);
4865                          solutionIn[iColumn]=solValue;
4866                        }
4867                        if (priValue!=1000000) {
4868                          assert (prioritiesIn);
4869                          prioritiesIn[iColumn]=priValue;
4870                        }
4871                      } else {
4872                        nBadName++;
4873                      }
4874                    }
4875                    if (!noPrinting) {
4876                      printf("%d fields and %d records",nAcross,nLine);
4877                      if (nBadPseudo)
4878                        printf(" %d bad pseudo costs",nBadPseudo);
4879                      if (nBadDir)
4880                        printf(" %d bad directions",nBadDir);
4881                      if (nBadPri)
4882                        printf(" %d bad priorities",nBadPri);
4883                      if (nBadName)
4884                        printf(" ** %d records did not match on name/sequence",nBadName);
4885                      printf("\n");
4886                    }
4887                    for (iColumn=0;iColumn<numberColumns;iColumn++) {
4888                      free(columnNames[iColumn]);
4889                    }
4890                    delete [] columnNames;
4891                  } else {
4892                    std::cout<<"Duplicate or unknown keyword - or name/number fields wrong"<<line<<std::endl;
4893                  }
4894                }
4895                fclose(fp);
4896              } else {
4897                std::cout<<"Unable to open file "<<fileName<<std::endl;
4898              }
4899            } else {
4900              std::cout<<"** Current model not valid"<<std::endl;
4901            }
4902            break;
4903          case DEBUG:
4904            if (goodModel) {
4905              delete [] debugValues;
4906              debugValues=NULL;
4907              // get next field
4908              field = CoinReadGetString(argc,argv);
4909              if (field=="$") {
4910                field = parameters[iParam].stringValue();
4911              } else if (field=="EOL") {
4912                parameters[iParam].printString();
4913                break;
4914              } else {
4915                parameters[iParam].setStringValue(field);
4916                debugFile=field;
4917                if (debugFile=="create"||
4918                    debugFile=="createAfterPre") {
4919                  printf("Will create a debug file so this run should be a good one\n");
4920                  break;
4921                }
4922              }
4923              std::string fileName;
4924              if (field[0]=='/'||field[0]=='\\') {
4925                fileName = field;
4926              } else if (field[0]=='~') {
4927                char * environVar = getenv("HOME");
4928                if (environVar) {
4929                  std::string home(environVar);
4930                  field=field.erase(0,1);
4931                  fileName = home+field;
4932                } else {
4933                  fileName=field;
4934                }
4935              } else {
4936                fileName = directory+field;
4937              }
4938              FILE *fp=fopen(fileName.c_str(),"rb");
4939              if (fp) {
4940                // can open - lets go for it
4941                int numRows;
4942                double obj;
4943                fread(&numRows,sizeof(int),1,fp);
4944                fread(&numberDebugValues,sizeof(int),1,fp);
4945                fread(&obj,sizeof(double),1,fp);
4946                debugValues = new double[numberDebugValues+numRows];
4947                fread(debugValues,sizeof(double),numRows,fp);
4948                fread(debugValues,sizeof(double),numRows,fp);
4949                fread(debugValues,sizeof(double),numberDebugValues,fp);
4950                printf("%d doubles read into debugValues\n",numberDebugValues);
4951                fclose(fp);
4952              } else {
4953                std::cout<<"Unable to open file "<<fileName<<std::endl;
4954              }
4955            } else {
4956              std::cout<<"** Current model not valid"<<std::endl;
4957            }
4958            break;
4959          case PRINTMASK:
4960            // get next field
4961            {
4962              std::string name = CoinReadGetString(argc,argv);
4963              if (name!="EOL") {
4964                parameters[iParam].setStringValue(name);
4965                printMask = name;
4966              } else {
4967                parameters[iParam].printString();
4968              }
4969            }
4970            break;
4971          case BASISOUT:
4972            if (goodModel) {
4973              // get next field
4974              field = CoinReadGetString(argc,argv);
4975              if (field=="$") {
4976                field = parameters[iParam].stringValue();
4977              } else if (field=="EOL") {
4978                parameters[iParam].printString();
4979                break;
4980              } else {
4981                parameters[iParam].setStringValue(field);
4982              }
4983              std::string fileName;
4984              bool canOpen=false;
4985              if (field[0]=='/'||field[0]=='\\') {
4986                fileName = field;
4987              } else if (field[0]=='~') {
4988                char * environVar = getenv("HOME");
4989                if (environVar) {
4990                  std::string home(environVar);
4991                  field=field.erase(0,1);
4992                  fileName = home+field;
4993                } else {
4994                  fileName=field;
4995                }
4996              } else {
4997                fileName = directory+field;
4998              }
4999              FILE *fp=fopen(fileName.c_str(),"w");
5000              if (fp) {
5001                // can open - lets go for it
5002                fclose(fp);
5003                canOpen=true;
5004              } else {
5005                std::cout<<"Unable to open file "<<fileName<<std::endl;
5006              }
5007              if (canOpen) {
5008