source: trunk/Cbc/src/CbcSymmetry.cpp @ 2080

Last change on this file since 2080 was 2051, checked in by forrest, 5 years ago

fix stupid bug in orbital branching

File size: 36.5 KB
Line 
1/* $Id: CbcSymmetry.cpp 925 2012-11-27 19:11:04Z stefan $
2 *
3 * hacked from CouenneSymmetry.cpp
4 * Name:    CbcSymmetry.cpp
5 * Author:  Jim Ostrowski (the good bits - rest JJHF)
6 * Purpose: methods for exploiting symmetry
7 * Date:    October 13, 2010
8 *
9 * This file is licensed under the Eclipse Public License (EPL)
10 */
11//#define PRINT_MORE 1
12#include <stdio.h>
13
14#ifdef COIN_HAS_NTY
15
16#include <cassert>
17#include <vector>
18#include <algorithm>
19#include <ostream>
20#include <iterator>
21
22#include "CbcSymmetry.hpp"
23#include "CbcBranchingObject.hpp"
24#include "CoinTime.hpp"
25/* Deliberately not threadsafe to save effort
26   Just for statistics
27   and not worth gathering across threads
28   can redo later
29 */
30static int nautyBranchCalls_ = 0;
31static int nautyBranchSucceeded_ = 0;
32static int nautyFixCalls_ = 0;
33static int nautyFixSucceeded_ = 0;
34static double nautyTime_ = 0.0;
35static double nautyFixes_= 0.0; 
36static double nautyOtherBranches_ = 0.0;
37
38void Node::node(int i, double c , double l, double u, int cod, int s){
39  index = i;
40  coeff = c;
41  lb = l;
42  ub = u;
43  color = -1;
44  code = cod;
45  sign = s;
46}
47
48inline bool CbcSymmetry::compare (register Node &a, register Node &b) const {
49  if(a.get_code() == b.get_code() )
50    if(a.get_coeff() == b.get_coeff() )
51      if(a.get_sign() == b.get_sign() )
52        if( fabs ( a.get_lb() - b.get_lb() ) <= COUENNE_HACKED_EPS )
53          if( fabs ( a.get_ub() - b.get_ub() ) <= COUENNE_HACKED_EPS )
54            return 1; 
55  return 0;
56}
57
58void CbcSymmetry::Compute_Symmetry() const{
59
60
61  std::sort(node_info_. begin (), node_info_. end (), node_sort);
62
63  for (std::vector <Node>:: iterator i = node_info_. begin (); i != node_info_. end (); ++i) 
64    (*i).color_vertex(-1);
65 
66  int color = 1;
67  for (std::vector <Node>:: iterator i = node_info_. begin (); i != node_info_. end (); ++i) {
68    if( (*i).get_color() == -1){
69      (*i).color_vertex(color);
70#ifdef PRINT_MORE
71      printf ("Graph vertex %d is given color %d\n", (*i).get_index(), color);
72#endif
73      nauty_info_ -> color_node((*i).get_index(), color);
74      for (std::vector <Node>:: iterator j = i+1; j != node_info_. end (); ++j)
75        if( compare( (*i) , (*j) ) ==1){
76          (*j).color_vertex(color);
77          nauty_info_ -> color_node((*j).get_index(),color);
78#ifdef PRINT_MORE
79          printf ("Graph vertex %d is given color %d, the same as vertex %d\n", (*j).get_index(), color, (*i).get_index());
80#endif
81        }
82      //       else
83      // j = node_info_. end();
84      color++;
85    }
86  }
87
88  //Print_Orbits ();
89  nauty_info_ -> computeAuto();
90  //Print_Orbits ();
91}
92
93int
94CbcSymmetry::statsOrbits(CbcModel * model, int type) const
95{
96  char general[200];
97  int returnCode=0;
98  if (type) {
99    double branchSuccess=0.0;
100    if (nautyBranchSucceeded_) 
101      branchSuccess = nautyOtherBranches_/nautyBranchSucceeded_;
102    double fixSuccess=0.0;
103    if (nautyFixSucceeded_) 
104      fixSuccess = nautyFixes_/nautyFixSucceeded_;
105    sprintf(general,"Orbital branching tried %d times, succeeded %d times - average extra %7.3f, fixing %d times (%d, %7.3f) - %.2f seconds",
106            nautyBranchCalls_,nautyBranchSucceeded_,branchSuccess,
107            nautyFixCalls_,nautyFixSucceeded_,fixSuccess,nautyTime_);
108  } else {
109    returnCode = nauty_info_->getNumGenerators();
110    if (returnCode) {
111      sprintf (general,"Nauty: %d orbits, %d generators, group size: %g - dense size %d, sparse %d - going %s",
112               nauty_info_->getNumOrbits(),
113               nauty_info_ -> getNumGenerators () ,
114               nauty_info_ -> getGroupSize (),
115               whichOrbit_[0],whichOrbit_[1],nauty_info_->isSparse() ? "sparse" : "dense");
116    } else {
117      if ((model->moreSpecialOptions2()&(128|256))!=(128|256)) 
118        sprintf(general,"Nauty did not find any useful orbits");
119      else
120        sprintf(general,"Nauty did not find any useful orbits - but keeping Nauty on");
121    }
122  }
123  model->messageHandler()->message(CBC_GENERAL,
124                                   model->messages())
125    << general << CoinMessageEol ;
126  return returnCode;
127}
128 
129void CbcSymmetry::Print_Orbits () const {
130
131  //printf ("num gens = %d, num orbits = %d \n", nauty_info_ -> getNumGenerators(), nauty_info_ -> getNumOrbits() );
132
133  std::vector<std::vector<int> > *new_orbits = nauty_info_->getOrbits();
134
135  printf ("Nauty: %d generators, group size: %.0g",
136          //  nauty_info_->getNumOrbits(),
137          nauty_info_ -> getNumGenerators () ,
138          nauty_info_ -> getGroupSize ());
139
140  int nNonTrivialOrbits = 0;
141
142  for (unsigned int i = 0; i < new_orbits -> size(); i++) {
143    if ((*new_orbits)[i].size() > 1) 
144      nNonTrivialOrbits++;
145    else continue;
146
147    // int orbsize = (*new_orbits)[i].size();
148    // printf( "Orbit %d [size: %d] [", i, orbsize);
149    // copy ((*new_orbits)[i].begin(), (*new_orbits)[i].end(),
150    //    std::ostream_iterator<int>(std::cout, " "));
151    // printf("] \n");
152  }
153
154  printf (" (%d non-trivial orbits).\n", nNonTrivialOrbits);
155
156#if 1
157  if (nNonTrivialOrbits) {
158
159    int orbCnt = 0;
160
161    std::vector<std::vector<int> > *orbits = nauty_info_ -> getOrbits ();
162
163    for   (std::vector <std::vector<int> >::iterator i = orbits -> begin ();  i != orbits -> end (); ++i) {
164
165      printf ("Orbit %d: ", orbCnt++);
166
167      for (std::vector<int>::iterator j = i -> begin (); j != i -> end (); ++j)
168        printf (" %d", *j);
169
170      printf ("\n");
171    }
172  }
173#endif
174
175
176#if 0
177  if (nNonTrivialOrbits)
178    for (int i=0; i< nVars (); i++) {
179
180      std::vector< int > *branch_orbit = Find_Orbit (i);
181
182      if (branch_orbit -> size () > 1) {
183        printf ("x%04d: ", i);
184
185        for (std::vector<int>::iterator it = branch_orbit -> begin (); it != branch_orbit -> end (); ++it)
186          printf ("%d ", *it);
187        printf ("\n");
188      }
189    }
190#endif
191  delete new_orbits;
192}
193void 
194CbcSymmetry::fillOrbits()
195{
196  for (int i=0;i<numberColumns_;i++)
197    whichOrbit_[i]=-1;
198  numberUsefulOrbits_=0;
199
200  std::vector<std::vector<int> > *orbits = nauty_info_ -> getOrbits ();
201
202  for   (std::vector <std::vector<int> >::iterator i = orbits -> begin ();  i != orbits -> end (); ++i) {
203    int nUseful=0;
204    int jColumn=-2;
205    for (std::vector<int>::iterator j = i -> begin (); j != i -> end (); ++j) {
206      int iColumn=*j;
207      if (iColumn<numberColumns_) {
208        whichOrbit_[iColumn]=numberUsefulOrbits_;
209        nUseful++;
210        jColumn=iColumn;
211      }
212    }
213    if (nUseful>1) {
214      numberUsefulOrbits_++;
215    } else if (jColumn>=0) {
216      assert (nUseful);
217      whichOrbit_[jColumn]=-2;
218    }
219  }
220  delete orbits;
221}
222int 
223CbcSymmetry::largestOrbit(const double * lower, const double * upper) const
224{
225  int * counts = new int[numberUsefulOrbits_];
226  memset(counts,0,numberUsefulOrbits_*sizeof(int));
227  for (int i=0;i<numberColumns_;i++) {
228    int iOrbit=whichOrbit_[i];
229    if (iOrbit>=0) {
230      if (lower[i]==0.0&&upper[i]==1.0)
231        counts[iOrbit]++;
232    }
233  }
234  int iOrbit=-1;
235  int maxOrbit=0;
236  for (int i=0;i<numberUsefulOrbits_;i++) {
237    if (counts[i]>maxOrbit) {
238      maxOrbit=counts[i];
239      iOrbit=i;
240    }
241  }
242  delete [] counts;
243  return iOrbit;
244}
245
246std::vector<int> *CbcSymmetry::Find_Orbit(int index) const{
247
248  std::vector<int> *orbit = new std::vector <int>;
249  int which_orbit = -1;
250  std::vector<std::vector<int> > *new_orbits = nauty_info_->getOrbits();
251
252  for (unsigned int i = 0; i < new_orbits -> size(); i++) {
253    for (unsigned int j = 0; j < (*new_orbits)[i].size(); j++) {
254      //   for (std::vector <int>:: iterator j = new_orbits[i].begin(); new_orbits[i].end(); ++j){
255      if( (*new_orbits)[i][j] ==  index)
256        which_orbit = i;
257    }
258  }
259 
260  //  for (std::vector <int>:: iterator j = new_orbits[which_orbit].begin(); new_orbits[which_orbit].end(), ++j)
261  for (unsigned int j = 0; j < (*new_orbits)[which_orbit].size(); j++) 
262    orbit -> push_back ((*new_orbits)[which_orbit][j]);
263
264  delete new_orbits;
265
266  return orbit;
267}
268
269
270void CbcSymmetry::ChangeBounds (const double * new_lb, const double * new_ub, 
271                                int num_cols, bool justFixedAtOne) const {
272  if (justFixedAtOne)
273    nautyFixCalls_++;
274  else
275    nautyBranchCalls_++;
276  std::sort(node_info_. begin (), node_info_. end (), index_sort);
277
278  for (int  i = 0; i < num_cols; i++) {
279    //   printf("Var %d  lower bound: %f   upper bound %f \n", i, new_lb[i], new_ub[i]);
280   
281    assert (node_info_[i].get_index () == i);
282    double newLower = new_lb[i];
283    double newUpper = new_ub[i];
284    if (justFixedAtOne) {
285      // free up all fixed at zero
286      if (!newLower) 
287        newUpper = 1.0;
288    }
289    node_info_[i ].bounds ( newLower , newUpper );
290    //printf("Var %d  INPUT lower bound: %f   upper bound %f \n", i, node_info_[i].get_lb(), node_info_[i].get_ub());
291  }
292}
293void CbcSymmetry::setupSymmetry (const OsiSolverInterface & solver) {
294  const double *objective = solver.getObjCoefficients() ;
295  const double *columnLower = solver.getColLower() ;
296  const double *columnUpper = solver.getColUpper() ;
297  int numberColumns = solver.getNumCols() ;
298  int numberRows = solver.getNumRows();
299  int iRow, iColumn;
300 
301  // Row copy
302  CoinPackedMatrix matrixByRow(*solver.getMatrixByRow());
303  const double * elementByRow = matrixByRow.getElements();
304  const int * column = matrixByRow.getIndices();
305  const CoinBigIndex * rowStart = matrixByRow.getVectorStarts();
306  const int * rowLength = matrixByRow.getVectorLengths();
307
308  const double * rowLower = solver.getRowLower();
309  const double * rowUpper = solver.getRowUpper();
310  //  // Find Coefficients
311
312  /// initialize nauty
313
314  int num_affine = 0;
315
316  for (iColumn = 0; iColumn < numberColumns; iColumn++) {
317    if (objective[iColumn] && objective[iColumn]!=1.0)
318      num_affine++;
319  }
320  for (iRow = 0; iRow < numberRows; iRow++) {
321    for (CoinBigIndex j = rowStart[iRow]; 
322         j < rowStart[iRow] + rowLength[iRow]; j++) {
323      int jColumn = column[j];
324      double value = elementByRow[j];
325      if (value!=1.0)
326        num_affine++;
327    }
328  }
329
330  // Create Nauty object
331
332  int coef_count= numberRows + numberColumns +1;
333  int nc = num_affine + coef_count;
334  // create graph (part 1)
335
336  for (iColumn = 0; iColumn < numberColumns; iColumn++) {
337    Node var_vertex;
338    var_vertex.node(iColumn,0.0,columnLower[iColumn],columnUpper[iColumn],-1,-1 );
339    node_info_.push_back(var_vertex);
340  }
341  // objective
342  int index = numberColumns;
343  {
344    Node vertex;
345    vertex.node( index , 0.0 , -COIN_DBL_MAX, COIN_DBL_MAX, 
346                 COUENNE_HACKED_EXPRGROUP, 0);
347    node_info_.push_back( vertex);
348  }
349
350  // compute space for sparse
351  size_t * v = NULL;
352  int * d = NULL;
353  int * e = NULL;
354  bool sparse=false;
355  double spaceDense = ((nc+WORDSIZE-1)*(nc+WORDSIZE-1))/WORDSIZE;
356  int spaceSparse = 0;
357  {
358    size_t numberElements = 0;
359    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
360      double value = objective[iColumn];
361      if (value) {
362        if (value==1.0) {
363          numberElements+=2;
364        } else {
365          numberElements+=4;
366          coef_count ++;
367        }
368      }
369    }
370    for (iRow = 0; iRow < numberRows; iRow++) {
371      for (CoinBigIndex j = rowStart[iRow]; 
372           j < rowStart[iRow] + rowLength[iRow]; j++) {
373        int jColumn = column[j];
374        double value = elementByRow[j];
375        if (value==1.0) {
376          numberElements+=2;
377        } else {
378          numberElements+=4;
379          coef_count ++;
380        }
381      }
382    }
383    spaceSparse = 2*nc+numberElements;
384    //printf("Space for sparse is %d for dense %g\n",
385    //     spaceSparse,spaceDense);
386#ifdef NTY_TRACES
387    bool goSparse = true;
388#else
389    bool goSparse = (spaceSparse<spaceDense);
390#endif
391    // for now always sparse
392    goSparse=true;
393    if (goSparse) {
394      sparse=true;
395      v = new size_t [nc+1];
396      d = new int [nc];
397      e = new int [numberElements];
398      size_t * counts = new size_t [coef_count+1];
399      memset(counts,0,coef_count*sizeof(size_t));
400      coef_count= numberRows + numberColumns +1;
401      // create graph (part 2)
402      for (iColumn = 0; iColumn < numberColumns; iColumn++) {
403        double value = objective[iColumn];
404        if (value) {
405          if (value==1.0) {
406            counts[index]++;
407            counts[iColumn]++;
408          } else {
409            counts[index]++;
410            counts[coef_count]+=2;
411            counts[iColumn]++;
412            coef_count ++;
413          }
414        }
415      }
416      index++;
417      for (iRow = 0; iRow < numberRows; iRow++) {
418        for (CoinBigIndex j = rowStart[iRow]; 
419             j < rowStart[iRow] + rowLength[iRow]; j++) {
420          int jColumn = column[j];
421          double value = elementByRow[j];
422          if (value==1.0) {
423            counts[index]++;
424            counts[jColumn]++;
425          } else {
426            counts[index]++;
427            counts[coef_count]+=2;
428            counts[jColumn]++;
429            coef_count ++;
430          }
431        }
432        index++;
433      }
434      // create graph (part 3)
435      assert (nc==coef_count);
436      numberElements=0;
437      v[0]=0;
438      for (int i=0;i<nc;i++) {
439        int count=counts[i];
440        d[i]=count;
441        counts[i]=v[i];
442        numberElements+=count;;
443        v[i+1]=numberElements;
444      }
445      index = numberColumns;
446      coef_count= numberRows + numberColumns +1;
447      for (iColumn = 0; iColumn < numberColumns; iColumn++) {
448        double value = objective[iColumn];
449        if (value) {
450          int where;
451          if (value==1.0) {
452            where = counts[index];
453            counts[index]++;
454            e[where]=iColumn;
455            where = counts[iColumn];
456            counts[iColumn]++;
457            e[where]=index;
458          } else {
459            Node coef_vertex;
460            coef_vertex.node( coef_count, value, value, value, -2, 0 );
461            node_info_.push_back(coef_vertex);
462            where = counts[index];
463            counts[index]++;
464            e[where]=coef_count;
465            where = counts[coef_count];
466            counts[coef_count]++;
467            e[where]=index;
468            where = counts[coef_count];
469            counts[coef_count]++;
470            e[where]=iColumn;
471            where = counts[iColumn];
472            counts[iColumn]++;
473            e[where]=coef_count;
474            coef_count ++;
475          }
476      }
477      }
478      index++;
479      for (iRow = 0; iRow < numberRows; iRow++) {
480        Node vertex;
481        vertex.node( index , 0.0 , rowLower[iRow], rowUpper[iRow],
482                     COUENNE_HACKED_EXPRGROUP, 0);
483        node_info_.push_back( vertex);
484        for (CoinBigIndex j = rowStart[iRow]; 
485             j < rowStart[iRow] + rowLength[iRow]; j++) {
486          int jColumn = column[j];
487          double value = elementByRow[j];
488          int where;
489          if (value==1.0) {
490            where = counts[index];
491            counts[index]++;
492            e[where]=jColumn;
493            where = counts[jColumn];
494            counts[jColumn]++;
495            e[where]=index;
496          } else {
497            Node coef_vertex;
498            coef_vertex.node( coef_count, value, value, value, -2, 0 );
499            node_info_.push_back(coef_vertex);
500            where = counts[index];
501            counts[index]++;
502            e[where]=coef_count;
503            where = counts[coef_count];
504            counts[coef_count]++;
505            e[where]=index;
506            where = counts[coef_count];
507            counts[coef_count]++;
508            e[where]=jColumn;
509            where = counts[jColumn];
510            counts[jColumn]++;
511            e[where]=coef_count;
512            coef_count ++;
513          }
514        }
515        index++;
516      }
517      delete [] counts;
518    }
519  }
520
521  nauty_info_ = new CbcNauty(nc,v,d,e); 
522  delete [] v;
523  delete [] d;
524  delete [] e;
525  if (!sparse) {
526    // create graph (part 2)
527    coef_count= numberRows + numberColumns +1;
528    for (iColumn = 0; iColumn < numberColumns; iColumn++) {
529      double value = objective[iColumn];
530      if (value) {
531        if (value==1.0) {
532          nauty_info_->addElement(index,iColumn);
533          nauty_info_->addElement(iColumn,index);
534        } else {
535          Node coef_vertex;
536          coef_vertex.node( coef_count, value, value, value, -2, 0 );
537          node_info_.push_back(coef_vertex);
538          nauty_info_->addElement(index,  coef_count);
539          nauty_info_->addElement(coef_count, index);
540          nauty_info_->addElement(coef_count,  iColumn);
541          nauty_info_->addElement(iColumn, coef_count);
542          coef_count ++;
543        }
544      }
545    }
546    index++;
547    for (iRow = 0; iRow < numberRows; iRow++) {
548      Node vertex;
549      vertex.node( index , 0.0 , rowLower[iRow], rowUpper[iRow],
550                   COUENNE_HACKED_EXPRGROUP, 0);
551      node_info_.push_back( vertex);
552      for (CoinBigIndex j = rowStart[iRow]; 
553           j < rowStart[iRow] + rowLength[iRow]; j++) {
554        int jColumn = column[j];
555        double value = elementByRow[j];
556        if (value==1.0) {
557          nauty_info_->addElement(index,jColumn);
558          nauty_info_->addElement(jColumn,index);
559        } else {
560          Node coef_vertex;
561          coef_vertex.node( coef_count, value, value, value, -2, 0 );
562          node_info_.push_back(coef_vertex);
563          nauty_info_->addElement(index,  coef_count);
564          nauty_info_->addElement(coef_count, index);
565          nauty_info_->addElement(coef_count,  jColumn);
566          nauty_info_->addElement(jColumn, coef_count);
567          coef_count ++;
568        }
569      }
570      index++;
571    }
572  }
573  numberColumns_ = numberColumns;
574  whichOrbit_ = new int [2*numberColumns_];
575  nautyBranchCalls_ = 0;
576  nautyBranchSucceeded_ = 0;
577  nautyFixCalls_ = 0;
578  nautyFixSucceeded_ = 0;
579  nautyTime_ = 0.0;
580  nautyFixes_= 0.0; 
581  nautyOtherBranches_ = 0.0;
582  Compute_Symmetry ();
583  //Print_Orbits ();
584  // stats in array
585  whichOrbit_[0]=spaceDense;
586  whichOrbit_[1]=spaceSparse;
587}
588// Fixes variables using orbits (returns number fixed)
589int 
590CbcSymmetry::orbitalFixing(OsiSolverInterface * solver)
591{
592  int numberColumns = solver->getNumCols();
593  char * status = new char [numberColumns];
594  ChangeBounds(solver->getColLower(),
595               solver->getColUpper(),
596               solver->getNumCols(),true);
597  Compute_Symmetry();
598  fillOrbits();
599  int n=0;
600  //#define PRINT_MORE 1
601  const int * alternativeOrbits = whichOrbit();
602  if (alternativeOrbits) {
603    for (int i=0;i<numberColumns;i++) {
604      char type='0';
605      if (solver->getColUpper()[i]) {
606        if (solver->getColLower()[i]) {
607          type='1';
608        } else {
609          double value=solver->getColSolution()[i];
610          if (value<0.0001)
611            type='L';
612          else if (value>0.9999)
613            type='U';
614          else
615            type='X';
616        }
617      }
618      status[i]=type;
619    }
620    n=0;
621    for (int i=0;i<numberColumns;i++) {
622      if (status[i]!='0'&&status[i]!='1') {
623        int iOrbit=alternativeOrbits[i];
624        if (iOrbit<0)
625          continue;
626        for (int j=i+1;j<numberColumns;j++) {
627          if (status[j]=='0'&&alternativeOrbits[j]==iOrbit) {
628#if PRINT_MORE>1
629            printf("In alternative orbit %d - %d free (%c), %d fixed to 0\n",
630                   iOrbit,i,status[i],j);
631#endif
632            status[i]='0'; // can fix on both branches
633            solver->setColUpper(i,0.0);
634            n++;
635            break;
636          }
637        }
638      }
639    }
640  }
641  delete [] status;
642  if (n) {
643    nautyFixSucceeded_++;
644    nautyFixes_ += n;
645#if PRINT_MORE
646    printf("%d orbital fixes\n",n);
647#endif
648  }
649  return n;
650}
651// Default Constructor
652CbcSymmetry::CbcSymmetry ()
653  : nauty_info_(NULL),
654    numberColumns_(0),
655    numberUsefulOrbits_(0),
656    whichOrbit_(NULL)
657{
658}
659// Copy constructor
660CbcSymmetry::CbcSymmetry ( const CbcSymmetry & rhs)
661{
662  node_info_ = rhs.node_info_;
663  nauty_info_ = new CbcNauty(*rhs.nauty_info_);
664  numberUsefulOrbits_ = rhs.numberUsefulOrbits_;
665  numberColumns_ = rhs.numberColumns_;
666  if (rhs.whichOrbit_) 
667    whichOrbit_=CoinCopyOfArray(rhs.whichOrbit_,numberColumns_);
668  else
669    whichOrbit_ = NULL;
670}
671
672// Assignment operator
673CbcSymmetry &
674CbcSymmetry::operator=( const CbcSymmetry & rhs)
675{
676  if (this != &rhs) {
677    delete nauty_info_;
678    node_info_ = rhs.node_info_;
679    nauty_info_ = new CbcNauty(*rhs.nauty_info_);
680    delete [] whichOrbit_;
681    numberColumns_ = rhs.numberColumns_;
682    numberUsefulOrbits_ = rhs.numberUsefulOrbits_;
683    if (rhs.whichOrbit_) 
684      whichOrbit_=CoinCopyOfArray(rhs.whichOrbit_,numberColumns_);
685    else
686      whichOrbit_ = NULL;
687  }
688  return *this;
689}
690
691// Destructor
692CbcSymmetry::~CbcSymmetry ()
693{
694  delete nauty_info_; 
695  delete [] whichOrbit_;
696}
697
698CbcNauty::CbcNauty(int vertices, const size_t * v, const int * d, const int * e)
699{
700  //printf("Need sparse nauty - wordsize %d\n",WORDSIZE);
701  n_ = vertices;
702  m_ = (n_ + WORDSIZE - 1)/WORDSIZE;
703  if (v)
704    nel_ = v[n_];
705  else
706    nel_ = 0;
707
708  //printf ("size of long = %d (%d)\nwordsize = %d\nn,m = %d,%d\n",
709  //          SIZEOF_LONG, sizeof (long), WORDSIZE, n_, m_);
710
711  nauty_check (WORDSIZE, m_, n_, NAUTYVERSIONID);
712
713  /// Apparently sizes are skewed on 64bit machines
714
715#define MULTIPLIER 1
716
717  if (!nel_) {
718    G_ = (graph *) malloc(MULTIPLIER * m_ * n_ * sizeof(int));
719    GSparse_ = NULL;
720  } else {
721    G_ = NULL;
722    GSparse_ = (sparsegraph *) malloc(sizeof(sparsegraph));
723    SG_INIT(*GSparse_);
724    SG_ALLOC(*GSparse_,n_,nel_,"malloc");
725    GSparse_->nv = n_; /* Number of vertices */
726    GSparse_->nde = nel_;
727  }
728  lab_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int)); 
729  ptn_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int));
730  active_ = NULL;
731  orbits_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int));
732#ifndef NTY_TRACES
733  options_ = (optionblk *) malloc(MULTIPLIER * sizeof(optionblk));
734  stats_ = (statsblk *) malloc(MULTIPLIER * sizeof(statsblk));
735#else
736  options_ = (TracesOptions *) malloc(MULTIPLIER * sizeof(TracesOptions));
737  stats_ = (TracesStats *) malloc(MULTIPLIER * sizeof(TracesStats));
738#endif
739  worksize_ = 100*m_;
740  workspace_ = (setword *) malloc(MULTIPLIER * worksize_*sizeof(setword));
741  canonG_ = NULL;
742  if ((G_ == 0&&GSparse_ == 0) || lab_ == 0 || ptn_ == 0 || 
743      orbits_ == 0 || options_ == 0 || stats_ == 0 ||
744      workspace_ == 0) assert(0);
745
746  // Zero allocated memory
747  if (G_) {
748    memset(G_, 0, m_*n_*sizeof(int));
749  } else {
750    //for (int i=0;i<n_;i++) {
751    //GSparse_->v[i]=v[i];
752    //}
753    memcpy(GSparse_->v,v,n_*sizeof(size_t));
754    memcpy(GSparse_->d,d,n_*sizeof(int));
755    memcpy(GSparse_->e,e,nel_*sizeof(int));
756  }
757  memset(lab_, 0, n_*sizeof(int));
758  memset(ptn_, 0, n_*sizeof(int));
759  memset(orbits_, 0, n_*sizeof(int));
760  memset(workspace_, 0, worksize_*sizeof(setword));
761#ifndef NTY_TRACES
762  memset(options_, 0,MULTIPLIER * sizeof(optionblk));
763#else
764  memset(options_, 0,MULTIPLIER * sizeof(TracesOptions));
765#endif
766
767  // Set the options you want
768#ifndef NTY_TRACES
769  options_->getcanon = FALSE;
770  options_->digraph = FALSE;
771  options_->writeautoms = FALSE;
772  options_->writemarkers = FALSE;
773  options_->defaultptn = TRUE;
774  options_->cartesian = FALSE;
775  options_->linelength = 78;
776  options_->outfile = NULL;
777  options_->userrefproc = NULL;
778  options_->userautomproc = NULL;
779  options_->userlevelproc = NULL;
780  options_->usernodeproc = NULL;
781  //  options_->usertcellproc = NULL;
782  options_->invarproc = NULL;
783  options_->tc_level = 100;
784  options_->mininvarlevel = 0;
785  options_->maxinvarlevel = 1;
786  options_->invararg = 0;
787  options_->dispatch = &dispatch_graph;
788#else
789  options_->getcanon = FALSE;
790  options_->writeautoms = FALSE;
791  options_->cartesian = FALSE;
792  options_->digraph = FALSE;
793  options_->defaultptn = TRUE;
794  options_->linelength = 78;
795#endif
796  if (G_) {
797    // Make an empty graph
798    for (int j = 0; j < n_; j++) {
799      set *gv = GRAPHROW(G_, j, m_);
800      EMPTYSET(gv, m_);
801    }
802  }
803
804  vstat_ = new int[n_];
805  clearPartitions();
806  afp_ = NULL;
807}
808
809CbcNauty::~CbcNauty()
810{
811  if (G_) free(G_);
812  if (GSparse_) {
813    SG_FREE(*GSparse_);
814    free(GSparse_);
815  }
816  if (lab_) free(lab_);
817  if (ptn_) free(ptn_);
818  if (active_) free(active_);
819  if (orbits_) free(orbits_);
820  if (options_) free(options_);
821  if (stats_) free(stats_);
822  if (workspace_) free(workspace_);
823  if (canonG_) free(canonG_);
824  if (vstat_) delete [] vstat_;
825}
826// Copy constructor
827CbcNauty::CbcNauty ( const CbcNauty & rhs)
828{
829  n_ = rhs.n_;
830  m_ = rhs.m_;
831  nel_ = rhs.nel_;
832  G_ = NULL;
833  GSparse_ = NULL;
834  if (!nel_) {
835    G_ = (graph *) malloc(MULTIPLIER * m_ * n_ * sizeof(int));
836  } else {
837    GSparse_ = (sparsegraph *) malloc(sizeof(sparsegraph));
838    SG_INIT(*GSparse_);
839    SG_ALLOC(*GSparse_,n_,nel_,"malloc");
840    GSparse_->nv = n_; /* Number of vertices */
841    GSparse_->nde = nel_;
842  }
843  lab_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int)); 
844  ptn_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int));
845  orbits_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int));
846#ifndef NTY_TRACES
847  options_ = (optionblk *) malloc(MULTIPLIER * sizeof(optionblk));
848  stats_ = (statsblk *) malloc(MULTIPLIER * sizeof(statsblk));
849#else
850  options_ = (TracesOptions *) malloc(MULTIPLIER * sizeof(TracesOptions));
851  stats_ = (TracesStats *) malloc(MULTIPLIER * sizeof(TracesStats));
852#endif
853  worksize_ = 100*m_;
854  workspace_ = (setword *) malloc(MULTIPLIER * worksize_*sizeof(setword));
855  vstat_ = new int[n_];
856  canonG_ = NULL;
857  if ((G_ == 0 && GSparse_ == 0) || lab_ == 0 || ptn_ == 0 || 
858      orbits_ == 0 || options_ == 0 || stats_ == 0 ||
859      workspace_ == 0) assert(0);
860
861  // Copy allocated memory
862  if (G_) {
863    memcpy(G_, rhs.G_, m_*n_*sizeof(int));
864  } else {
865    memcpy(GSparse_->v,rhs.GSparse_->v,n_*sizeof(size_t));
866    memcpy(GSparse_->d,rhs.GSparse_->d,n_*sizeof(int));
867    memcpy(GSparse_->e,rhs.GSparse_->e,nel_*sizeof(int));
868  }
869  memcpy(lab_, rhs.lab_, n_*sizeof(int));
870  memcpy(ptn_, rhs.ptn_, n_*sizeof(int));
871  memcpy(orbits_, rhs.orbits_, n_*sizeof(int));
872  memcpy(workspace_, rhs.workspace_, worksize_*sizeof(setword));
873#ifndef NTY_TRACES
874  memcpy(options_, rhs.options_, MULTIPLIER * sizeof(optionblk));
875  memcpy(stats_, rhs.stats_, MULTIPLIER * sizeof(statsblk));
876#else
877  memcpy(options_, rhs.options_,MULTIPLIER * sizeof(TracesOptions));
878  memcpy(stats_, rhs.stats_, MULTIPLIER * sizeof(TracesStats));
879#endif
880  memcpy(vstat_,rhs.vstat_,n_*sizeof(int));
881
882  // ? clearPartitions();
883  active_ = NULL;
884  afp_ = rhs.afp_; // ? no copy ?
885}
886
887// Assignment operator
888CbcNauty &
889CbcNauty::operator=( const CbcNauty & rhs)
890{
891  if (this != &rhs) {
892    if (G_) free(G_);
893    if (GSparse_) {
894      SG_FREE(*GSparse_);
895      free(GSparse_);
896    }
897    if (lab_) free(lab_);
898    if (ptn_) free(ptn_);
899    if (active_) free(active_);
900    if (orbits_) free(orbits_);
901    if (options_) free(options_);
902    if (stats_) free(stats_);
903    if (workspace_) free(workspace_);
904    if (canonG_) free(canonG_);
905    if (vstat_) delete [] vstat_;
906    {
907      n_ = rhs.n_;
908      m_ = rhs.m_;
909      nel_ = rhs.nel_;
910      G_ = NULL;
911      GSparse_ = NULL;
912      if (!nel_) {
913        G_ = (graph *) malloc(MULTIPLIER * m_ * n_ * sizeof(int));
914      } else {
915        GSparse_ = (sparsegraph *) malloc(sizeof(sparsegraph));
916        SG_INIT(*GSparse_);
917        SG_ALLOC(*GSparse_,n_,nel_,"malloc");
918        GSparse_->nv = n_; /* Number of vertices */
919        GSparse_->nde = nel_;
920      }
921      lab_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int)); 
922      ptn_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int));
923      orbits_ = (int *) malloc(MULTIPLIER * n_ * sizeof(int));
924#ifndef NTY_TRACES
925      options_ = (optionblk *) malloc(MULTIPLIER * sizeof(optionblk));
926      stats_ = (statsblk *) malloc(MULTIPLIER * sizeof(statsblk));
927#else
928      options_ = (TracesOptions *) malloc(MULTIPLIER * sizeof(TracesOptions));
929      stats_ = (TracesStats *) malloc(MULTIPLIER * sizeof(TracesStats));
930#endif
931      worksize_ = 100*m_;
932      workspace_ = (setword *) malloc(MULTIPLIER * worksize_*sizeof(setword));
933      vstat_ = new int[n_];
934      canonG_ = NULL;
935      if ((G_ == 0 && GSparse_ == 0) || lab_ == 0 || ptn_ == 0 || 
936          orbits_ == 0 || options_ == 0 || stats_ == 0 ||
937          workspace_ == 0) assert(0);
938     
939      // Copy allocated memory
940      if (!nel_) {
941        memcpy(G_, rhs.G_, m_*n_*sizeof(int));
942      } else {
943        memcpy(GSparse_->v,rhs.GSparse_->v,n_*sizeof(size_t));
944        memcpy(GSparse_->d,rhs.GSparse_->d,n_*sizeof(int));
945        memcpy(GSparse_->e,rhs.GSparse_->e,nel_*sizeof(int));
946      }
947      memcpy(lab_, rhs.lab_, n_*sizeof(int));
948      memcpy(ptn_, rhs.ptn_, n_*sizeof(int));
949      memcpy(orbits_, rhs.orbits_, n_*sizeof(int));
950      memcpy(workspace_, rhs.workspace_, worksize_*sizeof(setword));
951#ifndef NTY_TRACES
952      memcpy(options_, rhs.options_, MULTIPLIER * sizeof(optionblk));
953      memcpy(stats_, rhs.stats_, MULTIPLIER * sizeof(statsblk));
954#else
955      memcpy(options_, rhs.options_,MULTIPLIER * sizeof(TracesOptions));
956      memcpy(stats_, rhs.stats_, MULTIPLIER * sizeof(TracesStats));
957#endif
958      memcpy(vstat_,rhs.vstat_,n_*sizeof(int));
959     
960      // ? clearPartitions();
961      active_ = NULL;
962      afp_ = rhs.afp_; // ? no copy ?
963    }
964  }
965  return *this;
966}
967
968void
969CbcNauty::addElement(int ix, int jx)
970{
971  // Right now die if bad index.  Can throw exception later
972  //printf("addelement %d %d \n", ix, jx);
973  assert(ix < n_ && jx < n_);
974  if(ix != jx){  //No Loops
975    set *gv = GRAPHROW(G_, ix, m_);
976    ADDELEMENT(gv, jx);
977    set *gv2 = GRAPHROW(G_, jx, m_);
978    ADDELEMENT(gv2, ix);
979    autoComputed_ = false;
980  }
981}
982
983void 
984CbcNauty::clearPartitions()
985{
986  for (int j = 0; j < n_; j++) {
987    vstat_[j] = 1;
988    //printf("vstat %d = %d", j, vstat_[j]);
989  }
990
991  autoComputed_ = false;
992}
993
994void
995CbcNauty::computeAuto()
996{
997
998  //  if (autoComputed_) return;
999
1000  double startCPU = CoinCpuTime ();
1001
1002  options_->defaultptn = FALSE;
1003
1004  // Here we only implement the partitions
1005  // [ fix1 | fix0 (union) free | constraints ]
1006  int ix = 0;
1007 
1008  for( int color = 1; color <= n_; color++){
1009    for (int j = 0; j < n_; j++) {
1010      if (vstat_[j] == color) {
1011        lab_[ix] = j;
1012        ptn_[ix] = color;
1013        ix++;
1014      }
1015    }
1016     if (ix > 0) ptn_[ix-1] = 0;
1017  }
1018 
1019  /*
1020  for (int j = 0; j < n_; j++)
1021    printf("ptn %d = %d      lab = %d \n", j, ptn_[j], lab_[j]);
1022  */
1023
1024  // Should be number of columns
1025  assert(ix == n_);
1026  // Now the constraints if needed
1027 
1028  // Compute Partition
1029   
1030  if (G_) {
1031#ifndef NTY_TRACES
1032    nauty(G_, lab_, ptn_, active_, orbits_, options_, 
1033          stats_, workspace_, worksize_, m_, n_, canonG_);
1034#else
1035    abort();
1036#endif
1037  } else {
1038#ifndef NTY_TRACES
1039    options_->dispatch = &dispatch_sparse;
1040    sparsenauty(GSparse_, lab_, ptn_, orbits_, options_, 
1041          stats_, NULL);
1042#else
1043    //options_->dispatch = &dispatch_sparse;
1044    Traces(GSparse_, lab_, ptn_, orbits_, options_, 
1045          stats_, NULL);
1046#endif
1047  }
1048  autoComputed_ = true;
1049
1050  double endCPU = CoinCpuTime ();
1051
1052  nautyTime_ += endCPU - startCPU;
1053  // Need to make sure all generators are written
1054  if (afp_) fflush(afp_);   
1055}
1056
1057void
1058CbcNauty::deleteElement(int ix, int jx)
1059{
1060  // Right now die if bad index.  Can throw exception later
1061  assert(ix < n_ && jx < n_);
1062  set *gv = GRAPHROW(G_, ix, m_);
1063  if (ISELEMENT(gv, jx)) {
1064    DELELEMENT(gv, jx);
1065  } 
1066  autoComputed_ = false;
1067}
1068
1069double
1070CbcNauty::getGroupSize() const
1071{
1072  if (!autoComputed_) return -1.0;
1073  return( stats_->grpsize1 * pow(10.0, (double) stats_->grpsize2) );
1074}
1075
1076int
1077CbcNauty::getNumGenerators() const
1078{
1079  if (!autoComputed_) return -1;
1080  return(stats_->numgenerators);
1081}
1082
1083int
1084CbcNauty::getNumOrbits() const
1085{
1086  if (!autoComputed_) return -1;
1087  return(stats_->numorbits);
1088}
1089
1090std::vector<std::vector<int> >
1091*CbcNauty::getOrbits() const
1092{
1093  std::vector<std::vector<int> > *orb = new std::vector<std::vector<int> >;
1094  if (!autoComputed_) return orb;
1095  orb -> resize(getNumOrbits());
1096  std::multimap<int, int> orbmap;
1097  std::set<int> orbkeys;
1098  for (int j = 0; j < n_; j++) {
1099    orbkeys.insert(orbits_[j]);
1100    orbmap.insert(std::make_pair(orbits_[j], j));
1101  }
1102
1103  int orbix = 0;
1104  for (std::set<int>::iterator it = orbkeys.begin();
1105       it != orbkeys.end(); ++it) {
1106    std::multimap<int, int>::iterator pos;
1107    for (pos = orbmap.lower_bound(*it);
1108         pos != orbmap.upper_bound(*it); ++pos) {
1109      (*orb)[orbix].push_back(pos->second);
1110    }
1111    orbix++;
1112  }
1113
1114  assert(orbix == getNumOrbits());
1115  return orb;
1116}
1117
1118void
1119CbcNauty::getVstat(double *v, int nv)
1120{
1121  assert(nv == n_);
1122  memcpy(v, vstat_, nv * sizeof(VarStatus));
1123}
1124
1125/*
1126bool
1127CbcNauty::isAllFixOneOrbit(const std::vector<int> &orbit) const
1128{
1129
1130  for(std::vector<int>::const_iterator it = orbit.begin();
1131      it != orbit.end(); ++it) {
1132    if (*it >= n_) return false;
1133    if (vstat_[*it] != FIX_AT_ONE) return false;
1134  }
1135  return true;
1136}
1137
1138bool
1139CbcNauty::isAllFreeOrbit(const std::vector<int> &orbit) const
1140{
1141  for(std::vector<int>::const_iterator it = orbit.begin();
1142      it != orbit.end(); ++it) {
1143    if (*it >= n_) return false;
1144    if (vstat_[*it] != FREE) return false;
1145  }
1146  return true;
1147}
1148
1149bool
1150CbcNauty::isConstraintOrbit(const std::vector<int> &orbit) const
1151{
1152  for(std::vector<int>::const_iterator it = orbit.begin();
1153      it != orbit.end(); ++it) {
1154    if (*it >= n_) return true;
1155  }
1156  return false;
1157 
1158}
1159
1160bool
1161CbcNauty::isMixedFreeZeroOrbit(const std::vector<int> &orbit) const
1162{
1163  bool containsFree = false;
1164  bool containsZero = false;
1165
1166  for(std::vector<int>::const_iterator it = orbit.begin();
1167      it != orbit.end(); ++it) {
1168    if (*it >= n_) return false;   
1169    if (vstat_[*it] == FREE) containsFree = true;
1170    if (vstat_[*it] == FIX_AT_ZERO) containsZero = true;   
1171    if (containsFree && containsZero) break;
1172  } 
1173  return (containsFree && containsZero);
1174}
1175*/
1176
1177void 
1178CbcNauty::setWriteAutoms(const std::string &fname)
1179{
1180  afp_ = fopen(fname.c_str(), "w");
1181  options_->writeautoms = TRUE;
1182#ifndef NTY_TRACES
1183  options_->writemarkers = FALSE;
1184#endif
1185  options_->outfile = afp_;
1186
1187}
1188
1189void 
1190CbcNauty::unsetWriteAutoms()
1191{
1192  fclose(afp_);
1193  options_->writeautoms = FALSE;
1194}
1195
1196// Default Constructor
1197CbcOrbitalBranchingObject::CbcOrbitalBranchingObject()
1198        : CbcBranchingObject(),
1199          column_(-1),
1200          numberOther_(0),
1201          numberExtra_(0),
1202          fixToZero_(NULL)
1203{
1204}
1205
1206// Useful constructor
1207CbcOrbitalBranchingObject::CbcOrbitalBranchingObject (CbcModel * model, int column,
1208                                                      int way ,
1209                                                      int numberExtra,
1210                                                      const int * extraToZero)
1211  : CbcBranchingObject(model, -1, way, 0.5),
1212    column_(column),
1213    numberOther_(0),
1214    numberExtra_(0),
1215    fixToZero_(NULL)
1216{
1217  CbcSymmetry * symmetryInfo = model->symmetryInfo();
1218  assert (symmetryInfo);
1219  // Filled in (hopefully)
1220  const int * orbit = symmetryInfo->whichOrbit();
1221  int iOrbit=orbit[column];
1222  assert (iOrbit>=0);
1223  int numberColumns = model->getNumCols();
1224  numberOther_=-1;
1225  for (int i=0;i<numberColumns;i++) {
1226    if (orbit[i]==iOrbit)
1227      numberOther_++;
1228  }
1229  assert (numberOther_>0);
1230  nautyBranchSucceeded_++;
1231  nautyOtherBranches_ += numberOther_;
1232  numberExtra_ = numberExtra;
1233  fixToZero_ = new int [numberOther_+numberExtra_];
1234  int n=0;
1235  for (int i=0;i<numberColumns;i++) {
1236    if (orbit[i]==iOrbit && i != column)
1237      fixToZero_[n++]=i;
1238  }
1239  for (int i=0;i<numberExtra;i++) {
1240      fixToZero_[n++]=extraToZero[i];
1241  }
1242}
1243
1244// Copy constructor
1245CbcOrbitalBranchingObject::CbcOrbitalBranchingObject (const CbcOrbitalBranchingObject & rhs)
1246        : CbcBranchingObject(rhs),
1247          column_(rhs.column_),
1248          numberOther_(rhs.numberOther_),
1249          numberExtra_(rhs.numberExtra_)
1250{
1251  fixToZero_ = CoinCopyOfArray(rhs.fixToZero_,numberOther_+numberExtra_);
1252}
1253
1254// Assignment operator
1255CbcOrbitalBranchingObject &
1256CbcOrbitalBranchingObject::operator=( const CbcOrbitalBranchingObject & rhs)
1257{
1258    if (this != &rhs) {
1259        CbcBranchingObject::operator=(rhs);
1260        delete [] fixToZero_;
1261        column_=rhs.column_;
1262        numberOther_=rhs.numberOther_;
1263        numberExtra_=rhs.numberExtra_;
1264        fixToZero_ = CoinCopyOfArray(rhs.fixToZero_,numberOther_+numberExtra_);
1265    }
1266    return *this;
1267}
1268CbcBranchingObject *
1269CbcOrbitalBranchingObject::clone() const
1270{
1271    return (new CbcOrbitalBranchingObject(*this));
1272}
1273
1274
1275// Destructor
1276CbcOrbitalBranchingObject::~CbcOrbitalBranchingObject ()
1277{
1278  delete [] fixToZero_;
1279}
1280
1281double
1282CbcOrbitalBranchingObject::branch()
1283{
1284  decrementNumberBranchesLeft();
1285  if (model_->logLevel()>1)
1286    print();
1287  OsiSolverInterface * solver = model_->solver();
1288  if (way_ < 0) {
1289    solver->setColUpper(column_,0.0);
1290    for ( int i = 0; i < numberOther_+numberExtra_; i++) {
1291      solver->setColUpper(fixToZero_[i],0.0);
1292    }
1293    way_ = 1;     // Swap direction
1294 } else {
1295    solver->setColLower(column_,1.0);
1296    for (int i = numberOther_; i < numberOther_+numberExtra_; i++) {
1297      solver->setColUpper(fixToZero_[i],0.0);
1298    }
1299    way_ = -1;    // Swap direction
1300  }
1301  return 0.0;
1302}
1303/* Update bounds in solver as in 'branch' and update given bounds.
1304   branchState is -1 for 'down' +1 for 'up' */
1305void
1306CbcOrbitalBranchingObject::fix(OsiSolverInterface * solver,
1307                           double * lower, double * upper,
1308                           int branchState) const
1309{
1310  if (branchState < 0) {
1311    upper[column_]=0.0;
1312    for ( int i = 0; i < numberOther_+numberExtra_; i++) {
1313      upper[fixToZero_[i]]=0.0;;
1314    }
1315  } else {
1316    lower[column_]=1.0;
1317    for (int i = numberOther_; i < numberOther_+numberExtra_; i++) {
1318      upper[fixToZero_[i]]=0.0;;
1319    }
1320  }
1321}
1322// Print what would happen
1323void
1324CbcOrbitalBranchingObject::print()
1325{
1326    if (way_ < 0) {
1327      printf("Orbital Down - to zero %d",column_);
1328      for ( int i = 0; i < numberOther_+numberExtra_; i++) {
1329        printf(" %d",fixToZero_[i]);
1330      }
1331    } else {
1332      printf("Orbital Up - to one %d, to zero",column_);
1333      for (int i = numberOther_; i < numberOther_+numberExtra_; i++) {
1334        printf(" %d",fixToZero_[i]);
1335      }
1336    }
1337    printf("\n");
1338}
1339
1340/** Compare the original object of \c this with the original object of \c
1341    brObj. Assumes that there is an ordering of the original objects.
1342    This method should be invoked only if \c this and brObj are of the same
1343    type.
1344    Return negative/0/positive depending on whether \c this is
1345    smaller/same/larger than the argument.
1346*/
1347int
1348CbcOrbitalBranchingObject::compareOriginalObject
1349(const CbcBranchingObject* brObj) const
1350{
1351  const CbcOrbitalBranchingObject* br =
1352    dynamic_cast<const CbcOrbitalBranchingObject*>(brObj);
1353  assert(!br);
1354  abort();
1355  return 0;
1356}
1357
1358/** Compare the \c this with \c brObj. \c this and \c brObj must be os the
1359    same type and must have the same original object, but they may have
1360    different feasible regions.
1361    Return the appropriate CbcRangeCompare value (first argument being the
1362    sub/superset if that's the case). In case of overlap (and if \c
1363    replaceIfOverlap is true) replace the current branching object with one
1364    whose feasible region is the overlap.
1365*/
1366CbcRangeCompare
1367CbcOrbitalBranchingObject::compareBranchingObject
1368(const CbcBranchingObject* brObj, const bool replaceIfOverlap)
1369{
1370  const CbcOrbitalBranchingObject* br =
1371    dynamic_cast<const CbcOrbitalBranchingObject*>(brObj);
1372  assert(!br);
1373  abort();
1374  return CbcRangeDisjoint;
1375}
1376#endif
1377
Note: See TracBrowser for help on using the repository browser.