Improved clustergraph.h/cpp code and wrote unit tests
authorJoris Mooij <joris.mooij@tuebingen.mpg.de>
Wed, 7 Apr 2010 10:38:43 +0000 (12:38 +0200)
committerJoris Mooij <joris.mooij@tuebingen.mpg.de>
Wed, 7 Apr 2010 10:38:43 +0000 (12:38 +0200)
* Improved clustergraph.h/clustergraph.cpp:
  - Made (previously public) members G, vars and clusters private and added
    bipGraph(), vars() and clusters() methods which offer read-only access to
    these members.
  - Deprecated toVector()
  - Added nrVars() method
  - Renamed size() method into nrClusters()
  - Added var( size_t ) method
  - Added cluster( size_t ) method

ChangeLog
Makefile
include/dai/clustergraph.h
include/dai/prob.h
src/clustergraph.cpp
src/jtree.cpp
src/regiongraph.cpp
tests/unit/clustergraph.cpp [new file with mode: 0644]

index 1e3395d..3091044 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -3,6 +3,15 @@ git HEAD
 
 * Fixed some bugs in the MatLab interface build system
 * Fixed a bug in utils/fginfo.cpp
+* Improved clustergraph.h/clustergraph.cpp:
+  - Made (previously public) members G, vars and clusters private and added
+    bipGraph(), vars() and clusters() methods which offer read-only access to 
+    these members.
+  - Deprecated toVector()
+  - Added nrVars() method
+  - Renamed size() method into nrClusters()
+  - Added var( size_t ) method
+  - Added cluster( size_t ) method
 * Improved factorgraph.h/factorgraph.cpp:
   - FactorGraph::clamped() now contains a delta factor for the clamped variable
   - Renamed FactorGraph::Cliques() into FactorGraph::maximalFactorDomains()
@@ -10,7 +19,7 @@ git HEAD
   - Fixed bug in FactorGraph::clone()
   - FactorGraph::findVars( const VarSet& ) now returns a SmallSet<size_t>
     and its argument now has a const-qualifier (which fixes a small bug)
-  - Made previously public member G private and added the bipGraph() method,
+  - Made (previously public) member G private and added the bipGraph() method,
     which offers read-only access to it
 * Improved factor.h/cpp:
   - Fixed bug in min( const TFactor<T> &, const TFactor<T> & )
@@ -97,7 +106,7 @@ git HEAD
 * Compressed Makefile
 * Added unit tests for var, smallset, varset, graph, bipgraph,
   weightedgraph, enum, util, exceptions, properties, index, prob, 
-  factor, factorgraph
+  factor, factorgraph, clustergraph
 * Added unit testing framework
 * Added initialization of TRWBP weights by sampling spanning trees
 * Cleaned up MR code:
index d8e1f69..9a147e6 100644 (file)
--- a/Makefile
+++ b/Makefile
@@ -124,7 +124,7 @@ examples : examples/example$(EE) examples/example_bipgraph$(EE) examples/example
 
 matlabs : matlab/dai$(ME) matlab/dai_readfg$(ME) matlab/dai_writefg$(ME) matlab/dai_potstrength$(ME)
 
-unittests : tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE) tests/unit/graph$(EE) tests/unit/bipgraph$(EE) tests/unit/weightedgraph$(EE) tests/unit/enum$(EE) tests/unit/enum$(EE) tests/unit/util$(EE) tests/unit/exceptions$(EE) tests/unit/properties$(EE) tests/unit/index$(EE) tests/unit/prob$(EE) tests/unit/factor$(EE) tests/unit/factorgraph$(EE)
+unittests : tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE) tests/unit/graph$(EE) tests/unit/bipgraph$(EE) tests/unit/weightedgraph$(EE) tests/unit/enum$(EE) tests/unit/enum$(EE) tests/unit/util$(EE) tests/unit/exceptions$(EE) tests/unit/properties$(EE) tests/unit/index$(EE) tests/unit/prob$(EE) tests/unit/factor$(EE) tests/unit/factorgraph$(EE) tests/unit/clustergraph$(EE)
        echo Running unit tests...
        tests/unit/var$(EE)
        tests/unit/smallset$(EE)
@@ -140,6 +140,7 @@ unittests : tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE)
        tests/unit/prob$(EE)
        tests/unit/factor$(EE)
        tests/unit/factorgraph$(EE)
+       tests/unit/clustergraph$(EE)
 ifneq ($(OS),WINDOWS)
        rm factorgraph_test.fg
 else
@@ -292,7 +293,7 @@ clean :
        -rm matlab/*$(ME)
        -rm examples/example$(EE) examples/example_bipgraph$(EE) examples/example_varset$(EE) examples/example_permute$(EE) examples/example_sprinkler$(EE) examples/example_sprinkler_gibbs$(EE) examples/example_sprinkler_em$(EE)
        -rm tests/testdai$(EE) tests/testem/testem$(EE) tests/testbbp$(EE)
-       -rm tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE) tests/unit/graph$(EE) tests/unit/bipgraph$(EE) tests/unit/weightedgraph$(EE) tests/unit/enum$(EE) tests/unit/util$(EE) tests/unit/exceptions$(EE) tests/unit/properties$(EE) tests/unit/index$(EE) tests/unit/prob$(EE) tests/unit/factor$(EE) tests/unit/factorgraph$(EE)
+       -rm tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE) tests/unit/graph$(EE) tests/unit/bipgraph$(EE) tests/unit/weightedgraph$(EE) tests/unit/enum$(EE) tests/unit/util$(EE) tests/unit/exceptions$(EE) tests/unit/properties$(EE) tests/unit/index$(EE) tests/unit/prob$(EE) tests/unit/factor$(EE) tests/unit/factorgraph$(EE) tests/unit/clustergraph$(EE)
        -rm utils/fg2dot$(EE) utils/createfg$(EE) utils/fginfo$(EE)
        -rm -R doc
        -rm -R lib
@@ -337,6 +338,7 @@ clean :
        -del tests\unit\prob$(EE)
        -del tests\unit\factor$(EE)
        -del tests\unit\factorgraph$(EE)
+       -del tests\unit\clustergraph$(EE)
        -del $(LIB)\libdai$(LE)
        -rmdir lib
 endif
index 5f9845c..fb16303 100644 (file)
@@ -11,6 +11,8 @@
 
 /// \file
 /// \brief Defines class ClusterGraph, which is used by JTree, TreeEP and HAK
+/// \todo The "MinFill" and "WeightedMinFill" variable elimination heuristics seem designed for Markov graphs;
+/// add similar heuristics which are designed for factor graphs.
 
 
 #ifndef __defined_libdai_clustergraph_h
@@ -28,67 +30,105 @@ namespace dai {
 
     /// A ClusterGraph is a hypergraph with variables as nodes, and "clusters" (sets of variables) as hyperedges.
     /** It is implemented as a bipartite graph with variable (Var) nodes and cluster (VarSet) nodes.
+     *  One may think of a ClusterGraph as a FactorGraph without the actual factor values.
      */
     class ClusterGraph {
         public:
-            /// Stores the neighborhood structure
-            BipartiteGraph       G;
-
-            /// Stores the variables corresponding to the nodes
-            std::vector<Var>     vars;
-
-            /// Stores the clusters corresponding to the hyperedges
-            std::vector<VarSet>  clusters;
-
             /// Shorthand for BipartiteGraph::Neighbor
             typedef BipartiteGraph::Neighbor Neighbor;
 
             /// Shorthand for BipartiteGraph::Edge
             typedef BipartiteGraph::Edge     Edge;
 
+        private:
+            /// Stores the neighborhood structure
+            BipartiteGraph       _G;
+
+            /// Stores the variables corresponding to the nodes
+            std::vector<Var>     _vars;
+
+            /// Stores the clusters corresponding to the hyperedges
+            std::vector<VarSet>  _clusters;
+
         public:
         /// \name Constructors and destructors
         //@{
             /// Default constructor
-            ClusterGraph() : G(), vars(), clusters() {}
+            ClusterGraph() : _G(), _vars(), _clusters() {}
 
             /// Construct from vector of VarSet 's
-            ClusterGraph( const std::vector<VarSet> & cls );
+            ClusterGraph( const std::vector<VarSet>& cls );
         //@}
 
         /// \name Queries
         //@{
-            /// Returns a constant reference to the clusters
-            const std::vector<VarSet> & toVector() const { return clusters; }
+            /// Returns a constant reference to the graph structure
+            const BipartiteGraph& bipGraph() const { return _G; }
+
+            /// Returns number of variables
+            size_t nrVars() const { return _vars.size(); }
+
+            /// Returns a constant reference to the variables
+            const std::vector<Var>& vars() const { return _vars; }
+
+            /// Returns a constant reference to the \a i'th variable
+            const Var& var( size_t i ) const {
+                DAI_DEBASSERT( i < nrVars() );
+                return _vars[i]; 
+            }
 
             /// Returns number of clusters
-            size_t size() const {
-                return G.nrNodes2();
+            size_t nrClusters() const { return _clusters.size(); }
+
+            /// Returns a constant reference to the clusters
+            const std::vector<VarSet>& clusters() const { return _clusters; }
+
+            /// Returns a constant reference to the \a I'th cluster
+            const VarSet& cluster( size_t I ) const {
+                DAI_DEBASSERT( I < nrClusters() );
+                return _clusters[I]; 
             }
 
+            /// Returns a constant reference to the clusters
+            /** \deprecated Please use dai::ClusterGraph::clusters() instead
+             */
+            const std::vector<VarSet>& toVector() const { return _clusters; }
+
+            /// Returns number of clusters
+            /** \deprecated Please use dai::ClusterGraph::nrClusters() instead
+             */
+            size_t size() const { return _G.nrNodes2(); }
+
             /// Returns the index of variable \a n
-            size_t findVar( const Var &n ) const {
-                return find( vars.begin(), vars.end(), n ) - vars.begin();
+            /** \throw OBJECT_NOT_FOUND if the variable does not occur in the cluster graph
+             */
+            size_t findVar( const Var& n ) const {
+                size_t r = find( _vars.begin(), _vars.end(), n ) - _vars.begin();
+                if( r == _vars.size() )
+                    DAI_THROW(OBJECT_NOT_FOUND);
+                return r;
             }
 
             /// Returns union of clusters that contain the \a i 'th variable
             VarSet Delta( size_t i ) const {
                 VarSet result;
-                foreach( const Neighbor &I, G.nb1(i) )
-                    result |= clusters[I];
+                foreach( const Neighbor &I, _G.nb1(i) )
+                    result |= _clusters[I];
                 return result;
             }
 
             /// Returns union of clusters that contain the \a i 'th (except this variable itself)
             VarSet delta( size_t i ) const {
-                return Delta( i ) / vars[i];
+                return Delta( i ) / _vars[i];
             }
 
             /// Returns \c true if variables with indices \a i1 and \a i2 are adjacent, i.e., both contained in the same cluster
             bool adj( size_t i1, size_t i2 ) const {
+                if( i1 == i2 )
+                    return false;
                 bool result = false;
-                foreach( const Neighbor &I, G.nb1(i1) )
-                    if( find( G.nb2(I).begin(), G.nb2(I).end(), i2 ) != G.nb2(I).end() ) {
+                foreach( const Neighbor &I, _G.nb1(i1) )
+                    if( find( _G.nb2(I).begin(), _G.nb2(I).end(), i2 ) != _G.nb2(I).end() ) {
                         result = true;
                         break;
                     }
@@ -97,13 +137,13 @@ namespace dai {
 
             /// Returns \c true if cluster \a I is not contained in a larger cluster
             bool isMaximal( size_t I ) const {
-                DAI_DEBASSERT( I < G.nrNodes2() );
-                const VarSet & clI = clusters[I];
+                DAI_DEBASSERT( I < _G.nrNodes2() );
+                const VarSet & clI = _clusters[I];
                 bool maximal = true;
                 // The following may not be optimal, since it may repeatedly test the same cluster *J
-                foreach( const Neighbor &i, G.nb2(I) ) {
-                    foreach( const Neighbor &J, G.nb1(i) )
-                        if( (J != I) && (clI << clusters[J]) ) {
+                foreach( const Neighbor &i, _G.nb2(I) ) {
+                    foreach( const Neighbor &J, _G.nb1(i) )
+                        if( (J != I) && (clI << _clusters[J]) ) {
                             maximal = false;
                             break;
                         }
@@ -117,29 +157,29 @@ namespace dai {
         /// \name Operations
         //@{
             /// Inserts a cluster (if it does not already exist)
-            void insert( const VarSet &cl ) {
-                if( find( clusters.begin(), clusters.end(), cl ) == clusters.end() ) {
-                    clusters.push_back( cl );
+            void insert( const VarSetcl ) {
+                if( find( _clusters.begin(), _clusters.end(), cl ) == _clusters.end() ) {
+                    _clusters.push_back( cl );
                     // add variables (if necessary) and calculate neighborhood of new cluster
                     std::vector<size_t> nbs;
                     for( VarSet::const_iterator n = cl.begin(); n != cl.end(); n++ ) {
-                        size_t iter = find( vars.begin(), vars.end(), *n ) - vars.begin();
+                        size_t iter = find( _vars.begin(), _vars.end(), *n ) - _vars.begin();
                         nbs.push_back( iter );
-                        if( iter == vars.size() ) {
-                            G.addNode1();
-                            vars.push_back( *n );
+                        if( iter == _vars.size() ) {
+                            _G.addNode1();
+                            _vars.push_back( *n );
                         }
                     }
-                    G.addNode2( nbs.begin(), nbs.end(), nbs.size() );
+                    _G.addNode2( nbs.begin(), nbs.end(), nbs.size() );
                 }
             }
 
             /// Erases all clusters that are not maximal
             ClusterGraph& eraseNonMaximal() {
-                for( size_t I = 0; I < G.nrNodes2(); ) {
+                for( size_t I = 0; I < _G.nrNodes2(); ) {
                     if( !isMaximal(I) ) {
-                        clusters.erase( clusters.begin() + I );
-                        G.eraseNode2(I);
+                        _clusters.erase( _clusters.begin() + I );
+                        _G.eraseNode2(I);
                     } else
                         I++;
                 }
@@ -148,9 +188,9 @@ namespace dai {
 
             /// Erases all clusters that contain the \a i 'th variable
             ClusterGraph& eraseSubsuming( size_t i ) {
-                while( G.nb1(i).size() ) {
-                    clusters.erase( clusters.begin() + G.nb1(i)[0] );
-                    G.eraseNode2( G.nb1(i)[0] );
+                while( _G.nb1(i).size() ) {
+                    _clusters.erase( _clusters.begin() + _G.nb1(i)[0] );
+                    _G.eraseNode2( _G.nb1(i)[0] );
                 }
                 return *this;
             }
@@ -159,15 +199,15 @@ namespace dai {
         /// \name Input/Ouput
         //@{
             /// Writes a ClusterGraph to an output stream
-            friend std::ostream & operator << ( std::ostream & os, const ClusterGraph & cl ) {
-                os << cl.toVector();
+            friend std::ostream& operator << ( std::ostream& os, const ClusterGraph& cl ) {
+                os << cl.clusters();
                 return os;
             }
         //@}
 
         /// \name Variable elimination
         //@{
-            /// Performs Variable Elimination, only keeping track of the interactions that are created along the way.
+            /// Performs Variable Elimination, keeping track of the interactions that are created along the way.
             /** \tparam EliminationChoice should support "size_t operator()( const ClusterGraph &cl, const std::set<size_t> &remainingVars )"
              *  \param f function object which returns the next variable index to eliminate; for example, a dai::greedyVariableElimination object.
              *  \return A set of elimination "cliques".
@@ -182,13 +222,13 @@ namespace dai {
 
                 // Construct set of variable indices
                 std::set<size_t> varindices;
-                for( size_t i = 0; i < vars.size(); ++i )
+                for( size_t i = 0; i < _vars.size(); ++i )
                     varindices.insert( i );
 
                 // Do variable elimination
                 while( !varindices.empty() ) {
                     size_t i = f( cl, varindices );
-                    DAI_ASSERT( i < vars.size() );
+                    DAI_ASSERT( i < _vars.size() );
                     result.insert( cl.Delta( i ) );
                     cl.insert( cl.delta( i ) );
                     cl.eraseSubsuming( i );
@@ -243,7 +283,7 @@ namespace dai {
             /// Returns the best variable from \a remainingVars to eliminate in the cluster graph \a cl by greedily minimizing the cost function.
             /** This function calculates the cost for eliminating each variable in \a remaingVars and returns the variable which has lowest cost.
              */
-            size_t operator()( const ClusterGraph &cl, const std::set<size_t> &remainingVars );
+            size_t operator()( const ClusterGraph &cl, const std::set<size_t>remainingVars );
     };
 
 
@@ -252,7 +292,7 @@ namespace dai {
      *  where the adjacency graph has the variables as its nodes and connects
      *  nodes \a i1 and \a i2 iff \a i1 and \a i2 occur together in some common cluster.
      */
-    size_t eliminationCost_MinNeighbors( const ClusterGraph &cl, size_t i );
+    size_t eliminationCost_MinNeighbors( const ClusterGraphcl, size_t i );
 
 
     /// Calculates cost of eliminating the \a i 'th variable from cluster graph \a cl according to the "MinWeight" criterion.
@@ -261,7 +301,7 @@ namespace dai {
      *  nodes \a i1 and \a i2 iff \a i1 and \a i2 occur together in some common cluster.
      *  The weight of a node is the number of states of the corresponding variable.
      */
-    size_t eliminationCost_MinWeight( const ClusterGraph &cl, size_t i );
+    size_t eliminationCost_MinWeight( const ClusterGraphcl, size_t i );
 
 
     /// Calculates cost of eliminating the \a i 'th variable from cluster graph \a cl according to the "MinFill" criterion.
@@ -269,7 +309,7 @@ namespace dai {
      *  where the adjacency graph has the variables as its nodes and connects
      *  nodes \a i1 and \a i2 iff \a i1 and \a i2 occur together in some common cluster.
      */
-    size_t eliminationCost_MinFill( const ClusterGraph &cl, size_t i );
+    size_t eliminationCost_MinFill( const ClusterGraphcl, size_t i );
 
 
     /// Calculates cost of eliminating the \a i 'th variable from cluster graph \a cl according to the "WeightedMinFill" criterion.
@@ -278,7 +318,7 @@ namespace dai {
      *  nodes \a i1 and \a i2 iff \a i1 and \a i2 occur together in some common cluster.
      *  The weight of an edge is the product of the number of states of the variables corresponding with its nodes.
      */
-    size_t eliminationCost_WeightedMinFill( const ClusterGraph &cl, size_t i );
+    size_t eliminationCost_WeightedMinFill( const ClusterGraphcl, size_t i );
 
 
 } // end of namespace dai
index 7212f7a..9a72416 100644 (file)
@@ -198,6 +198,8 @@ class TProb {
     public:
         /// Type of data structure used for storing the values
         typedef std::vector<T> container_type;
+
+        /// Shorthand
         typedef TProb<T> this_type;
 
     private:
index 0953918..46afcc1 100644 (file)
@@ -22,26 +22,26 @@ namespace dai {
 using namespace std;
 
 
-ClusterGraph::ClusterGraph( const std::vector<VarSet> & cls ) : G(), vars(), clusters() {
+ClusterGraph::ClusterGraph( const std::vector<VarSet> & cls ) : _G(), _vars(), _clusters() {
     // construct vars, clusters and edge list
     vector<Edge> edges;
     foreach( const VarSet &cl, cls ) {
-        if( find( clusters.begin(), clusters.end(), cl ) == clusters.end() ) {
+        if( find( clusters().begin(), clusters().end(), cl ) == clusters().end() ) {
             // add cluster
-            size_t n2 = clusters.size();
-            clusters.push_back( cl );
+            size_t n2 = nrClusters();
+            _clusters.push_back( cl );
             for( VarSet::const_iterator n = cl.begin(); n != cl.end(); n++ ) {
-                size_t n1 = find( vars.begin(), vars.end(), *n ) - vars.begin();
-                if( n1 == vars.size() )
+                size_t n1 = find( vars().begin(), vars().end(), *n ) - vars().begin();
+                if( n1 == nrVars() )
                     // add variable
-                    vars.push_back( *n );
+                    _vars.push_back( *n );
                 edges.push_back( Edge( n1, n2 ) );
             }
         } // disregard duplicate clusters
     }
 
     // Create bipartite graph
-    G.construct( vars.size(), clusters.size(), edges.begin(), edges.end() );
+    _G.construct( nrVars(), nrClusters(), edges.begin(), edges.end() );
 }
 
 
@@ -65,23 +65,23 @@ size_t greedyVariableElimination::operator()( const ClusterGraph &cl, const std:
 
 
 size_t eliminationCost_MinNeighbors( const ClusterGraph &cl, size_t i ) {
-    return cl.G.delta1( i ).size();
+    return cl.bipGraph().delta1( i ).size();
 }
 
 
 size_t eliminationCost_MinWeight( const ClusterGraph &cl, size_t i ) {
-    SmallSet<size_t> id_n = cl.G.delta1( i );
+    SmallSet<size_t> id_n = cl.bipGraph().delta1( i );
     
     size_t cost = 1;
     for( SmallSet<size_t>::const_iterator it = id_n.begin(); it != id_n.end(); it++ )
-        cost *= cl.vars[*it].states();
+        cost *= cl.vars()[*it].states();
 
     return cost;
 }
 
 
 size_t eliminationCost_MinFill( const ClusterGraph &cl, size_t i ) {
-    SmallSet<size_t> id_n = cl.G.delta1( i );
+    SmallSet<size_t> id_n = cl.bipGraph().delta1( i );
 
     size_t cost = 0;
     // for each unordered pair {i1,i2} adjacent to n
@@ -98,7 +98,7 @@ size_t eliminationCost_MinFill( const ClusterGraph &cl, size_t i ) {
 
 
 size_t eliminationCost_WeightedMinFill( const ClusterGraph &cl, size_t i ) {
-    SmallSet<size_t> id_n = cl.G.delta1( i );
+    SmallSet<size_t> id_n = cl.bipGraph().delta1( i );
 
     size_t cost = 0;
     // for each unordered pair {i1,i2} adjacent to n
@@ -107,7 +107,7 @@ size_t eliminationCost_WeightedMinFill( const ClusterGraph &cl, size_t i ) {
             if( it1 != it2 ) {
                 // if i1 and i2 are not adjacent, eliminating n would make them adjacent
                 if( !cl.adj(*it1, *it2) )
-                    cost += cl.vars[*it1].states() * cl.vars[*it2].states();
+                    cost += cl.vars()[*it1].states() * cl.vars()[*it2].states();
             }
 
     return cost;
index 5c034f8..b132302 100644 (file)
@@ -101,7 +101,7 @@ JTree::JTree( const FactorGraph &fg, const PropertySet &opts, bool automatic ) :
             default:
                 DAI_THROW(UNKNOWN_ENUM_VALUE);
         }
-        vector<VarSet> ElimVec = _cg.VarElim( greedyVariableElimination( ec ) ).eraseNonMaximal().toVector();
+        vector<VarSet> ElimVec = _cg.VarElim( greedyVariableElimination( ec ) ).eraseNonMaximal().clusters();
         if( props.verbose >= 3 )
             cerr << "VarElim result: " << ElimVec << endl;
 
@@ -531,7 +531,7 @@ std::pair<size_t,double> boundTreewidth( const FactorGraph &fg, greedyVariableEl
     _cg.eraseNonMaximal();
 
     // Obtain elimination sequence
-    vector<VarSet> ElimVec = _cg.VarElim( greedyVariableElimination( fn ) ).eraseNonMaximal().toVector();
+    vector<VarSet> ElimVec = _cg.VarElim( greedyVariableElimination( fn ) ).eraseNonMaximal().clusters();
 
     // Calculate treewidth
     size_t treewidth = 0;
index 8707354..0755fd8 100644 (file)
@@ -63,9 +63,9 @@ void RegionGraph::constructCVM( const FactorGraph &fg, const std::vector<VarSet>
 
     // Create inner regions - first pass
     set<VarSet> betas;
-    for( size_t alpha = 0; alpha < cg.clusters.size(); alpha++ )
-        for( size_t alpha2 = alpha; (++alpha2) != cg.clusters.size(); ) {
-            VarSet intersection = cg.clusters[alpha] & cg.clusters[alpha2];
+    for( size_t alpha = 0; alpha < cg.nrClusters(); alpha++ )
+        for( size_t alpha2 = alpha; (++alpha2) != cg.nrClusters(); ) {
+            VarSet intersection = cg.cluster(alpha) & cg.cluster(alpha2);
             if( intersection.size() > 0 )
                 betas.insert( intersection );
         }
@@ -92,12 +92,12 @@ void RegionGraph::constructCVM( const FactorGraph &fg, const std::vector<VarSet>
     // Create edges
     vector<pair<size_t,size_t> > edges;
     for( size_t beta = 0; beta < irs.size(); beta++ )
-        for( size_t alpha = 0; alpha < cg.clusters.size(); alpha++ )
-            if( cg.clusters[alpha] >> irs[beta] )
+        for( size_t alpha = 0; alpha < cg.nrClusters(); alpha++ )
+            if( cg.cluster(alpha) >> irs[beta] )
                 edges.push_back( pair<size_t,size_t>(alpha,beta) );
 
     // Construct region graph
-    construct( fg, cg.clusters, irs, edges );
+    construct( fg, cg.clusters(), irs, edges );
 
     // Calculate counting numbers
     calcCountingNumbers();
diff --git a/tests/unit/clustergraph.cpp b/tests/unit/clustergraph.cpp
new file mode 100644 (file)
index 0000000..dcdecf1
--- /dev/null
@@ -0,0 +1,427 @@
+/*  This file is part of libDAI - http://www.libdai.org/
+ *
+ *  libDAI is licensed under the terms of the GNU General Public License version
+ *  2, or (at your option) any later version. libDAI is distributed without any
+ *  warranty. See the file COPYING for more details.
+ *
+ *  Copyright (C) 2010  Joris Mooij      [joris dot mooij at libdai dot org]
+ */
+
+
+#define BOOST_TEST_DYN_LINK
+
+
+#include <dai/clustergraph.h>
+#include <vector>
+#include <strstream>
+
+
+using namespace dai;
+
+
+const double tol = 1e-8;
+
+
+#define BOOST_TEST_MODULE ClusterGraphTest
+
+
+#include <boost/test/unit_test.hpp>
+#include <boost/test/floating_point_comparison.hpp>
+
+
+BOOST_AUTO_TEST_CASE( ConstructorsTest ) {
+    ClusterGraph G;
+    BOOST_CHECK_EQUAL( G.clusters(), std::vector<VarSet>() );
+    BOOST_CHECK( G.bipGraph() == BipartiteGraph() );
+    BOOST_CHECK_EQUAL( G.nrVars(), 0 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 0 );
+    BOOST_CHECK_THROW( G.var( 0 ), Exception );
+    BOOST_CHECK_THROW( G.cluster( 0 ), Exception );
+    BOOST_CHECK_THROW( G.findVar( Var( 0, 2 ) ), Exception );
+
+    Var v0( 0, 2 );
+    Var v1( 1, 3 );
+    Var v2( 2, 2 );
+    Var v3( 3, 4 );
+    std::vector<Var> vs;
+    vs.push_back( v0 );
+    vs.push_back( v1 );
+    vs.push_back( v2 );
+    vs.push_back( v3 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v03( v0, v3 );
+    VarSet v12( v1, v2 );
+    VarSet v13( v1, v3 );
+    VarSet v23( v2, v3 );
+    std::vector<VarSet> cl;
+    cl.push_back( v01 );
+    cl.push_back( v12 );
+    cl.push_back( v23 );
+    cl.push_back( v13 );
+    ClusterGraph G2( cl );
+    BOOST_CHECK_EQUAL( G2.nrVars(), 4 );
+    BOOST_CHECK_EQUAL( G2.nrClusters(), 4 );
+    BOOST_CHECK_EQUAL( G2.vars(), vs );
+    BOOST_CHECK_EQUAL( G2.clusters(), cl );
+    BOOST_CHECK_EQUAL( G2.findVar( v0 ), 0 );
+    BOOST_CHECK_EQUAL( G2.findVar( v1 ), 1 );
+    BOOST_CHECK_EQUAL( G2.findVar( v2 ), 2 );
+    BOOST_CHECK_EQUAL( G2.findVar( v3 ), 3 );
+
+    ClusterGraph Gb( G );
+    BOOST_CHECK( G.bipGraph() == Gb.bipGraph() );
+    BOOST_CHECK( G.vars() == Gb.vars() );
+    BOOST_CHECK( G.clusters() == Gb.clusters() );
+
+    ClusterGraph Gc = G;
+    BOOST_CHECK( G.bipGraph() == Gc.bipGraph() );
+    BOOST_CHECK( G.vars() == Gc.vars() );
+    BOOST_CHECK( G.clusters() == Gc.clusters() );
+
+    ClusterGraph G2b( G2 );
+    BOOST_CHECK( G2.bipGraph() == G2b.bipGraph() );
+    BOOST_CHECK( G2.vars() == G2b.vars() );
+    BOOST_CHECK( G2.clusters() == G2b.clusters() );
+
+    ClusterGraph G2c = G2;
+    BOOST_CHECK( G2.bipGraph() == G2c.bipGraph() );
+    BOOST_CHECK( G2.vars() == G2c.vars() );
+    BOOST_CHECK( G2.clusters() == G2c.clusters() );
+}
+
+
+BOOST_AUTO_TEST_CASE( QueriesTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 3 );
+    Var v2( 2, 2 );
+    Var v3( 3, 4 );
+    Var v4( 4, 2 );
+    std::vector<Var> vs;
+    vs.push_back( v0 );
+    vs.push_back( v1 );
+    vs.push_back( v2 );
+    vs.push_back( v3 );
+    vs.push_back( v4 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v03( v0, v3 );
+    VarSet v04( v0, v4 );
+    VarSet v12( v1, v2 );
+    VarSet v13( v1, v3 );
+    VarSet v14( v1, v4 );
+    VarSet v23( v2, v3 );
+    VarSet v24( v2, v4 );
+    VarSet v34( v3, v4 );
+    VarSet v123 = v12 | v3;
+    std::vector<VarSet> cl;
+    cl.push_back( v01 );
+    cl.push_back( v12 );
+    cl.push_back( v123 );
+    cl.push_back( v34 );
+    cl.push_back( v04 );
+    ClusterGraph G( cl );
+
+    BOOST_CHECK_EQUAL( G.nrVars(), 5 );
+    BOOST_CHECK_EQUAL( G.vars(), vs );
+    BOOST_CHECK_EQUAL( G.var(0), v0 );
+    BOOST_CHECK_EQUAL( G.var(1), v1 );
+    BOOST_CHECK_EQUAL( G.var(2), v2 );
+    BOOST_CHECK_EQUAL( G.var(3), v3 );
+    BOOST_CHECK_EQUAL( G.var(4), v4 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 5 );
+    BOOST_CHECK_EQUAL( G.clusters(), cl );
+    BOOST_CHECK_EQUAL( G.cluster(0), v01 );
+    BOOST_CHECK_EQUAL( G.cluster(1), v12 );
+    BOOST_CHECK_EQUAL( G.cluster(2), v123 );
+    BOOST_CHECK_EQUAL( G.cluster(3), v34 );
+    BOOST_CHECK_EQUAL( G.cluster(4), v04 );
+    BOOST_CHECK_EQUAL( G.findVar( v0 ), 0 );
+    BOOST_CHECK_EQUAL( G.findVar( v1 ), 1 );
+    BOOST_CHECK_EQUAL( G.findVar( v2 ), 2 );
+    BOOST_CHECK_EQUAL( G.findVar( v3 ), 3 );
+    BOOST_CHECK_EQUAL( G.findVar( v4 ), 4 );
+    BipartiteGraph H( 5, 5 );
+    H.addEdge( 0, 0 );
+    H.addEdge( 1, 0 );
+    H.addEdge( 1, 1 );
+    H.addEdge( 2, 1 );
+    H.addEdge( 1, 2 );
+    H.addEdge( 2, 2 );
+    H.addEdge( 3, 2 );
+    H.addEdge( 3, 3 );
+    H.addEdge( 4, 3 );
+    H.addEdge( 0, 4 );
+    H.addEdge( 4, 4 );
+    BOOST_CHECK( G.bipGraph() == H );
+
+    BOOST_CHECK_EQUAL( G.delta( 0 ), v14 );
+    BOOST_CHECK_EQUAL( G.delta( 1 ), v02 | v3 );
+    BOOST_CHECK_EQUAL( G.delta( 2 ), v13 );
+    BOOST_CHECK_EQUAL( G.delta( 3 ), v12 | v4 );
+    BOOST_CHECK_EQUAL( G.delta( 4 ), v03 );
+    BOOST_CHECK_EQUAL( G.Delta( 0 ), v14 | v0 );
+    BOOST_CHECK_EQUAL( G.Delta( 1 ), v01 | v23 );
+    BOOST_CHECK_EQUAL( G.Delta( 2 ), v13 | v2 );
+    BOOST_CHECK_EQUAL( G.Delta( 3 ), v12 | v34 );
+    BOOST_CHECK_EQUAL( G.Delta( 4 ), v03 | v4 );
+
+    BOOST_CHECK( !G.adj( 0, 0 ) );
+    BOOST_CHECK(  G.adj( 0, 1 ) );
+    BOOST_CHECK( !G.adj( 0, 2 ) );
+    BOOST_CHECK( !G.adj( 0, 3 ) );
+    BOOST_CHECK(  G.adj( 0, 4 ) );
+    BOOST_CHECK(  G.adj( 1, 0 ) );
+    BOOST_CHECK( !G.adj( 1, 1 ) );
+    BOOST_CHECK(  G.adj( 1, 2 ) );
+    BOOST_CHECK(  G.adj( 1, 3 ) );
+    BOOST_CHECK( !G.adj( 1, 4 ) );
+    BOOST_CHECK( !G.adj( 2, 0 ) );
+    BOOST_CHECK(  G.adj( 2, 1 ) );
+    BOOST_CHECK( !G.adj( 2, 2 ) );
+    BOOST_CHECK(  G.adj( 2, 3 ) );
+    BOOST_CHECK( !G.adj( 2, 4 ) );
+    BOOST_CHECK( !G.adj( 3, 0 ) );
+    BOOST_CHECK(  G.adj( 3, 1 ) );
+    BOOST_CHECK(  G.adj( 3, 2 ) );
+    BOOST_CHECK( !G.adj( 3, 3 ) );
+    BOOST_CHECK(  G.adj( 3, 4 ) );
+    BOOST_CHECK(  G.adj( 4, 0 ) );
+    BOOST_CHECK( !G.adj( 4, 1 ) );
+    BOOST_CHECK( !G.adj( 4, 2 ) );
+    BOOST_CHECK(  G.adj( 4, 3 ) );
+    BOOST_CHECK( !G.adj( 4, 4 ) );
+
+    BOOST_CHECK(  G.isMaximal( 0 ) );
+    BOOST_CHECK( !G.isMaximal( 1 ) );
+    BOOST_CHECK(  G.isMaximal( 2 ) );
+    BOOST_CHECK(  G.isMaximal( 3 ) );
+    BOOST_CHECK(  G.isMaximal( 4 ) );
+}
+
+
+BOOST_AUTO_TEST_CASE( OperationsTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 3 );
+    Var v2( 2, 2 );
+    Var v3( 3, 4 );
+    Var v4( 4, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v03( v0, v3 );
+    VarSet v04( v0, v4 );
+    VarSet v12( v1, v2 );
+    VarSet v13( v1, v3 );
+    VarSet v14( v1, v4 );
+    VarSet v23( v2, v3 );
+    VarSet v24( v2, v4 );
+    VarSet v34( v3, v4 );
+    VarSet v123 = v12 | v3;
+    std::vector<VarSet> cl;
+    cl.push_back( v01 );
+    cl.push_back( v12 );
+    cl.push_back( v123 );
+    cl.push_back( v34 );
+    cl.push_back( v04 );
+    ClusterGraph G( cl );
+
+    BipartiteGraph H( 5, 5 );
+    H.addEdge( 0, 0 );
+    H.addEdge( 1, 0 );
+    H.addEdge( 1, 1 );
+    H.addEdge( 2, 1 );
+    H.addEdge( 1, 2 );
+    H.addEdge( 2, 2 );
+    H.addEdge( 3, 2 );
+    H.addEdge( 3, 3 );
+    H.addEdge( 4, 3 );
+    H.addEdge( 0, 4 );
+    H.addEdge( 4, 4 );
+    BOOST_CHECK( G.bipGraph() == H );
+
+    G.eraseNonMaximal();
+    BOOST_CHECK_EQUAL( G.nrClusters(), 4 );
+    H.eraseNode2( 1 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.eraseSubsuming( 4 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 2 );
+    H.eraseNode2( 2 );
+    H.eraseNode2( 2 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.insert( v34 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 3 );
+    G.insert( v123 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 3 );
+    H.addNode2();
+    H.addEdge( 3, 2 );
+    H.addEdge( 4, 2 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.insert( v12 );
+    G.insert( v23 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 5 );
+    H.addNode2();
+    H.addNode2();
+    H.addEdge( 1, 3 );
+    H.addEdge( 2, 3 );
+    H.addEdge( 2, 4 );
+    H.addEdge( 3, 4 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.eraseNonMaximal();
+    BOOST_CHECK_EQUAL( G.nrClusters(), 3 );
+    H.eraseNode2( 3 );
+    H.eraseNode2( 3 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.eraseSubsuming( 2 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 2 );
+    H.eraseNode2( 1 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.eraseNonMaximal();
+    BOOST_CHECK_EQUAL( G.nrClusters(), 2 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.eraseSubsuming( 0 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 1 );
+    H.eraseNode2( 0 );
+    BOOST_CHECK( G.bipGraph() == H );
+    G.eraseSubsuming( 4 );
+    BOOST_CHECK_EQUAL( G.nrClusters(), 0 );
+    H.eraseNode2( 0 );
+    BOOST_CHECK( G.bipGraph() == H );
+}
+
+
+BOOST_AUTO_TEST_CASE( VarElimTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 3 );
+    Var v2( 2, 2 );
+    Var v3( 3, 4 );
+    Var v4( 4, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v03( v0, v3 );
+    VarSet v04( v0, v4 );
+    VarSet v12( v1, v2 );
+    VarSet v13( v1, v3 );
+    VarSet v14( v1, v4 );
+    VarSet v23( v2, v3 );
+    VarSet v24( v2, v4 );
+    VarSet v34( v3, v4 );
+    VarSet v123 = v12 | v3;
+    std::vector<VarSet> cl;
+    cl.push_back( v01 );
+    cl.push_back( v12 );
+    cl.push_back( v123 );
+    cl.push_back( v34 );
+    cl.push_back( v04 );
+    ClusterGraph G( cl );
+    ClusterGraph Gorg = G;
+
+    BipartiteGraph H( 5, 5 );
+    H.addEdge( 0, 0 );
+    H.addEdge( 1, 0 );
+    H.addEdge( 1, 1 );
+    H.addEdge( 2, 1 );
+    H.addEdge( 1, 2 );
+    H.addEdge( 2, 2 );
+    H.addEdge( 3, 2 );
+    H.addEdge( 3, 3 );
+    H.addEdge( 4, 3 );
+    H.addEdge( 0, 4 );
+    H.addEdge( 4, 4 );
+    BOOST_CHECK( G.bipGraph() == H );
+
+    BOOST_CHECK_EQUAL( eliminationCost_MinFill( G, 0 ), 1 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinFill( G, 1 ), 2 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinFill( G, 2 ), 0 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinFill( G, 3 ), 2 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinFill( G, 4 ), 1 );
+    cl.clear();
+    cl.push_back( v123 );
+    cl.push_back( v01 | v4 );
+    cl.push_back( v13 | v4 );
+    cl.push_back( v34 );
+    cl.push_back( v4 );
+    BOOST_CHECK_EQUAL( G.VarElim( greedyVariableElimination( eliminationCost_MinFill ) ).clusters(), cl );
+
+    G = Gorg;
+    BOOST_CHECK_EQUAL( eliminationCost_WeightedMinFill( G, 0 ), 2*3 );
+    BOOST_CHECK_EQUAL( eliminationCost_WeightedMinFill( G, 1 ), 2*2+2*4 );
+    BOOST_CHECK_EQUAL( eliminationCost_WeightedMinFill( G, 2 ), 0 );
+    BOOST_CHECK_EQUAL( eliminationCost_WeightedMinFill( G, 3 ), 3*2+2*2 );
+    BOOST_CHECK_EQUAL( eliminationCost_WeightedMinFill( G, 4 ), 2*4 );
+    cl.clear();
+    cl.push_back( v123 );
+    cl.push_back( v01 | v4 );
+    cl.push_back( v13 | v4 );
+    cl.push_back( v34 );
+    cl.push_back( v4 );
+    BOOST_CHECK_EQUAL( G.VarElim( greedyVariableElimination( eliminationCost_WeightedMinFill ) ).clusters(), cl );
+
+    G = Gorg;
+    BOOST_CHECK_EQUAL( eliminationCost_MinNeighbors( G, 0 ), 2 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinNeighbors( G, 1 ), 3 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinNeighbors( G, 2 ), 2 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinNeighbors( G, 3 ), 3 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinNeighbors( G, 4 ), 2 );
+    cl.clear();
+    cl.push_back( v01 | v4 );
+    cl.push_back( v123 );
+    cl.push_back( v13 | v4 );
+    cl.push_back( v34 );
+    cl.push_back( v4 );
+    BOOST_CHECK_EQUAL( G.VarElim( greedyVariableElimination( eliminationCost_MinNeighbors ) ).clusters(), cl );
+
+    G = Gorg;
+    BOOST_CHECK_EQUAL( eliminationCost_MinWeight( G, 0 ), 3*2 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinWeight( G, 1 ), 2*2*4 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinWeight( G, 2 ), 3*4 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinWeight( G, 3 ), 3*2*2 );
+    BOOST_CHECK_EQUAL( eliminationCost_MinWeight( G, 4 ), 2*4 );
+    cl.clear();
+    cl.push_back( v01 | v4 );
+    cl.push_back( v123 );
+    cl.push_back( v13 | v4 );
+    cl.push_back( v14 );
+    cl.push_back( v4 );
+    BOOST_CHECK_EQUAL( G.VarElim( greedyVariableElimination( eliminationCost_MinWeight ) ).clusters(), cl );
+
+    G = Gorg;
+    std::vector<Var> vs;
+    vs.push_back( v4 );
+    vs.push_back( v3 );
+    vs.push_back( v2 );
+    vs.push_back( v1 );
+    vs.push_back( v0 );
+    cl.clear();
+    cl.push_back( v03 | v4 );
+    cl.push_back( v01 | v23 );
+    cl.push_back( v01 | v2 );
+    cl.push_back( v01 );
+    cl.push_back( v0 );
+    BOOST_CHECK_EQUAL( G.VarElim( sequentialVariableElimination( vs ) ).clusters(), cl );
+}
+
+
+BOOST_AUTO_TEST_CASE( IOTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 3 );
+    Var v2( 2, 2 );
+    Var v3( 3, 4 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v03( v0, v3 );
+    VarSet v12( v1, v2 );
+    VarSet v13( v1, v3 );
+    VarSet v23( v2, v3 );
+    std::vector<VarSet> cl;
+    cl.push_back( v01 );
+    cl.push_back( v12 );
+    cl.push_back( v23 );
+    cl.push_back( v13 );
+    ClusterGraph G( cl );
+
+    std::stringstream ss;
+    ss << G;
+    std::string s;
+    getline( ss, s );
+    BOOST_CHECK_EQUAL( s, "({x0, x1}, {x1, x2}, {x2, x3}, {x1, x3})" );
+}