1 /* $Id: CbcSymmetry.hpp 1033 2013-12-14 19:34:28Z pbelotti $
2  *
3  * Name:    Hacked from CouenneProblem.hpp
4  * Author:  Pietro Belotti, Lehigh University
5  *          Andreas Waechter, IBM
6  * Purpose: define the class CouenneProblem
7  *
8  * (C) Carnegie-Mellon University, 2006-11.
9  * This file is licensed under the Eclipse Public License (EPL)
10  */
11 /*
12   If this is much used then we could improve build experience
13   Download nauty - say to /disk/nauty25r9
14   In that directory ./configure --enable-tls --enable-wordsize=32
15   make
16   copy nauty.a to libnauty.a
17 
18   In Cbc's configure
19   add -DCOIN_HAS_NTY to CXXDEFS
20   add -I/disk/nauty25r9 to CXXDEFS or ADD_CXXFLAGS
21   add -L/disk/nauty25r9 -lnauty to LDFLAGS
22 
23   If you wish to use Traces rather than nauty then add -DNTY_TRACES
24 
25   To use it is -orbit on
26 
27   Nauty has this -
28 *   Permission
29 *   is hereby given for use and/or distribution with the exception of        *
30 *   sale for profit or application with nontrivial military significance.    *
31   so you can use it internally even if you are a company.
32  */
33 #ifndef CBC_SYMMETRY_HPP
34 #define CBC_SYMMETRY_HPP
35 extern "C" {
36 #include "nauty.h"
37 #include "nausparse.h"
38 #ifdef NTY_TRACES
39 #include "traces.h"
40 #endif
41 }
42 
43 #include <vector>
44 #include <map>
45 #include <string.h>
46 
47 #include "CbcModel.hpp"
48 
49 class OsiObject;
50 // when to give up (depth since last success)
51 #ifndef NTY_BAD_DEPTH
52 #define NTY_BAD_DEPTH 10
53 #endif
54 class CbcNauty;
55 
56 #define COUENNE_HACKED_EPS 1.e-07
57 #define COUENNE_HACKED_EPS_SYMM 1e-8
58 #define COUENNE_HACKED_EXPRGROUP 8
59 
60 /** Class to deal with symmetry
61  *
62  *  Hacked from Couenne
63  *  Thanks, but it had been nice to make sure that there are no symbol collisions when building Couenne with this Cbc.
64  */
65 
66 class CbcSymmetry {
67   class Node {
68     int index;
69     double coeff;
70     double lb;
71     double ub;
72     int color;
73     int code;
74     int sign;
75 
76   public:
77     void node(int, double, double, double, int, int);
color_vertex(register int k)78     inline void color_vertex(register int k) { color = k; }
get_index() const79     inline int get_index() const { return index; }
get_coeff() const80     inline double get_coeff() const { return coeff; }
get_lb() const81     inline double get_lb() const { return lb; }
get_ub() const82     inline double get_ub() const { return ub; }
get_color() const83     inline int get_color() const { return color; }
get_code() const84     inline int get_code() const { return code; }
get_sign() const85     inline int get_sign() const { return sign; }
bounds(register double a,register double b)86     inline void bounds(register double a, register double b)
87     {
88       lb = a;
89       ub = b;
90     }
91   };
92 
93   struct myclass0 {
operator ()CbcSymmetry::myclass094     inline bool operator()(register const Node &a, register const Node &b)
95     {
96 
97       return ((a.get_code() < b.get_code()) || ((a.get_code() == b.get_code() && ((a.get_coeff() < b.get_coeff() - COUENNE_HACKED_EPS_SYMM) || ((fabs(a.get_coeff() - b.get_coeff()) < COUENNE_HACKED_EPS_SYMM) && ((a.get_lb() < b.get_lb() - COUENNE_HACKED_EPS_SYMM) || ((fabs(a.get_lb() - b.get_lb()) < COUENNE_HACKED_EPS_SYMM) && ((a.get_ub() < b.get_ub() - COUENNE_HACKED_EPS_SYMM) || ((fabs(a.get_ub() - b.get_ub()) < COUENNE_HACKED_EPS_SYMM) && ((a.get_index() < b.get_index())))))))))));
98     }
99   };
100 
101   struct myclass {
operator ()CbcSymmetry::myclass102     inline bool operator()(register const Node &a, register const Node &b)
103     {
104       return (a.get_index() < b.get_index());
105     }
106   };
107 
108   struct less_than_str {
operator ()CbcSymmetry::less_than_str109     inline bool operator()(register const char *a, register const char *b) const
110     {
111       return strcmp(a, b) < 0;
112     }
113   };
114 
115 public:
116   /**@name Constructors and destructors */
117   //@{
118   /// Default constructor
119   CbcSymmetry();
120 
121   /// Copy constructor
122   CbcSymmetry(const CbcSymmetry &);
123 
124   /// Assignment operator
125   CbcSymmetry &operator=(const CbcSymmetry &rhs);
126 
127   /// Destructor
128   ~CbcSymmetry();
129   //@}
130 
131   // Symmetry Info
132 
133   std::vector< int > *Find_Orbit(int) const;
134 
135   myclass0 node_sort;
136   myclass index_sort;
137 
138   void Compute_Symmetry() const;
139   int statsOrbits(CbcModel *model, int type) const;
140   //double timeNauty () const;
141   void Print_Orbits() const;
142   void fillOrbits();
143   /// Fixes variables using orbits (returns number fixed)
144   int orbitalFixing(OsiSolverInterface *solver);
whichOrbit()145   inline int *whichOrbit()
146   {
147     return numberUsefulOrbits_ ? whichOrbit_ : NULL;
148   }
numberUsefulOrbits() const149   inline int numberUsefulOrbits() const
150   {
151     return numberUsefulOrbits_;
152   }
numberUsefulObjects() const153   inline int numberUsefulObjects() const
154   {
155     return numberUsefulObjects_;
156   }
157   int largestOrbit(const double *lower, const double *upper) const;
158   void ChangeBounds(const double *lower, const double *upper,
159     int numberColumns, bool justFixedAtOne) const;
160   inline bool compare(register Node &a, register Node &b) const;
getNtyInfo()161   CbcNauty *getNtyInfo() { return nauty_info_; }
162 
163   // bool node_sort (  Node  a, Node  b);
164   // bool index_sort (  Node  a, Node  b);
165 
166   /// empty if no NTY, symmetry data structure setup otherwise
167   void setupSymmetry(CbcModel * model);
168 
169 private:
170   mutable std::vector< Node > node_info_;
171   mutable CbcNauty *nauty_info_;
172   int numberColumns_;
173   int numberUsefulOrbits_;
174   int numberUsefulObjects_;
175   int *whichOrbit_;
176 };
177 
178 class CbcNauty {
179 
180 public:
181   enum VarStatus { FIX_AT_ZERO,
182     FIX_AT_ONE,
183     FREE };
184 
185   /**@name Constructors and destructors */
186   //@{
187 private:
188   /// Default constructor
189   CbcNauty();
190 
191 public:
192   /// Normal constructor (if dense - NULLS)
193   CbcNauty(int n, const size_t *v, const int *d, const int *e);
194 
195   /// Copy constructor
196   CbcNauty(const CbcNauty &);
197 
198   /// Assignment operator
199   CbcNauty &operator=(const CbcNauty &rhs);
200 
201   /// Destructor
202   ~CbcNauty();
203   //@}
204 
205   void addElement(int ix, int jx);
206   void clearPartitions();
207   void computeAuto();
208   void deleteElement(int ix, int jx);
color_node(int ix,int color)209   void color_node(int ix, int color) { vstat_[ix] = color; }
insertRHS(int rhs,int cons)210   void insertRHS(int rhs, int cons) { constr_rhs.insert(std::pair< int, int >(rhs, cons)); }
211 
212   double getGroupSize() const;
213   //int getNautyCalls() const { return nautyCalls_; }
214   //double getNautyTime() const { return nautyTime_; }
215 
getN() const216   int getN() const { return n_; }
217 
218   int getNumGenerators() const;
219   int getNumOrbits() const;
220 
221   /// Returns the orbits in a "convenient" form
222   std::vector< std::vector< int > > *getOrbits() const;
223 
224   void getVstat(double *v, int nv);
isSparse() const225   inline bool isSparse() const
226   {
227     return GSparse_ != NULL;
228   }
errorStatus() const229   inline int errorStatus() const
230 #ifndef NTY_TRACES
231   {
232     return stats_->errstatus;
233   }
234 #else
235   {
236     return 0;
237   }
238 #endif
239   /// Pointer to options
options() const240   inline optionblk *options() const
241   {
242     return options_;
243   }
244   /**
245    * Methods to classify orbits.  Not horribly efficient, but gets the job done
246    */
247   //  bool isAllFixOneOrbit(const std::vector<int> &orbit) const;
248   // bool isAllFreeOrbit(const std::vector<int> &orbit) const;
249   //bool isAutoComputed() const { return autoComputed_; }
250   //bool isConstraintOrbit(const std::vector<int> &orbit) const;
251   //bool isMixedFreeZeroOrbit(const std::vector<int> &orbit) const;
252   //void makeFree(int ix) { vstat_[ix] = FREE; }
253 
254   void setWriteAutoms(const std::string &afilename);
255   void unsetWriteAutoms();
256 
257 private:
258   // The base nauty stuff
259   graph *G_;
260   sparsegraph *GSparse_;
261   int *lab_;
262   int *ptn_;
263   set *active_;
264   int *orbits_;
265 #ifndef NTY_TRACES
266   optionblk *options_;
267   statsblk *stats_;
268 #else
269   TracesOptions *options_;
270   TracesStats *stats_;
271 #endif
272   setword *workspace_;
273   int worksize_;
274   int m_;
275   int n_;
276   size_t nel_;
277   graph *canonG_;
278 
279   bool autoComputed_;
280 
281   int *vstat_;
282 
283   //static int nautyCalls_;
284   //static double nautyTime_;
285 
286   std::multimap< int, int > constr_rhs;
287   std::multimap< int, int >::iterator it;
288 
289   std::pair< std::multimap< int, int >::iterator,
290     std::multimap< int, int >::iterator >
291     ret;
292 
293   // File pointer for automorphism group
294   FILE *afp_;
295 };
296 
297 /** Branching object for Orbital branching
298 
299     Variable_ is the set id number (redundant, as the object also holds a
300     pointer to the set.
301  */
302 class CbcOrbitalBranchingObject : public CbcBranchingObject {
303 
304 public:
305   // Default Constructor
306   CbcOrbitalBranchingObject();
307 
308   // Useful constructor
309   CbcOrbitalBranchingObject(CbcModel *model, int column,
310     int way,
311     int numberExtra, const int *extraToZero);
312 
313   // Copy constructor
314   CbcOrbitalBranchingObject(const CbcOrbitalBranchingObject &);
315 
316   // Assignment operator
317   CbcOrbitalBranchingObject &operator=(const CbcOrbitalBranchingObject &rhs);
318 
319   /// Clone
320   virtual CbcBranchingObject *clone() const;
321 
322   // Destructor
323   virtual ~CbcOrbitalBranchingObject();
324 
325   using CbcBranchingObject::branch;
326   /// Does next branch and updates state
327   virtual double branch();
328   /** Update bounds in solver as in 'branch' and update given bounds.
329         branchState is -1 for 'down' +1 for 'up' */
330   virtual void fix(OsiSolverInterface *solver,
331     double *lower, double *upper,
332     int branchState) const;
333 
334   /** Reset every information so that the branching object appears to point to
335         the previous child. This method does not need to modify anything in any
336         solver. */
previousBranch()337   virtual void previousBranch()
338   {
339     CbcBranchingObject::previousBranch();
340   }
341 
342   using CbcBranchingObject::print;
343   /** \brief Print something about branch - only if log level high
344     */
345   virtual void print();
346 
347   /** Return the type (an integer identifier) of \c this */
type() const348   virtual CbcBranchObjType type() const
349   {
350     return SoSBranchObj;
351   }
352 
353   /** Compare the original object of \c this with the original object of \c
354         brObj. Assumes that there is an ordering of the original objects.
355         This method should be invoked only if \c this and brObj are of the same
356         type.
357         Return negative/0/positive depending on whether \c this is
358         smaller/same/larger than the argument.
359     */
360   virtual int compareOriginalObject(const CbcBranchingObject *brObj) const;
361 
362   /** Compare the \c this with \c brObj. \c this and \c brObj must be os the
363         same type and must have the same original object, but they may have
364         different feasible regions.
365         Return the appropriate CbcRangeCompare value (first argument being the
366         sub/superset if that's the case). In case of overlap (and if \c
367         replaceIfOverlap is true) replace the current branching object with one
368         whose feasible region is the overlap.
369      */
370   virtual CbcRangeCompare compareBranchingObject(const CbcBranchingObject *brObj, const bool replaceIfOverlap = false);
371 
372 private:
373   /// Column to go to 1
374   int column_;
375   /// Number (without column) going to zero on down branch
376   int numberOther_;
377   /// Number extra
378   int numberExtra_;
379   /// Fix to zero
380   int *fixToZero_;
381 };
382 
383 #endif
384 
385 /* vi: softtabstop=2 shiftwidth=2 expandtab tabstop=2
386 */
387