1 // $Id$
2 // Copyright (C) 2000, International Business Machines
3 // Corporation and others.  All Rights Reserved.
4 // This code is licensed under the terms of the Eclipse Public License (EPL).
5 
6 #include <cstdlib>
7 #include <cstdio>
8 #include <cmath>
9 #include <cfloat>
10 #include <cassert>
11 #include <iostream>
12 
13 #include "CoinPragma.hpp"
14 #include "CglTreeInfo.hpp"
15 #include "CoinHelperFunctions.hpp"
16 #include "CoinSort.hpp"
17 #include "CoinPackedMatrix.hpp"
18 #include "CglStored.hpp"
19 #include "OsiRowCut.hpp"
20 
21 // Default constructor
CglTreeInfo()22 CglTreeInfo::CglTreeInfo()
23   : level(-1)
24   , pass(-1)
25   , formulation_rows(-1)
26   , options(0)
27   , inTree(false)
28   , hasParent(0)
29   , parentSolver(NULL)
30   , originalColumns(NULL)
31   , strengthenRow(NULL)
32   , randomNumberGenerator(NULL)
33 {
34 }
35 
36 // Copy constructor
CglTreeInfo(const CglTreeInfo & rhs)37 CglTreeInfo::CglTreeInfo(const CglTreeInfo &rhs)
38   : level(rhs.level)
39   , pass(rhs.pass)
40   , formulation_rows(rhs.formulation_rows)
41   , options(rhs.options)
42   , inTree(rhs.inTree)
43   , hasParent(rhs.hasParent)
44   , parentSolver(rhs.parentSolver)
45   , originalColumns(rhs.originalColumns)
46   , strengthenRow(rhs.strengthenRow)
47   , randomNumberGenerator(rhs.randomNumberGenerator)
48 {
49 }
50 // Clone
51 CglTreeInfo *
clone() const52 CglTreeInfo::clone() const
53 {
54   return new CglTreeInfo(*this);
55 }
56 
57 // Assignment operator
58 CglTreeInfo &
operator =(const CglTreeInfo & rhs)59 CglTreeInfo::operator=(const CglTreeInfo &rhs)
60 {
61   if (this != &rhs) {
62     //CglCutGenerator::operator=(rhs);
63     level = rhs.level;
64     pass = rhs.pass;
65     formulation_rows = rhs.formulation_rows;
66     options = rhs.options;
67     inTree = rhs.inTree;
68     hasParent = rhs.hasParent;
69     parentSolver = rhs.parentSolver;
70     originalColumns = rhs.originalColumns;
71     strengthenRow = rhs.strengthenRow;
72     randomNumberGenerator = rhs.randomNumberGenerator;
73   }
74   return *this;
75 }
76 
77 // Destructor
~CglTreeInfo()78 CglTreeInfo::~CglTreeInfo()
79 {
80 }
81 
82 // Default constructor
CglTreeProbingInfo()83 CglTreeProbingInfo::CglTreeProbingInfo()
84   : CglTreeInfo()
85   , fixEntry_(NULL)
86   , toZero_(NULL)
87   , toOne_(NULL)
88   , integerVariable_(NULL)
89   , backward_(NULL)
90   , fixingEntry_(NULL)
91   , numberVariables_(0)
92   , numberIntegers_(0)
93   , maximumEntries_(0)
94   , numberEntries_(-1)
95 {
96 }
97 // Constructor from model
CglTreeProbingInfo(const OsiSolverInterface * model)98 CglTreeProbingInfo::CglTreeProbingInfo(const OsiSolverInterface *model)
99   : CglTreeInfo()
100   , fixEntry_(NULL)
101   , toZero_(NULL)
102   , toOne_(NULL)
103   , integerVariable_(NULL)
104   , backward_(NULL)
105   , fixingEntry_(NULL)
106   , numberVariables_(0)
107   , numberIntegers_(0)
108   , maximumEntries_(0)
109   , numberEntries_(-1)
110 {
111   numberVariables_ = model->getNumCols();
112   // Too many ... but
113   integerVariable_ = new int[numberVariables_];
114   backward_ = new int[numberVariables_];
115   int i;
116   // Get integer types
117   const char *columnType = model->getColType(true);
118   for (i = 0; i < numberVariables_; i++) {
119     backward_[i] = -1;
120     if (columnType[i]) {
121       if (columnType[i] == 1) {
122         backward_[i] = numberIntegers_;
123         integerVariable_[numberIntegers_++] = i;
124       } else {
125         backward_[i] = -2;
126       }
127     }
128   }
129   // Set up to arrays
130   toOne_ = new int[numberIntegers_];
131   toZero_ = new int[numberIntegers_ + 1];
132   // zero out
133   CoinZeroN(toOne_, numberIntegers_);
134   CoinZeroN(toZero_, numberIntegers_ + 1);
135 }
136 
137 // Copy constructor
CglTreeProbingInfo(const CglTreeProbingInfo & rhs)138 CglTreeProbingInfo::CglTreeProbingInfo(const CglTreeProbingInfo &rhs)
139   : CglTreeInfo(rhs)
140   , fixEntry_(NULL)
141   , toZero_(NULL)
142   , toOne_(NULL)
143   , integerVariable_(NULL)
144   , backward_(NULL)
145   , fixingEntry_(NULL)
146   , numberVariables_(rhs.numberVariables_)
147   , numberIntegers_(rhs.numberIntegers_)
148   , maximumEntries_(rhs.maximumEntries_)
149   , numberEntries_(rhs.numberEntries_)
150 {
151   if (numberVariables_) {
152     fixEntry_ = new CliqueEntry[maximumEntries_];
153     memcpy(fixEntry_, rhs.fixEntry_, maximumEntries_ * sizeof(CliqueEntry));
154     if (numberEntries_ < 0) {
155       // in order
156       toZero_ = CoinCopyOfArray(rhs.toZero_, numberIntegers_ + 1);
157       toOne_ = CoinCopyOfArray(rhs.toOne_, numberIntegers_);
158     } else {
159       // not in order
160       fixingEntry_ = CoinCopyOfArray(rhs.fixingEntry_, maximumEntries_);
161     }
162     integerVariable_ = CoinCopyOfArray(rhs.integerVariable_, numberIntegers_);
163     backward_ = CoinCopyOfArray(rhs.backward_, numberVariables_);
164   }
165 }
166 // Clone
167 CglTreeInfo *
clone() const168 CglTreeProbingInfo::clone() const
169 {
170   return new CglTreeProbingInfo(*this);
171 }
172 
173 // Assignment operator
174 CglTreeProbingInfo &
operator =(const CglTreeProbingInfo & rhs)175 CglTreeProbingInfo::operator=(const CglTreeProbingInfo &rhs)
176 {
177   if (this != &rhs) {
178     CglTreeInfo::operator=(rhs);
179     delete[] fixEntry_;
180     delete[] toZero_;
181     delete[] toOne_;
182     delete[] integerVariable_;
183     delete[] backward_;
184     delete[] fixingEntry_;
185     numberVariables_ = rhs.numberVariables_;
186     numberIntegers_ = rhs.numberIntegers_;
187     maximumEntries_ = rhs.maximumEntries_;
188     numberEntries_ = rhs.numberEntries_;
189     if (numberVariables_) {
190       fixEntry_ = new CliqueEntry[maximumEntries_];
191       memcpy(fixEntry_, rhs.fixEntry_, maximumEntries_ * sizeof(CliqueEntry));
192       if (numberEntries_ < 0) {
193         // in order
194         toZero_ = CoinCopyOfArray(rhs.toZero_, numberIntegers_ + 1);
195         toOne_ = CoinCopyOfArray(rhs.toOne_, numberIntegers_);
196         fixingEntry_ = NULL;
197       } else {
198         // not in order
199         fixingEntry_ = CoinCopyOfArray(rhs.fixingEntry_, maximumEntries_);
200         toZero_ = NULL;
201         toOne_ = NULL;
202       }
203       toZero_ = CoinCopyOfArray(rhs.toZero_, numberIntegers_ + 1);
204       toOne_ = CoinCopyOfArray(rhs.toOne_, numberIntegers_);
205       integerVariable_ = CoinCopyOfArray(rhs.integerVariable_, numberIntegers_);
206       backward_ = CoinCopyOfArray(rhs.backward_, numberVariables_);
207     } else {
208       fixEntry_ = NULL;
209       toZero_ = NULL;
210       toOne_ = NULL;
211       integerVariable_ = NULL;
212       backward_ = NULL;
213       fixingEntry_ = NULL;
214     }
215   }
216   return *this;
217 }
218 
219 // Destructor
~CglTreeProbingInfo()220 CglTreeProbingInfo::~CglTreeProbingInfo()
221 {
222   delete[] fixEntry_;
223   delete[] toZero_;
224   delete[] toOne_;
225   delete[] integerVariable_;
226   delete[] backward_;
227   delete[] fixingEntry_;
228 }
outDupsEtc(int numberIntegers,int & numberCliques,int & numberMatrixCliques,CoinBigIndex * & cliqueStart,char * & cliqueType,CliqueEntry * & entry,int numberLastTime,int printit)229 static int outDupsEtc(int numberIntegers, int &numberCliques, int &numberMatrixCliques,
230   CoinBigIndex *&cliqueStart, char *&cliqueType, CliqueEntry *&entry,
231   int numberLastTime, int printit)
232 {
233   bool allNew = false;
234   int *whichP = new int[numberIntegers];
235   int iClique;
236   assert(sizeof(int) == 4);
237   assert(sizeof(CliqueEntry) == 4);
238   // If lots then get rid of short ones
239 #define KEEP_CLIQUES 10000
240   if (numberCliques - numberMatrixCliques > KEEP_CLIQUES) {
241     int *sort = new int[numberCliques];
242     for (iClique = numberMatrixCliques; iClique < numberCliques; iClique++) {
243       CoinBigIndex j = cliqueStart[iClique];
244       int n = static_cast< int >(cliqueStart[iClique + 1] - j);
245       sort[iClique] = n;
246     }
247     std::sort(sort + numberMatrixCliques, sort + numberCliques);
248     int allow = sort[numberCliques - KEEP_CLIQUES];
249     int nEqual = 0;
250     for (iClique = numberCliques - KEEP_CLIQUES; iClique < numberCliques; iClique++) {
251       if (sort[iClique] > allow)
252         break;
253       else
254         nEqual++;
255     }
256     delete[] sort;
257     CoinBigIndex j = cliqueStart[numberMatrixCliques];
258     CoinBigIndex put = j;
259     int nClique = numberMatrixCliques;
260     for (iClique = numberMatrixCliques; iClique < numberCliques; iClique++) {
261       CoinBigIndex end = cliqueStart[iClique + 1];
262       int n = static_cast< int >(end - j);
263       bool copy = false;
264       if (n > allow) {
265         copy = true;
266       } else if (n == allow && nEqual) {
267         copy = true;
268         nEqual--;
269       }
270       if (copy) {
271         cliqueType[nClique++] = cliqueType[iClique];
272         for (; j < end; j++)
273           entry[put++] = entry[j];
274       }
275       j = cliqueStart[iClique + 1];
276       cliqueStart[nClique] = put;
277     }
278     numberCliques = nClique;
279   }
280   // sort
281   for (iClique = 0; iClique < numberCliques; iClique++) {
282     CoinBigIndex j = cliqueStart[iClique];
283     int n = static_cast< int >(cliqueStart[iClique + 1] - j);
284     for (int i = 0; i < n; i++)
285       whichP[i] = sequenceInCliqueEntry(entry[i + j]);
286     CoinSort_2(whichP, whichP + n, (reinterpret_cast< int * >(entry)) + j);
287   }
288   // lexicographic sort
289   int *which = new int[numberCliques];
290   CoinBigIndex *position = new CoinBigIndex[numberCliques];
291   int *sort = new int[numberCliques];
292   int *value = new int[numberCliques];
293   for (iClique = 0; iClique < numberCliques; iClique++) {
294     which[iClique] = iClique;
295     sort[iClique] = sequenceInCliqueEntry(entry[cliqueStart[iClique]]);
296     value[iClique] = sort[iClique];
297     position[iClique] = 0;
298   }
299   CoinSort_2(sort, sort + numberCliques, which);
300   int lastDone = -1;
301   int nDup = 0;
302   CoinBigIndex nSave = 0;
303   while (lastDone < numberCliques - 1) {
304     int jClique = lastDone + 1;
305     int jFirst = jClique;
306     int iFirst = which[jFirst];
307     int iValue = value[iFirst];
308     CoinBigIndex iPos = position[iFirst];
309     jClique++;
310     for (; jClique < numberCliques; jClique++) {
311       int kClique = which[jClique];
312       int jValue = value[kClique];
313       if (jValue > iValue || position[kClique] < iPos)
314         break;
315     }
316     if (jClique == jFirst + 1) {
317       // done that bit
318       lastDone++;
319     } else {
320       // use next bit to sort and then repeat
321       int jLast = jClique;
322       for (jClique = jFirst; jClique < jLast; jClique++) {
323         int kClique = which[jClique];
324         int iValue = value[kClique];
325         // put at end if finished
326         if (iValue < numberIntegers) {
327           CoinBigIndex kPos = position[kClique] + 1;
328           position[kClique] = kPos;
329           kPos += cliqueStart[kClique];
330           if (kPos == cliqueStart[kClique + 1]) {
331             iValue = numberIntegers;
332           } else {
333             iValue = sequenceInCliqueEntry(entry[kPos]);
334           }
335           value[kClique] = iValue;
336         }
337         sort[jClique] = iValue;
338       }
339       CoinSort_2(sort + jFirst, sort + jLast, which + jFirst);
340       // if duplicate mark and move on
341       int iLowest = numberCliques;
342       for (jClique = jFirst; jClique < jLast; jClique++) {
343         int kClique = which[jClique];
344         int iValue = value[kClique];
345         if (iValue < numberIntegers)
346           break;
347         iLowest = CoinMin(iLowest, kClique);
348       }
349       if (jClique > jFirst) {
350         // mark all apart from lowest number as duplicate and move on
351         lastDone = jClique - 1;
352         for (jClique = jFirst; jClique <= lastDone; jClique++) {
353           int kClique = which[jClique];
354           if (kClique != iLowest) {
355             value[kClique] = -2;
356             nDup++;
357             nSave += cliqueStart[kClique + 1] - cliqueStart[kClique];
358           }
359         }
360       }
361     }
362   }
363   if (printit)
364     printf("%d duplicates\n", nDup);
365   // Now see if any subset
366   int nOut = 0;
367   for (int jClique = 0; jClique < numberCliques; jClique++) {
368     if (value[jClique] != -2) {
369       position[jClique] = cliqueStart[jClique];
370       value[jClique] = sequenceInCliqueEntry(entry[cliqueStart[jClique]]);
371     }
372   }
373   nSave = 0;
374   int startLooking = 0;
375   for (int jClique = 0; jClique < numberCliques; jClique++) {
376     int kClique = which[jClique];
377     if (value[kClique] == -2) {
378       nOut++;
379       nSave += cliqueStart[kClique + 1] - cliqueStart[kClique];
380       if (jClique == startLooking)
381         startLooking++;
382       continue;
383     }
384     int kValue = value[kClique];
385     for (int iiClique = startLooking; iiClique < jClique; iiClique++) {
386       int iClique = which[iiClique];
387       int iValue = value[iClique];
388       if (iValue == -2 || iValue == numberIntegers) {
389         if (iiClique == startLooking)
390           startLooking++;
391         continue;
392       } else {
393         if (kValue > static_cast< int >(sequenceInCliqueEntry(entry[cliqueStart[iClique + 1] - 1]))) {
394           value[iClique] = numberIntegers;
395           continue;
396         }
397       }
398       if (iValue < kValue) {
399         while (iValue < kValue) {
400           CoinBigIndex iPos = position[iClique] + 1;
401           position[iClique] = iPos;
402           if (iPos == cliqueStart[iClique + 1]) {
403             iValue = numberIntegers;
404           } else {
405             iValue = sequenceInCliqueEntry(entry[iPos]);
406           }
407           value[iClique] = iValue;
408         }
409       }
410       if (iValue > kValue)
411         continue; // not a candidate
412       // See if subset (remember duplicates have gone)
413       if (cliqueStart[iClique + 1] - position[iClique] > cliqueStart[kClique + 1] - cliqueStart[kClique]) {
414         // could be subset ?
415         CoinBigIndex offset = cliqueStart[iClique] - position[kClique];
416         CoinBigIndex j;
417         bool subset = true;
418         // what about different fixes bool odd=false;
419         for (j = cliqueStart[kClique] + 1; j < cliqueStart[kClique + 1]; j++) {
420           int kColumn = sequenceInCliqueEntry(entry[j]);
421           int iColumn = sequenceInCliqueEntry(entry[j + offset]);
422           if (iColumn > kColumn) {
423             subset = false;
424           } else {
425             while (iColumn < kColumn) {
426               offset++;
427               if (j + offset < cliqueStart[iClique + 1]) {
428                 iColumn = sequenceInCliqueEntry(entry[j + offset]);
429               } else {
430                 subset = false;
431                 break;
432               }
433             }
434           }
435           if (!subset)
436             break;
437         }
438         if (subset) {
439           value[kClique] = -2;
440           if (printit > 1)
441             printf("clique %d is subset of %d\n", kClique, iClique);
442           nOut++;
443           break;
444         }
445       }
446     }
447   }
448   if (nOut) {
449     if (printit)
450       printf("Can get rid of %d cliques\n", nOut);
451     // make new copy
452     int nNewC = numberCliques - nOut;
453     int size = static_cast< int >(cliqueStart[numberCliques] - nSave);
454     int n = 0;
455     CoinBigIndex *start = new CoinBigIndex[nNewC + 1];
456     char *type = new char[nNewC];
457     start[0] = 0;
458     CliqueEntry *entryC = new CliqueEntry[size];
459     CoinBigIndex nel = 0;
460     allNew = true;
461     for (int jClique = 0; jClique < numberCliques; jClique++) {
462       int kClique = which[jClique];
463       if (value[kClique] != -2 && kClique < numberMatrixCliques) {
464         if (kClique >= numberLastTime)
465           allNew = false;
466         int nn = static_cast< int >(cliqueStart[kClique + 1] - cliqueStart[kClique]);
467         memcpy(entryC + nel, entry + cliqueStart[kClique], nn * sizeof(CliqueEntry));
468         nel += nn;
469         type[n++] = cliqueType[kClique];
470         start[n] = nel;
471       }
472     }
473     int nM = n;
474     for (int jClique = 0; jClique < numberCliques; jClique++) {
475       int kClique = which[jClique];
476       if (value[kClique] != -2 && kClique >= numberMatrixCliques) {
477         if (kClique >= numberLastTime)
478           allNew = false;
479         int nn = static_cast< int >(cliqueStart[kClique + 1] - cliqueStart[kClique]);
480         memcpy(entryC + nel, entry + cliqueStart[kClique], nn * sizeof(CliqueEntry));
481         nel += nn;
482         type[n++] = cliqueType[kClique];
483         start[n] = nel;
484       }
485     }
486     // move
487     numberCliques = n;
488     numberMatrixCliques = nM;
489     delete[] cliqueStart;
490     cliqueStart = start;
491     delete[] entry;
492     entry = entryC;
493     delete[] cliqueType;
494     cliqueType = type;
495     if (printit > 1) {
496       for (int jClique = 0; jClique < numberCliques; jClique++) {
497         printf("%d [ ", jClique);
498         for (CoinBigIndex i = cliqueStart[jClique]; i < cliqueStart[jClique + 1]; i++)
499           printf("%d(%d) ", sequenceInCliqueEntry(entry[i]), oneFixesInCliqueEntry(entry[i]));
500         printf("]\n");
501       }
502     }
503     if (printit)
504       printf("%d matrix cliques and %d found by probing\n", numberMatrixCliques, numberCliques - numberMatrixCliques);
505   }
506   delete[] value;
507   delete[] sort;
508   delete[] which;
509   delete[] position;
510   delete[] whichP;
511   if (!allNew)
512     return nOut;
513   else
514     return -1;
515 }
516 OsiSolverInterface *
analyze(const OsiSolverInterface & si,int createSolver,int numberExtraCliques,const CoinBigIndex * starts,const CliqueEntry * entries,const char * type)517 CglTreeProbingInfo::analyze(const OsiSolverInterface &si, int createSolver,
518   int numberExtraCliques, const CoinBigIndex *starts,
519   const CliqueEntry *entries, const char *type)
520 {
521   if (!createSolver)
522     return NULL;
523   convert();
524   if (!numberIntegers_)
525     return NULL;
526   bool alwaysDo = false;
527   if (numberExtraCliques < 0) {
528     alwaysDo = true;
529     numberExtraCliques = 0;
530   }
531   bool printit = false;
532   int numberCliques = 0;
533   CoinBigIndex numberEntries = 0;
534   CoinBigIndex *cliqueStart = NULL;
535   CliqueEntry *entry = NULL;
536   char *cliqueType = NULL;
537   int *whichP = new int[numberIntegers_];
538   int *whichM = new int[numberIntegers_];
539   int *whichClique = NULL;
540   int numberRows = si.getNumRows();
541   int numberMatrixCliques = 0;
542   const CoinPackedMatrix *rowCopy = si.getMatrixByRow();
543   assert(numberRows && si.getNumCols());
544   int iRow;
545   const int *column = rowCopy->getIndices();
546   const double *elementByRow = rowCopy->getElements();
547   const CoinBigIndex *rowStart = rowCopy->getVectorStarts();
548   const int *rowLength = rowCopy->getVectorLengths();
549   const double *lower = si.getColLower();
550   const double *upper = si.getColUpper();
551   const double *rowLower = si.getRowLower();
552   const double *rowUpper = si.getRowUpper();
553   for (int iPass = 0; iPass < 2; iPass++) {
554     if (iPass) {
555       CoinBigIndex numberExtraEntries = 0;
556       if (numberExtraCliques)
557         numberExtraEntries = starts[numberExtraCliques];
558       cliqueStart = new CoinBigIndex[numberCliques + 1 + numberExtraCliques];
559       cliqueStart[0] = 0;
560       entry = new CliqueEntry[numberEntries + numberExtraEntries];
561       cliqueType = new char[numberCliques + numberExtraCliques];
562       whichClique = new int[numberEntries + numberExtraEntries];
563       numberCliques = 0;
564       numberEntries = 0;
565     }
566 #if 1
567     for (iRow = 0; iRow < numberRows; iRow++) {
568       int numberP1 = 0, numberM1 = 0;
569       int numberTotal = 0;
570       CoinBigIndex j;
571       double upperValue = rowUpper[iRow];
572       double lowerValue = rowLower[iRow];
573       bool good = true;
574       for (j = rowStart[iRow]; j < rowStart[iRow] + rowLength[iRow]; j++) {
575         int iColumn = column[j];
576         double value = elementByRow[j];
577         if (upper[iColumn] - lower[iColumn] < 1.0e-8) {
578           // fixed
579           upperValue -= lower[iColumn] * value;
580           lowerValue -= lower[iColumn] * value;
581           continue;
582         } else if (backward_[iColumn] < 0) {
583           good = false;
584           break;
585         } else {
586           iColumn = backward_[iColumn];
587           numberTotal++;
588         }
589         if (fabs(value) != 1.0) {
590           good = false;
591         } else if (value > 0.0) {
592           assert(numberP1 < numberIntegers_);
593           whichP[numberP1++] = iColumn;
594           ;
595         } else {
596           assert(numberM1 < numberIntegers_);
597           whichM[numberM1++] = iColumn;
598         }
599       }
600       int iUpper = static_cast< int >(floor(upperValue + 1.0e-5));
601       int iLower = static_cast< int >(ceil(lowerValue - 1.0e-5));
602       int state = 0;
603       if (upperValue < 1.0e6) {
604         if (iUpper == 1 - numberM1)
605           state = 1;
606         else if (iUpper == -numberM1)
607           state = 2;
608         else if (iUpper < -numberM1)
609           state = 3;
610         if (fabs(static_cast< double >(iUpper) - upperValue) > 1.0e-9)
611           state = -1;
612       }
613       if (!state && lowerValue > -1.0e6) {
614         if (-iLower == 1 - numberP1)
615           state = -1;
616         else if (-iLower == -numberP1)
617           state = -2;
618         else if (-iLower < -numberP1)
619           state = -3;
620         if (fabs(static_cast< double >(iLower) - lowerValue) > 1.0e-9)
621           state = -1;
622       }
623       if (numberP1 + numberM1 < 2)
624         state = -1;
625       if (good && state > 0) {
626         if (abs(state) == 3) {
627           // infeasible
628           printf("FFF Infeasible\n");
629           ;
630           //feasible=false;
631           break;
632         } else if (abs(state) == 2) {
633           // we can fix all
634           //numberFixed += numberP1+numberM1;
635           printf("FFF can fix %d\n", numberP1 + numberM1);
636         } else {
637           for (j = 0; j < numberP1; j++) {
638             int iColumn = whichP[j];
639             if (iPass) {
640               CliqueEntry temp;
641               setOneFixesInCliqueEntry(temp, true);
642               setSequenceInCliqueEntry(temp, iColumn);
643               entry[numberEntries] = temp;
644             }
645             numberEntries++;
646           }
647           for (j = 0; j < numberM1; j++) {
648             int iColumn = whichM[j];
649             if (iPass) {
650               CliqueEntry temp;
651               setOneFixesInCliqueEntry(temp, false);
652               setSequenceInCliqueEntry(temp, iColumn);
653               entry[numberEntries] = temp;
654             }
655             numberEntries++;
656           }
657           if (iPass) {
658             if (iLower != iUpper) {
659               // slack
660               cliqueType[numberCliques] = 'S';
661             } else {
662               cliqueType[numberCliques] = 'E';
663             }
664             cliqueStart[numberCliques + 1] = numberEntries;
665           }
666           numberCliques++;
667         }
668       }
669     }
670 #endif
671     if (numberExtraCliques) {
672       CoinBigIndex numberExtraEntries = starts[numberExtraCliques];
673       memcpy(entry + numberEntries, entries, numberExtraEntries * sizeof(CliqueEntry));
674       for (int iClique = 0; iClique < numberExtraCliques; iClique++) {
675         cliqueType[numberCliques] = type[iClique];
676         numberCliques++;
677         cliqueStart[numberCliques] = starts[iClique] + numberEntries;
678       }
679       numberEntries += numberExtraEntries;
680     }
681     numberMatrixCliques = numberCliques;
682     //int size = toZero_[numberIntegers_];
683     //char * used = new char [size];
684     //memset(used,0,size);
685     // find two cliques
686     int nFix = 0;
687     for (int iColumn = 0; iColumn < static_cast< int >(numberIntegers_); iColumn++) {
688       int j;
689       for (j = toZero_[iColumn]; j < toOne_[iColumn]; j++) {
690         int jColumn = sequenceInCliqueEntry(fixEntry_[j]);
691         // just look at ones beore (this also skips non 0-1)
692         if (jColumn < iColumn) {
693           int k;
694           for (k = toZero_[jColumn]; k < toOne_[jColumn]; k++) {
695             if (sequenceInCliqueEntry(fixEntry_[k]) == (iColumn)) {
696               if (oneFixesInCliqueEntry(fixEntry_[j])) {
697                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
698                   if (printit && !iPass)
699                     printf("%d to zero implies %d to one and %d to zero implies %d to one\n",
700                       iColumn, jColumn, jColumn, iColumn);
701                   //0-0 illegal
702                   if (iPass) {
703                     CliqueEntry temp;
704                     setOneFixesInCliqueEntry(temp, false);
705                     setSequenceInCliqueEntry(temp, iColumn);
706                     entry[numberEntries] = temp;
707                   }
708                   numberEntries++;
709                   if (iPass) {
710                     CliqueEntry temp;
711                     setOneFixesInCliqueEntry(temp, false);
712                     setSequenceInCliqueEntry(temp, jColumn);
713                     entry[numberEntries] = temp;
714                   }
715                   numberEntries++;
716                   if (iPass) {
717                     // slack
718                     cliqueType[numberCliques] = 'S';
719                     cliqueStart[numberCliques + 1] = numberEntries;
720                   }
721                   numberCliques++;
722                 } else {
723                   if (printit && !iPass)
724                     printf("%d to zero implies %d to one and %d to zero implies %d to zero\n",
725                       iColumn, jColumn, jColumn, iColumn);
726                   nFix++; // jColumn is 1
727                 }
728               } else {
729                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
730                   if (printit && !iPass)
731                     printf("%d to zero implies %d to zero and %d to zero implies %d to one\n",
732                       iColumn, jColumn, jColumn, iColumn);
733                   nFix++; // iColumn is 1
734                 } else {
735                   if (printit && !iPass)
736                     printf("%d to zero implies %d to zero and %d to zero implies %d to zero\n",
737                       iColumn, jColumn, jColumn, iColumn);
738                   nFix++; // jColumn=iColumn
739                 }
740               }
741             }
742           }
743           for (k = toOne_[jColumn]; k < toZero_[jColumn + 1]; k++) {
744             if (sequenceInCliqueEntry(fixEntry_[k]) == (iColumn)) {
745               if (oneFixesInCliqueEntry(fixEntry_[j])) {
746                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
747                   if (printit && !iPass)
748                     printf("%d to zero implies %d to one and %d to one implies %d to one\n",
749                       iColumn, jColumn, jColumn, iColumn);
750                   nFix++; //iColumn is 1
751                 } else {
752                   if (printit && !iPass)
753                     printf("%d to zero implies %d to one and %d to one implies %d to zero\n",
754                       iColumn, jColumn, jColumn, iColumn);
755                   nFix++; // iColumn+jcolumn=1
756                 }
757               } else {
758                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
759                   if (printit && !iPass)
760                     printf("%d to zero implies %d to zero and %d to one implies %d to one\n",
761                       iColumn, jColumn, jColumn, iColumn);
762                   // 0-1 illegal
763                   if (iPass) {
764                     CliqueEntry temp;
765                     setOneFixesInCliqueEntry(temp, false);
766                     setSequenceInCliqueEntry(temp, iColumn);
767                     entry[numberEntries] = temp;
768                   }
769                   numberEntries++;
770                   if (iPass) {
771                     CliqueEntry temp;
772                     setOneFixesInCliqueEntry(temp, true);
773                     setSequenceInCliqueEntry(temp, jColumn);
774                     entry[numberEntries] = temp;
775                   }
776                   numberEntries++;
777                   if (iPass) {
778                     // slack
779                     cliqueType[numberCliques] = 'S';
780                     cliqueStart[numberCliques + 1] = numberEntries;
781                   }
782                   numberCliques++;
783                 } else {
784                   if (printit && !iPass)
785                     printf("%d to zero implies %d to zero and %d to one implies %d to zero\n",
786                       iColumn, jColumn, jColumn, iColumn);
787                   nFix++; // jColumn is 0
788                 }
789               }
790             }
791           }
792         }
793       }
794       for (j = toOne_[iColumn]; j < toZero_[iColumn + 1]; j++) {
795         int jColumn = sequenceInCliqueEntry(fixEntry_[j]);
796         if (jColumn < iColumn) {
797           int k;
798           for (k = toZero_[jColumn]; k < toOne_[jColumn]; k++) {
799             if (sequenceInCliqueEntry(fixEntry_[k]) == (iColumn)) {
800               if (oneFixesInCliqueEntry(fixEntry_[j])) {
801                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
802                   if (printit && !iPass)
803                     printf("%d to one implies %d to one and %d to zero implies %d to one\n",
804                       iColumn, jColumn, jColumn, iColumn);
805                   nFix++; // jColumn is 1
806                 } else {
807                   if (printit && !iPass)
808                     printf("%d to one implies %d to one and %d to zero implies %d to zero\n",
809                       iColumn, jColumn, jColumn, iColumn);
810                   // 1-0 illegal
811                   if (iPass) {
812                     CliqueEntry temp;
813                     setOneFixesInCliqueEntry(temp, true);
814                     setSequenceInCliqueEntry(temp, iColumn);
815                     entry[numberEntries] = temp;
816                   }
817                   numberEntries++;
818                   if (iPass) {
819                     CliqueEntry temp;
820                     setOneFixesInCliqueEntry(temp, false);
821                     setSequenceInCliqueEntry(temp, jColumn);
822                     entry[numberEntries] = temp;
823                   }
824                   numberEntries++;
825                   if (iPass) {
826                     // slack
827                     cliqueType[numberCliques] = 'S';
828                     cliqueStart[numberCliques + 1] = numberEntries;
829                   }
830                   numberCliques++;
831                 }
832               } else {
833                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
834                   if (printit && !iPass)
835                     printf("%d to one implies %d to zero and %d to zero implies %d to one\n",
836                       iColumn, jColumn, jColumn, iColumn);
837                   nFix++; // iColumn+jColumn=1
838                 } else {
839                   if (printit && !iPass)
840                     printf("%d to one implies %d to zero and %d to zero implies %d to zero\n",
841                       iColumn, jColumn, jColumn, iColumn);
842                   nFix++; // iColumn is 0
843                 }
844               }
845             }
846           }
847           for (k = toOne_[jColumn]; k < toZero_[jColumn + 1]; k++) {
848             if (sequenceInCliqueEntry(fixEntry_[k]) == (iColumn)) {
849               if (oneFixesInCliqueEntry(fixEntry_[j])) {
850                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
851                   if (printit && !iPass)
852                     printf("%d to one implies %d to one and %d to one implies %d to one\n",
853                       iColumn, jColumn, jColumn, iColumn);
854                   nFix++; // iColumn == jColumn
855                 } else {
856                   if (printit && !iPass)
857                     printf("%d to one implies %d to one and %d to one implies %d to zero\n",
858                       iColumn, jColumn, jColumn, iColumn);
859                   nFix++; // iColumn is 0
860                 }
861               } else {
862                 if (oneFixesInCliqueEntry(fixEntry_[k])) {
863                   if (printit && !iPass)
864                     printf("%d to one implies %d to zero and %d to one implies %d to one\n",
865                       iColumn, jColumn, jColumn, iColumn);
866                   nFix++; // jColumn is 0
867                 } else {
868                   if (printit && !iPass)
869                     printf("%d to one implies %d to zero and %d to one implies %d to zero\n",
870                       iColumn, jColumn, jColumn, iColumn);
871                   // 1-1 illegal
872                   if (iPass) {
873                     CliqueEntry temp;
874                     setOneFixesInCliqueEntry(temp, true);
875                     setSequenceInCliqueEntry(temp, iColumn);
876                     entry[numberEntries] = temp;
877                   }
878                   numberEntries++;
879                   if (iPass) {
880                     CliqueEntry temp;
881                     setOneFixesInCliqueEntry(temp, true);
882                     setSequenceInCliqueEntry(temp, jColumn);
883                     entry[numberEntries] = temp;
884                   }
885                   numberEntries++;
886                   if (iPass) {
887                     // slack
888                     cliqueType[numberCliques] = 'S';
889                     cliqueStart[numberCliques + 1] = numberEntries;
890                   }
891                   numberCliques++;
892                 }
893               }
894             }
895           }
896         }
897       }
898     }
899     if (!iPass)
900       printf("%d cliques and %d fixed (%d already from matrix))\n",
901         numberCliques - numberMatrixCliques, nFix, numberMatrixCliques);
902   }
903   int iClique;
904   outDupsEtc(numberIntegers_, numberCliques, numberMatrixCliques,
905     cliqueStart, cliqueType, entry, -1, printit ? 2 : 1);
906   printf("%d matrix cliques and %d found by probing\n", numberMatrixCliques, numberCliques - numberMatrixCliques);
907   int *zeroStart = new int[numberIntegers_ + 1];
908   int *oneStart = new int[numberIntegers_];
909   int *zeroCount = new int[numberIntegers_];
910   int *oneCount = new int[numberIntegers_];
911   char *mark = new char[numberIntegers_];
912   memset(mark, 0, numberIntegers_);
913   int nStrengthen = -1;
914   int iPass = 0;
915   while (nStrengthen && iPass < 20) {
916     iPass++;
917     int numberLastTime = numberCliques;
918     int *count = new int[numberCliques];
919     int i, iColumn;
920     for (i = 0; i < numberCliques; i++) {
921       count[i] = 0;
922     }
923     int *whichCount = new int[numberCliques];
924     CoinZeroN(zeroCount, numberIntegers_);
925     CoinZeroN(oneCount, numberIntegers_);
926     for (iClique = 0; iClique < numberCliques; iClique++) {
927       for (CoinBigIndex j = cliqueStart[iClique]; j < cliqueStart[iClique + 1]; j++) {
928         int iColumn = static_cast< int >(sequenceInCliqueEntry(entry[j]));
929         if (oneFixesInCliqueEntry(entry[j])) {
930           oneCount[iColumn]++;
931         } else {
932           zeroCount[iColumn]++;
933         }
934       }
935     }
936     int j;
937     zeroStart[0] = 0;
938     cliqueStart[0] = 0;
939     for (j = 0; j < numberIntegers_; j++) {
940       int n;
941       n = zeroCount[j];
942       zeroCount[j] = 0;
943       oneStart[j] = zeroStart[j] + n;
944       n = oneCount[j];
945       oneCount[j] = 0;
946       zeroStart[j + 1] = oneStart[j] + n;
947     }
948     for (iClique = 0; iClique < numberCliques; iClique++) {
949       for (CoinBigIndex j = cliqueStart[iClique]; j < cliqueStart[iClique + 1]; j++) {
950         int iColumn = static_cast< int >(sequenceInCliqueEntry(entry[j]));
951         if (oneFixesInCliqueEntry(entry[j])) {
952           int k = oneCount[iColumn];
953           oneCount[iColumn]++;
954           int put = oneStart[iColumn] + k;
955           whichClique[put] = iClique;
956         } else {
957           int k = zeroCount[iColumn];
958           zeroCount[iColumn]++;
959           int put = zeroStart[iColumn] + k;
960           whichClique[put] = iClique;
961         }
962       }
963     }
964     nStrengthen = 0;
965     CoinBigIndex numberEntries = cliqueStart[numberCliques];
966     CoinBigIndex maximumEntries = numberEntries;
967     int maximumCliques = numberCliques;
968     for (iColumn = 0; iColumn < numberIntegers_; iColumn++) {
969       int i;
970       int n;
971       int nCount = 0;
972       n = 0;
973       for (i = zeroStart[iColumn]; i < oneStart[iColumn]; i++) {
974         int jClique = whichClique[i];
975         //if (jClique<numberMatrixCliques)
976         //continue;
977         CoinBigIndex j = cliqueStart[jClique];
978         //assert (cliqueStart[jClique+1]==j+2);
979         for (; j < cliqueStart[jClique + 1]; j++) {
980           CliqueEntry eJ = entry[j];
981           int jColumn = sequenceInCliqueEntry(eJ);
982           if (jColumn > iColumn && !mark[jColumn]) {
983             mark[jColumn] = 1;
984             whichP[n++] = jColumn;
985             assert(n < numberIntegers_);
986             if (oneFixesInCliqueEntry(eJ)) {
987               for (int k = oneStart[jColumn]; k < zeroStart[jColumn + 1]; k++) {
988                 int jClique = whichClique[k];
989                 if (!count[jClique])
990                   whichCount[nCount++] = jClique;
991                 count[jClique]++;
992               }
993             } else {
994               for (int k = zeroStart[jColumn]; k < oneStart[jColumn]; k++) {
995                 int jClique = whichClique[k];
996                 if (!count[jClique])
997                   whichCount[nCount++] = jClique;
998                 count[jClique]++;
999               }
1000             }
1001           }
1002         }
1003       }
1004       std::sort(whichP, whichP + n);
1005       for (i = 0; i < nCount; i++) {
1006         int jClique = whichCount[i];
1007         int jCount = count[jClique];
1008         count[jClique] = 0;
1009         if (jCount == cliqueStart[jClique + 1] - cliqueStart[jClique]) {
1010           printf("Zero can extend %d [ ", jClique);
1011           for (CoinBigIndex i = cliqueStart[jClique]; i < cliqueStart[jClique + 1]; i++)
1012             printf("%d(%d) ", sequenceInCliqueEntry(entry[i]), oneFixesInCliqueEntry(entry[i]));
1013           printf("] by %d(0)\n", iColumn);
1014           nStrengthen++;
1015           if (numberEntries + jCount + 1 > maximumEntries) {
1016             maximumEntries = CoinMax(numberEntries + jCount + 1, (maximumEntries * 12) / 10 + 100);
1017             CliqueEntry *temp = new CliqueEntry[maximumEntries];
1018             memcpy(temp, entry, numberEntries * sizeof(CliqueEntry));
1019             delete[] entry;
1020             entry = temp;
1021             int *tempI = new int[maximumEntries];
1022             memcpy(tempI, whichClique, numberEntries * sizeof(int));
1023             delete[] whichClique;
1024             whichClique = tempI;
1025           }
1026           if (numberCliques == maximumCliques) {
1027             maximumCliques = (maximumCliques * 12) / 10 + 100;
1028             CoinBigIndex *temp = new CoinBigIndex[maximumCliques + 1];
1029             memcpy(temp, cliqueStart, (numberCliques + 1) * sizeof(CoinBigIndex));
1030             delete[] cliqueStart;
1031             cliqueStart = temp;
1032             char *tempT = new char[maximumCliques];
1033             memcpy(tempT, cliqueType, numberCliques);
1034             delete[] cliqueType;
1035             cliqueType = tempT;
1036           }
1037           CliqueEntry eI;
1038           eI.fixes = 0;
1039           setSequenceInCliqueEntry(eI, iColumn);
1040           setOneFixesInCliqueEntry(eI, false);
1041           entry[numberEntries++] = eI;
1042           whichM[0] = iColumn;
1043           int n = 1;
1044           for (CoinBigIndex i = cliqueStart[jClique]; i < cliqueStart[jClique + 1]; i++) {
1045             entry[numberEntries++] = entry[i];
1046             whichM[n++] = sequenceInCliqueEntry(entry[i]);
1047           }
1048           CoinSort_2(whichM, whichM + n, (reinterpret_cast< int * >(entry)) + numberEntries - n);
1049           cliqueType[numberCliques] = 'S';
1050           numberCliques++;
1051           cliqueStart[numberCliques] = numberEntries;
1052         }
1053       }
1054       for (i = 0; i < n; i++)
1055         mark[whichP[i]] = 0;
1056       nCount = 0;
1057       n = 0;
1058       for (i = oneStart[iColumn]; i < zeroStart[iColumn + 1]; i++) {
1059         int jClique = whichClique[i];
1060         //if (jClique<numberMatrixCliques)
1061         //continue;
1062         CoinBigIndex j = cliqueStart[jClique];
1063         //assert (cliqueStart[jClique+1]==j+2);
1064         for (; j < cliqueStart[jClique + 1]; j++) {
1065           CliqueEntry eJ = entry[j];
1066           int jColumn = sequenceInCliqueEntry(eJ);
1067           if (jColumn > iColumn && !mark[jColumn]) {
1068             mark[jColumn] = 1;
1069             whichP[n++] = jColumn;
1070             assert(n < numberIntegers_);
1071             if (oneFixesInCliqueEntry(eJ)) {
1072               for (int k = oneStart[jColumn]; k < zeroStart[jColumn + 1]; k++) {
1073                 int jClique = whichClique[k];
1074                 if (!count[jClique])
1075                   whichCount[nCount++] = jClique;
1076                 count[jClique]++;
1077               }
1078             } else {
1079               for (int k = zeroStart[jColumn]; k < oneStart[jColumn]; k++) {
1080                 int jClique = whichClique[k];
1081                 if (!count[jClique])
1082                   whichCount[nCount++] = jClique;
1083                 count[jClique]++;
1084               }
1085             }
1086           }
1087         }
1088       }
1089       std::sort(whichP, whichP + n);
1090       for (i = 0; i < nCount; i++) {
1091         int jClique = whichCount[i];
1092         int jCount = count[jClique];
1093         count[jClique] = 0;
1094         if (jCount == cliqueStart[jClique + 1] - cliqueStart[jClique]) {
1095 #if 1
1096           if (printit) {
1097             printf("One can extend %d [ ", jClique);
1098             for (CoinBigIndex i = cliqueStart[jClique]; i < cliqueStart[jClique + 1]; i++)
1099               printf("%d(%d) ", sequenceInCliqueEntry(entry[i]), oneFixesInCliqueEntry(entry[i]));
1100             printf("] by %d(1)\n", iColumn);
1101           }
1102 #endif
1103           nStrengthen++;
1104           if (numberEntries + jCount + 1 > maximumEntries) {
1105             maximumEntries = CoinMax(numberEntries + jCount + 1, (maximumEntries * 12) / 10 + 100);
1106             CliqueEntry *temp = new CliqueEntry[maximumEntries];
1107             memcpy(temp, entry, numberEntries * sizeof(CliqueEntry));
1108             delete[] entry;
1109             entry = temp;
1110             int *tempI = new int[maximumEntries];
1111             memcpy(tempI, whichClique, numberEntries * sizeof(int));
1112             delete[] whichClique;
1113             whichClique = tempI;
1114           }
1115           if (numberCliques == maximumCliques) {
1116             maximumCliques = (maximumCliques * 12) / 10 + 100;
1117             CoinBigIndex *temp = new CoinBigIndex[maximumCliques + 1];
1118             memcpy(temp, cliqueStart, (numberCliques + 1) * sizeof(CoinBigIndex));
1119             delete[] cliqueStart;
1120             cliqueStart = temp;
1121             char *tempT = new char[maximumCliques];
1122             memcpy(tempT, cliqueType, numberCliques);
1123             delete[] cliqueType;
1124             cliqueType = tempT;
1125           }
1126           CliqueEntry eI;
1127           eI.fixes = 0;
1128           setSequenceInCliqueEntry(eI, iColumn);
1129           setOneFixesInCliqueEntry(eI, true);
1130           entry[numberEntries++] = eI;
1131           whichM[0] = iColumn;
1132           int n = 1;
1133           for (CoinBigIndex i = cliqueStart[jClique]; i < cliqueStart[jClique + 1]; i++) {
1134             entry[numberEntries++] = entry[i];
1135             whichM[n++] = sequenceInCliqueEntry(entry[i]);
1136           }
1137           CoinSort_2(whichM, whichM + n, (reinterpret_cast< int * >(entry)) + numberEntries - n);
1138           cliqueType[numberCliques] = 'S';
1139           numberCliques++;
1140           cliqueStart[numberCliques] = numberEntries;
1141         }
1142       }
1143       for (i = 0; i < n; i++)
1144         mark[whichP[i]] = 0;
1145     }
1146     if (nStrengthen) {
1147       int numberDeleted = outDupsEtc(numberIntegers_, numberCliques, numberMatrixCliques,
1148         cliqueStart, cliqueType, entry, numberLastTime, printit ? 2 : 1);
1149       if (numberDeleted < 0 || (iPass > 1 && numberCliques - numberDeleted > 5000))
1150         nStrengthen = 0;
1151     }
1152     delete[] count;
1153     delete[] whichCount;
1154   }
1155 #if 0
1156   if (numberCliques>numberMatrixCliques) {
1157     // should keep as cliques and also use in branching ??
1158     double * element = new double [numberIntegers_];
1159     for (iClique=numberMatrixCliques;iClique<numberCliques;iClique++) {
1160       int n=0;
1161       double rhs=1.0;
1162       for (int i=cliqueStart[iClique];i<cliqueStart[iClique+1];i++) {
1163 	CliqueEntry eI=entry[i];
1164 	int iColumn = integerVariable_[sequenceInCliqueEntry(eI)];
1165 	whichP[n]=iColumn;
1166 	if (oneFixesInCliqueEntry(eI)) {
1167 	  element[n++]=1.0;
1168 	} else {
1169 	  element[n++]=-1.0;
1170 	  rhs -= 1.0;
1171 	}
1172       }
1173       stored->addCut(-COIN_DBL_MAX,rhs,n,whichP,element);
1174     }
1175     delete [] element;
1176   }
1177 #endif
1178   OsiSolverInterface *newSolver = NULL;
1179   if (numberCliques > numberMatrixCliques || alwaysDo) {
1180     newSolver = si.clone();
1181     // Delete all rows
1182     CoinBigIndex *start = new CoinBigIndex[CoinMax(numberRows, numberCliques + 1)];
1183     int i;
1184     int *start2 = reinterpret_cast< int * >(start);
1185     for (i = 0; i < numberRows; i++)
1186       start2[i] = i;
1187     newSolver->deleteRows(numberRows, start2);
1188     start[0] = 0;
1189     CoinBigIndex numberElements = cliqueStart[numberCliques];
1190     int *column = new int[numberElements];
1191     double *element = new double[numberElements];
1192     double *lower = new double[numberCliques];
1193     double *upper = new double[numberCliques];
1194     numberElements = 0;
1195     for (iClique = 0; iClique < numberCliques; iClique++) {
1196       double rhs = 1.0;
1197       for (CoinBigIndex i = cliqueStart[iClique]; i < cliqueStart[iClique + 1]; i++) {
1198         CliqueEntry eI = entry[i];
1199         int iColumn = integerVariable_[sequenceInCliqueEntry(eI)];
1200         column[numberElements] = iColumn;
1201         if (oneFixesInCliqueEntry(eI)) {
1202           element[numberElements++] = 1.0;
1203         } else {
1204           element[numberElements++] = -1.0;
1205           rhs -= 1.0;
1206         }
1207       }
1208       start[iClique + 1] = numberElements;
1209       assert(cliqueType[iClique] == 'S' || cliqueType[iClique] == 'E');
1210       if (cliqueType[iClique] == 'S')
1211         lower[iClique] = -COIN_DBL_MAX;
1212       else
1213         lower[iClique] = rhs;
1214       upper[iClique] = rhs;
1215     }
1216     newSolver->addRows(numberCliques, start, column, element, lower, upper);
1217     delete[] start;
1218     delete[] column;
1219     delete[] element;
1220     delete[] lower;
1221     delete[] upper;
1222   }
1223   delete[] mark;
1224   delete[] whichP;
1225   delete[] whichM;
1226   delete[] cliqueStart;
1227   delete[] entry;
1228   delete[] cliqueType;
1229   delete[] zeroStart;
1230   delete[] oneStart;
1231   delete[] zeroCount;
1232   delete[] oneCount;
1233   delete[] whichClique;
1234   return newSolver;
1235 }
1236 // Take action if cut generator can fix a variable (toValue -1 for down, +1 for up)
fixes(int variable,int toValue,int fixedVariable,bool fixedToLower)1237 bool CglTreeProbingInfo::fixes(int variable, int toValue, int fixedVariable, bool fixedToLower)
1238 {
1239   //printf("%d going to %d fixes %d at %g\n",variable,toValue,fixedVariable,fixedToValue);
1240   int intVariable = backward_[variable];
1241   if (intVariable < 0) // off as no longer in order FIX
1242     return true; // not 0-1 (well wasn't when constructor was called)
1243   int intFix = backward_[fixedVariable];
1244   if (intFix < 0)
1245     intFix = numberIntegers_ + fixedVariable; // not 0-1
1246   int fixedTo = fixedToLower ? 0 : 1;
1247   if (numberEntries_ == maximumEntries_) {
1248     // See if taking too much memory
1249     if (maximumEntries_ >= CoinMax(1000000, 10 * numberIntegers_))
1250       return false;
1251     maximumEntries_ += 100 + maximumEntries_ / 2;
1252     CliqueEntry *temp1 = new CliqueEntry[maximumEntries_];
1253     memcpy(temp1, fixEntry_, numberEntries_ * sizeof(CliqueEntry));
1254     delete[] fixEntry_;
1255     fixEntry_ = temp1;
1256     int *temp2 = new int[maximumEntries_];
1257     memcpy(temp2, fixingEntry_, numberEntries_ * sizeof(int));
1258     delete[] fixingEntry_;
1259     fixingEntry_ = temp2;
1260   }
1261   CliqueEntry entry1;
1262   entry1.fixes = 0;
1263   setOneFixesInCliqueEntry(entry1, fixedTo != 0);
1264   setSequenceInCliqueEntry(entry1, intFix);
1265   fixEntry_[numberEntries_] = entry1;
1266   assert(toValue == -1 || toValue == 1);
1267   assert(fixedTo == 0 || fixedTo == 1);
1268   if (toValue < 0)
1269     fixingEntry_[numberEntries_++] = intVariable << 1;
1270   else
1271     fixingEntry_[numberEntries_++] = (intVariable << 1) | 1;
1272   return true;
1273 }
1274 // Initalizes fixing arrays etc - returns true if we want to save info
initializeFixing(const OsiSolverInterface * model)1275 int CglTreeProbingInfo::initializeFixing(const OsiSolverInterface *model)
1276 {
1277   if (numberEntries_ >= 0)
1278     return 2; // already got arrays
1279   else if (numberEntries_ == -2)
1280     return numberEntries_;
1281   delete[] fixEntry_;
1282   delete[] toZero_;
1283   delete[] toOne_;
1284   delete[] integerVariable_;
1285   delete[] backward_;
1286   delete[] fixingEntry_;
1287   numberVariables_ = model->getNumCols();
1288   // Too many ... but
1289   integerVariable_ = new int[numberVariables_];
1290   backward_ = new int[numberVariables_];
1291   numberIntegers_ = 0;
1292   int i;
1293   // Get integer types
1294   const char *columnType = model->getColType(true);
1295   for (i = 0; i < numberVariables_; i++) {
1296     backward_[i] = -1;
1297     if (columnType[i]) {
1298       if (columnType[i] == 1) {
1299         backward_[i] = numberIntegers_;
1300         integerVariable_[numberIntegers_++] = i;
1301       } else {
1302         backward_[i] = -2;
1303       }
1304     }
1305   }
1306   toZero_ = NULL;
1307   toOne_ = NULL;
1308   fixEntry_ = NULL;
1309   fixingEntry_ = NULL;
1310   maximumEntries_ = 0;
1311   numberEntries_ = 0;
1312   return 1;
1313 }
1314 // Converts to ordered and takes out duplicates
convert()1315 void CglTreeProbingInfo::convert()
1316 {
1317   if (numberEntries_ >= 0) {
1318     CoinSort_2(fixingEntry_, fixingEntry_ + numberEntries_, fixEntry_);
1319     assert(!toZero_);
1320     toZero_ = new int[numberIntegers_ + 1];
1321     toOne_ = new int[numberIntegers_];
1322     toZero_[0] = 0;
1323     int n = 0;
1324     int put = 0;
1325     for (int intVariable = 0; intVariable < numberIntegers_; intVariable++) {
1326       int last = n;
1327       for (; n < numberEntries_; n++) {
1328         int value = fixingEntry_[n];
1329         int iVar = value >> 1;
1330         int way = value & 1;
1331         if (intVariable != iVar || way)
1332           break;
1333       }
1334       if (n > last) {
1335         // sort
1336         assert(sizeof(int) == 4);
1337         std::sort(reinterpret_cast< unsigned int * >(fixEntry_) + last,
1338           reinterpret_cast< unsigned int * >(fixEntry_) + n);
1339         CliqueEntry temp2;
1340         temp2.fixes = 0;
1341         setSequenceInCliqueEntry(temp2, numberVariables_ + 1);
1342         for (int i = last; i < n; i++) {
1343           if (sequenceInCliqueEntry(temp2) != sequenceInCliqueEntry(fixEntry_[i]) || oneFixesInCliqueEntry(temp2) || oneFixesInCliqueEntry(fixEntry_[i])) {
1344             temp2 = fixEntry_[i];
1345             fixEntry_[put++] = temp2;
1346           }
1347         }
1348       }
1349       toOne_[intVariable] = put;
1350       last = n;
1351       for (; n < numberEntries_; n++) {
1352         int value = fixingEntry_[n];
1353         int iVar = value >> 1;
1354         if (intVariable != iVar)
1355           break;
1356       }
1357       if (n > last) {
1358         // sort
1359         assert(sizeof(int) == 4);
1360         std::sort(reinterpret_cast< unsigned int * >(fixEntry_) + last,
1361           reinterpret_cast< unsigned int * >(fixEntry_) + n);
1362         CliqueEntry temp2;
1363         temp2.fixes = 0;
1364         setSequenceInCliqueEntry(temp2, numberVariables_ + 1);
1365         for (int i = last; i < n; i++) {
1366           if (sequenceInCliqueEntry(temp2) != sequenceInCliqueEntry(fixEntry_[i]) || oneFixesInCliqueEntry(temp2) || oneFixesInCliqueEntry(fixEntry_[i])) {
1367             temp2 = fixEntry_[i];
1368             fixEntry_[put++] = temp2;
1369           }
1370         }
1371         last = n;
1372       }
1373       toZero_[intVariable + 1] = put;
1374     }
1375     delete[] fixingEntry_;
1376     fixingEntry_ = NULL;
1377     numberEntries_ = -2;
1378   }
1379 }
1380 // Fix entries in a solver using implications
fixColumns(OsiSolverInterface & si) const1381 int CglTreeProbingInfo::fixColumns(OsiSolverInterface &si) const
1382 {
1383   int nFix = 0;
1384   const double *lower = si.getColLower();
1385   const double *upper = si.getColUpper();
1386   bool feasible = true;
1387   for (int jColumn = 0; jColumn < static_cast< int >(numberIntegers_); jColumn++) {
1388     int iColumn = integerVariable_[jColumn];
1389     if (upper[iColumn] == 0.0) {
1390       int j;
1391       for (j = toZero_[jColumn]; j < toOne_[jColumn]; j++) {
1392         int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1393         kColumn = integerVariable_[kColumn];
1394         bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1395         if (fixToOne) {
1396           if (lower[kColumn] == 0.0) {
1397             if (upper[kColumn] == 1.0) {
1398               si.setColLower(kColumn, 1.0);
1399               nFix++;
1400             } else {
1401               // infeasible!
1402               feasible = false;
1403             }
1404           }
1405         } else {
1406           if (upper[kColumn] == 1.0) {
1407             if (lower[kColumn] == 0.0) {
1408               si.setColUpper(kColumn, 0.0);
1409               nFix++;
1410             } else {
1411               // infeasible!
1412               feasible = false;
1413             }
1414           }
1415         }
1416       }
1417     } else if (lower[iColumn] == 1.0) {
1418       int j;
1419       for (j = toOne_[jColumn]; j < toZero_[jColumn + 1]; j++) {
1420         int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1421         kColumn = integerVariable_[kColumn];
1422         bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1423         if (fixToOne) {
1424           if (lower[kColumn] == 0.0) {
1425             if (upper[kColumn] == 1.0) {
1426               si.setColLower(kColumn, 1.0);
1427               nFix++;
1428             } else {
1429               // infeasible!
1430               feasible = false;
1431             }
1432           }
1433         } else {
1434           if (upper[kColumn] == 1.0) {
1435             if (lower[kColumn] == 0.0) {
1436               si.setColUpper(kColumn, 0.0);
1437               nFix++;
1438             } else {
1439               // infeasible!
1440               feasible = false;
1441             }
1442           }
1443         }
1444       }
1445     }
1446   }
1447   if (!feasible) {
1448 #ifdef COIN_DEVELOP
1449     printf("treeprobing says infeasible!\n");
1450 #endif
1451     nFix = -1;
1452   }
1453   return nFix;
1454 }
1455 // Fix entries in a solver using implications for one variable
fixColumns(int iColumn,int value,OsiSolverInterface & si) const1456 int CglTreeProbingInfo::fixColumns(int iColumn, int value, OsiSolverInterface &si) const
1457 {
1458   assert(value == 0 || value == 1);
1459   int nFix = 0;
1460   const double *lower = si.getColLower();
1461   const double *upper = si.getColUpper();
1462   bool feasible = true;
1463   int jColumn = backward_[iColumn];
1464   if (jColumn < 0 || !toZero_)
1465     return 0;
1466   if (!value) {
1467     int j;
1468     for (j = toZero_[jColumn]; j < toOne_[jColumn]; j++) {
1469       int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1470       kColumn = integerVariable_[kColumn];
1471       bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1472       if (fixToOne) {
1473         if (lower[kColumn] == 0.0) {
1474           if (upper[kColumn] == 1.0) {
1475             si.setColLower(kColumn, 1.0);
1476             nFix++;
1477           } else {
1478             // infeasible!
1479             feasible = false;
1480           }
1481         }
1482       } else {
1483         if (upper[kColumn] == 1.0) {
1484           if (lower[kColumn] == 0.0) {
1485             si.setColUpper(kColumn, 0.0);
1486             nFix++;
1487           } else {
1488             // infeasible!
1489             feasible = false;
1490           }
1491         }
1492       }
1493     }
1494   } else {
1495     int j;
1496     for (j = toOne_[jColumn]; j < toZero_[jColumn + 1]; j++) {
1497       int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1498       kColumn = integerVariable_[kColumn];
1499       bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1500       if (fixToOne) {
1501         if (lower[kColumn] == 0.0) {
1502           if (upper[kColumn] == 1.0) {
1503             si.setColLower(kColumn, 1.0);
1504             nFix++;
1505           } else {
1506             // infeasible!
1507             feasible = false;
1508           }
1509         }
1510       } else {
1511         if (upper[kColumn] == 1.0) {
1512           if (lower[kColumn] == 0.0) {
1513             si.setColUpper(kColumn, 0.0);
1514             nFix++;
1515           } else {
1516             // infeasible!
1517             feasible = false;
1518           }
1519         }
1520       }
1521     }
1522   }
1523   if (!feasible) {
1524 #ifdef COIN_DEVELOP
1525     printf("treeprobing says infeasible!\n");
1526 #endif
1527     nFix = -1;
1528   }
1529   return nFix;
1530 }
1531 // Packs down entries
packDown()1532 int CglTreeProbingInfo::packDown()
1533 {
1534   convert();
1535   int iPut = 0;
1536   int iLast = 0;
1537   for (int jColumn = 0; jColumn < static_cast< int >(numberIntegers_); jColumn++) {
1538     int j;
1539     for (j = iLast; j < toOne_[jColumn]; j++) {
1540       int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1541       if (kColumn < numberIntegers_)
1542         fixEntry_[iPut++] = fixEntry_[j];
1543     }
1544     iLast = toOne_[jColumn];
1545     toOne_[jColumn] = iPut;
1546     for (j = iLast; j < toZero_[jColumn + 1]; j++) {
1547       int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1548       if (kColumn < numberIntegers_)
1549         fixEntry_[iPut++] = fixEntry_[j];
1550     }
1551     iLast = toZero_[jColumn + 1];
1552     toZero_[jColumn + 1] = iPut;
1553   }
1554   return iPut;
1555 }
generateCuts(const OsiSolverInterface & si,OsiCuts & cs,const CglTreeInfo) const1556 void CglTreeProbingInfo::generateCuts(const OsiSolverInterface &si, OsiCuts &cs,
1557   const CglTreeInfo /*info*/) const
1558 {
1559   const double *lower = si.getColLower();
1560   const double *upper = si.getColUpper();
1561   const double *colsol = si.getColSolution();
1562   CoinPackedVector lbs(false);
1563   CoinPackedVector ubs(false);
1564   int numberFixed = 0;
1565   char *fixed = NULL;
1566   for (int jColumn = 0; jColumn < static_cast< int >(numberIntegers_); jColumn++) {
1567     int iColumn = integerVariable_[jColumn];
1568     assert(iColumn >= 0 && iColumn < si.getNumCols());
1569     if (lower[iColumn] == 0.0 && upper[iColumn] == 1.0) {
1570       double value1 = colsol[iColumn];
1571       int j;
1572       for (j = toZero_[jColumn]; j < toOne_[jColumn]; j++) {
1573         int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1574         kColumn = integerVariable_[kColumn];
1575         assert(kColumn >= 0 && kColumn < si.getNumCols());
1576         assert(kColumn != iColumn);
1577         if (lower[kColumn] == 0.0 && upper[kColumn] == 1.0) {
1578           double value2 = colsol[kColumn];
1579           bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1580           if (fixToOne) {
1581             if (value1 + value2 < 0.99999) {
1582               OsiRowCut rc;
1583               int index[2];
1584               double element[2];
1585               index[0] = iColumn;
1586               element[0] = 1.0;
1587               index[1] = kColumn;
1588               element[1] = 1.0;
1589               rc.setLb(1.0);
1590               rc.setUb(COIN_DBL_MAX);
1591               rc.setRow(2, index, element, false);
1592               cs.insertIfNotDuplicate(rc);
1593             }
1594           } else {
1595             if (value1 - value2 < -0.00001) {
1596               OsiRowCut rc;
1597               int index[2];
1598               double element[2];
1599               index[0] = iColumn;
1600               element[0] = 1.0;
1601               index[1] = kColumn;
1602               element[1] = -1.0;
1603               rc.setLb(0.0);
1604               rc.setUb(COIN_DBL_MAX);
1605               rc.setRow(2, index, element, false);
1606               cs.insertIfNotDuplicate(rc);
1607             }
1608           }
1609         }
1610       }
1611       for (j = toOne_[jColumn]; j < toZero_[jColumn + 1]; j++) {
1612         int kColumn = sequenceInCliqueEntry(fixEntry_[j]);
1613         kColumn = integerVariable_[kColumn];
1614         assert(kColumn >= 0 && kColumn < si.getNumCols());
1615         assert(kColumn != iColumn);
1616         if (lower[kColumn] == 0.0 && upper[kColumn] == 1.0) {
1617           double value2 = colsol[kColumn];
1618           bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1619           if (fixToOne) {
1620             if (value1 - value2 > 0.00001) {
1621               OsiRowCut rc;
1622               int index[2];
1623               double element[2];
1624               index[0] = iColumn;
1625               element[0] = 1.0;
1626               index[1] = kColumn;
1627               element[1] = -1.0;
1628               rc.setLb(-COIN_DBL_MAX);
1629               rc.setUb(0.0);
1630               rc.setRow(2, index, element, false);
1631               cs.insertIfNotDuplicate(rc);
1632             }
1633           } else {
1634             if (value1 + value2 > 1.00001) {
1635               OsiRowCut rc;
1636               int index[2];
1637               double element[2];
1638               index[0] = iColumn;
1639               element[0] = 1.0;
1640               index[1] = kColumn;
1641               element[1] = 1.0;
1642               rc.setLb(-COIN_DBL_MAX);
1643               rc.setUb(1.0);
1644               rc.setRow(2, index, element, false);
1645               cs.insertIfNotDuplicate(rc);
1646             }
1647           }
1648         }
1649       }
1650     } else if (upper[iColumn] == 0.0) {
1651       for (int j = toZero_[jColumn]; j < toOne_[jColumn]; j++) {
1652         int kColumn01 = sequenceInCliqueEntry(fixEntry_[j]);
1653         int kColumn = integerVariable_[kColumn01];
1654         assert(kColumn >= 0 && kColumn < si.getNumCols());
1655         bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1656         if (lower[kColumn] == 0.0 && upper[kColumn] == 1.0) {
1657           if (!fixed) {
1658             fixed = new char[numberIntegers_];
1659             memset(fixed, 0, numberIntegers_);
1660           }
1661           if (fixToOne) {
1662             if ((fixed[kColumn01] & 1) == 0) {
1663               fixed[kColumn01] |= 1;
1664               lbs.insert(kColumn, 1.0);
1665             }
1666           } else {
1667             if ((fixed[kColumn01] & 2) == 0) {
1668               fixed[kColumn01] |= 2;
1669               ubs.insert(kColumn, 0.0);
1670             }
1671           }
1672           numberFixed++;
1673         } else if ((fixToOne && upper[kColumn] == 0.0) || (!fixToOne && lower[kColumn] == 1.0)) {
1674           // infeasible!
1675           // generate infeasible cut and return
1676           OsiRowCut rc;
1677           rc.setLb(COIN_DBL_MAX);
1678           rc.setUb(0.0);
1679           cs.insertIfNotDuplicate(rc);
1680           //printf("IMPINFEAS!\n");
1681           return;
1682         }
1683       }
1684     } else {
1685       for (int j = toOne_[jColumn]; j < toZero_[jColumn + 1]; j++) {
1686         int kColumn01 = sequenceInCliqueEntry(fixEntry_[j]);
1687         int kColumn = integerVariable_[kColumn01];
1688         assert(kColumn >= 0 && kColumn < si.getNumCols());
1689         bool fixToOne = oneFixesInCliqueEntry(fixEntry_[j]);
1690         if (lower[kColumn] == 0.0 && upper[kColumn] == 1.0) {
1691           if (!fixed) {
1692             fixed = new char[numberIntegers_];
1693             memset(fixed, 0, numberIntegers_);
1694           }
1695           if (fixToOne) {
1696             if ((fixed[kColumn01] & 1) == 0) {
1697               fixed[kColumn01] |= 1;
1698               lbs.insert(kColumn, 1.0);
1699             }
1700           } else {
1701             if ((fixed[kColumn01] & 2) == 0) {
1702               fixed[kColumn01] |= 2;
1703               ubs.insert(kColumn, 0.0);
1704             }
1705           }
1706           numberFixed++;
1707         } else if ((fixToOne && upper[kColumn] == 0.0) || (!fixToOne && lower[kColumn] == 1.0)) {
1708           // infeasible!
1709           // generate infeasible cut and return
1710           OsiRowCut rc;
1711           rc.setLb(COIN_DBL_MAX);
1712           rc.setUb(0.0);
1713           cs.insertIfNotDuplicate(rc);
1714           //printf("IMPINFEAS!\n");
1715           return;
1716         }
1717       }
1718     }
1719   }
1720   if (numberFixed) {
1721     int feasible = true;
1722     for (int i = 0; i < numberIntegers_; i++) {
1723       if (fixed[i] == 3) {
1724         feasible = false;
1725         break;
1726       }
1727     }
1728     delete[] fixed;
1729     fixed = NULL;
1730     if (feasible) {
1731       //printf("IMP fixed %d\n",numberFixed);
1732       OsiColCut cc;
1733       cc.setUbs(ubs);
1734       cc.setLbs(lbs);
1735       cc.setEffectiveness(1.0);
1736       cs.insert(cc);
1737     } else {
1738       // infeasible!
1739       // generate infeasible cut
1740       OsiRowCut rc;
1741       rc.setLb(COIN_DBL_MAX);
1742       rc.setUb(0.0);
1743       cs.insertIfNotDuplicate(rc);
1744       //printf("IMPINFEAS!\n");
1745     }
1746   }
1747 
1748   if (fixed)
1749     delete[] fixed;
1750 }
1751 
1752 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
1753 */
1754