Some improvements to jtree and regiongraph and started work on regiongraph unit tests
authorJoris Mooij <joris.mooij@tuebingen.mpg.de>
Wed, 7 Apr 2010 18:40:38 +0000 (20:40 +0200)
committerJoris Mooij <joris.mooij@tuebingen.mpg.de>
Wed, 7 Apr 2010 18:40:38 +0000 (20:40 +0200)
* Improved jtree.h/cpp:
  - changed JTree::construct( const std::vector<VarSet>&, bool ) into
    JTree::construct( const FactorGraph&, const std::vector<VarSet>&, bool )
  - changed JTree::GenerateJT( const std::vector<VarSet> & )
    into JTree::GenerateJT( const FactorGraph &, const std::vector<VarSet> & )
* Improved regiongraph.h/cpp:
  - Made (previously public) members RegionGraph::G, RegionGraph::ORs,
    RegionGraph::IRs and RegionGraph::fac2OR protected.
  - Removed partial constructor RegionGraph::RegionGraph( const FactorGraph& fg )
  - Added some error checking code

12 files changed:
ChangeLog
Makefile
include/dai/jtree.h
include/dai/regiongraph.h
include/dai/smallset.h
include/dai/treeep.h
src/hak.cpp
src/jtree.cpp
src/regiongraph.cpp
src/treeep.cpp
tests/unit/factor.cpp
tests/unit/regiongraph.cpp [new file with mode: 0644]

index 3091044..c9acbc3 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -3,24 +3,35 @@ 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:
+* Improved jtree.h/cpp:
+  - changed JTree::construct( const std::vector<VarSet>&, bool ) into
+    JTree::construct( const FactorGraph&, const std::vector<VarSet>&, bool )
+  - changed JTree::GenerateJT( const std::vector<VarSet> & )
+    into JTree::GenerateJT( const FactorGraph &, const std::vector<VarSet> & )
+* Improved regiongraph.h/cpp:
+  - Made (previously public) members RegionGraph::G, RegionGraph::ORs,
+    RegionGraph::IRs and RegionGraph::fac2OR protected.
+  - Removed partial constructor RegionGraph::RegionGraph( const FactorGraph& fg )
+  - Added some error checking code
+* Improved clustergraph.h/cpp:
+  - Made (previously public) members ClusterGraph::G, ClusterGraph::vars and 
+    ClusterGraph::clusters private and added ClusterGraph::bipGraph(), 
+    ClusterGraph::vars() and ClusterGraph::clusters() methods which offer 
+    read-only access to these members.
+  - Deprecated ClusterGraph::toVector()
+  - Added ClusterGraph::nrVars() method
+  - Renamed ClusterGraph::size() method into ClusterGraph::nrClusters()
+  - Added ClusterGraph::var( size_t ) method
+  - Added ClusterGraph::cluster( size_t ) method
+* Improved factorgraph.h/cpp:
   - FactorGraph::clamped() now contains a delta factor for the clamped variable
   - Renamed FactorGraph::Cliques() into FactorGraph::maximalFactorDomains()
   - Added FactorGraph::MarkovGraph()
   - 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,
-    which offers read-only access to it
+    and its argument now has a const-qualifier (which fixes a small bug).
+  - Made (previously public) member FactorGraph::G private and added the 
+    FactorGraph::bipGraph() method, which offers read-only access to it.
 * Improved factor.h/cpp:
   - Fixed bug in min( const TFactor<T> &, const TFactor<T> & )
   - Fixed bug in max( const TFactor<T> &, const TFactor<T> & )
@@ -31,8 +42,8 @@ git HEAD
   - Added TFactor<T>::sumAbs() const
   - Added TFactor<T>::operator==( const TFactor<T> &y )
   - Renamed TFactor<T>::states() into TFactor<T>::nrStates()
-  - Added get(size_t) and set(size_t,T) operators; get() is 
-    equivalent to "operator[](size_t) const" and set() should
+  - Added TFactor<T>::get(size_t) and TFactor<T>::set(size_t,T) operators; 
+    get() is equivalent to "operator[](size_t) const" and set() should
     be used instead of the non-const operator[], which has been deprecated
 * Improved prob.h/cpp:
   - Deprecated TProb::TProb( TIterator begin, TIterator end, size_t sizeHint=0 )
@@ -47,8 +58,8 @@ git HEAD
   - Fixed bug by renaming TProb<T>::operator<=() to TProb<T>::operator<()
   - TProb<T>& operator/= (T x) now yields 0 when dividing by 0
   - Changed format of TProb<T> when streamed to an ostream
-  - Added get(size_t) and set(size_t,T) operators; get() is 
-    equivalent to "operator[](size_t) const" and set() should
+  - Added TProb<T>::get(size_t) and TProb<T>::set(size_t,T) operators; 
+    get() is equivalent to "operator[](size_t) const" and set() should
     be used instead of the non-const operator[], which has been deprecated
   - Added TProb<T>::resize(size_t)
 * Improved index.h/cpp:
@@ -64,16 +75,16 @@ git HEAD
   - Fixed bug in PropertySet::setAsString<T>()
 * Improved util.h/cpp:
   - Fixed a bug in rnd_seed()
-  - Removed Real max( const std::vector<Real> &v )
+  - Removed max( const std::vector<Real> &v )
 * Improved weightedgraph.h/cpp:
   - Renamed MaxSpanningTreePrims into MaxSpanningTree
   - Renamed MinSpanningTreePrims into MinSpanningTree
   - Added option to MaxSpanningTree and MinSpanningTree for
     choosing between Prim's algorithm and Kruskal's algorithm
-  - More error checking in RootedTree constructor
+  - Better error checking in RootedTree constructor
   - The vertices in DEdge and UEdge are now named "first" and "second"
     for compatibility with a std::pair<size_t,size_t>
-  - Added GraphEL( const GraphAL& G ) constructor
+  - Added GraphEL::GraphEL( const GraphAL& G ) constructor
 * Improved bipgraph.h/cpp:
   - Fixed bug in BipartiteGraph::eraseNode1()
   - Fixed bug in BipartiteGraph::eraseNode2()
@@ -101,7 +112,7 @@ git HEAD
     no longer has a default value in order to avoid confusion with the
      SmallSet::SmallSet( const T &t1, const T &t2 )
     constructor.
-  - Added SmallSet& insert( const T& t )
+  - Added SmallSet::insert( const T& t ) method
 * Removed RandomDRegularGraph(), please use createGraphRegular() instead
 * Compressed Makefile
 * Added unit tests for var, smallset, varset, graph, bipgraph,
index 9a147e6..b84a36e 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) tests/unit/clustergraph$(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) tests/unit/regiongraph$(EE)
        echo Running unit tests...
        tests/unit/var$(EE)
        tests/unit/smallset$(EE)
@@ -141,6 +141,7 @@ unittests : tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE)
        tests/unit/factor$(EE)
        tests/unit/factorgraph$(EE)
        tests/unit/clustergraph$(EE)
+       tests/unit/regiongraph$(EE)
 ifneq ($(OS),WINDOWS)
        rm factorgraph_test.fg
 else
@@ -293,7 +294,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) tests/unit/clustergraph$(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) tests/unit/regiongraph$(EE)
        -rm utils/fg2dot$(EE) utils/createfg$(EE) utils/fginfo$(EE)
        -rm -R doc
        -rm -R lib
@@ -339,6 +340,7 @@ clean :
        -del tests\unit\factor$(EE)
        -del tests\unit\factorgraph$(EE)
        -del tests\unit\clustergraph$(EE)
+       -del tests\unit\regiongraph$(EE)
        -del $(LIB)\libdai$(LE)
        -rmdir lib
 endif
index f22e462..3c6d926 100644 (file)
@@ -152,13 +152,13 @@ class JTree : public DAIAlgRG {
          *  Finally, Beliefs are constructed.
          *  If \a verify == \c true, checks whether each factor is subsumed by a clique.
          */
-        void construct( const std::vector<VarSet> &cl, bool verify=false );
+        void construct( const FactorGraph &fg, const std::vector<VarSet> &cl, bool verify=false );
 
         /// Constructs a junction tree based on the cliques \a cl (corresponding to some elimination sequence).
         /** Invokes construct() and then constructs messages.
          *  \see construct()
          */
-        void GenerateJT( const std::vector<VarSet> &cl );
+        void GenerateJT( const FactorGraph &fg, const std::vector<VarSet> &cl );
 
         /// Returns constant reference to the message from outer region \a alpha to its \a _beta 'th neighboring inner region
         const Factor & message( size_t alpha, size_t _beta ) const { return _mes[alpha][_beta]; }
index 4fcf74a..130d89f 100644 (file)
@@ -72,7 +72,7 @@ class FRegion : public Factor {
 /** A RegionGraph inherits from a FactorGraph and adds additional structure in the form of a "region graph". Our definition of region graph
  *  is inspired by [\ref HAK03], which is less general than the definition given in [\ref YFW05].
  *
- *  The extra structure described by a RegionGraph over that described by a FactorGraph is:
+ *  The extra structure described by a RegionGraph compared with that described by a FactorGraph is:
  *  - a set of outer regions (indexed by \f$\alpha\f$), where each outer region consists of
  *    - a factor defined on a subset of variables
  *    - a counting number
@@ -85,33 +85,30 @@ class FRegion : public Factor {
  *  of an outer region would be the product of all the factors that belong to that region.
  */
 class RegionGraph : public FactorGraph {
-    public:
+    protected:
         /// Stores the neighborhood structure
-        BipartiteGraph          G;
+        BipartiteGraph          _G;
 
         /// The outer regions (corresponding to nodes of type 1)
-        std::vector<FRegion>    ORs;
+        std::vector<FRegion>    _ORs;
 
         /// The inner regions (corresponding to nodes of type 2)
-        std::vector<Region>     IRs;
+        std::vector<Region>     _IRs;
 
         /// Stores for each factor index the index of the outer region it belongs to
-        std::vector<size_t>     fac2OR;
+        std::vector<size_t>     _fac2OR;
 
 
     public:
     /// \name Constructors and destructors
     //@{
         /// Default constructor
-        RegionGraph() : FactorGraph(), G(), ORs(), IRs(), fac2OR() {}
-
-        /// Partially constructs a region graph from a factor graph
-        RegionGraph( const FactorGraph &fg ) : FactorGraph(fg), G(), ORs(), IRs(), fac2OR() {}
+        RegionGraph() : FactorGraph(), _G(), _ORs(), _IRs(), _fac2OR() {}
 
         /// Constructs a region graph from a factor graph, a vector of outer regions, a vector of inner regions and a vector of edges
         /** The counting numbers for the outer regions are set to 1.
          */
-        RegionGraph( const FactorGraph &fg, const std::vector<VarSet> &ors, const std::vector<Region> &irs, const std::vector<std::pair<size_t,size_t> > &edges ) : FactorGraph(), G(), ORs(), IRs(), fac2OR() {
+        RegionGraph( const FactorGraph& fg, const std::vector<VarSet>& ors, const std::vector<Region>& irs, const std::vector<std::pair<size_t,size_t> >& edges ) : FactorGraph(), _G(), _ORs(), _IRs(), _fac2OR() {
             construct( fg, ors, irs, edges );
 
             // Check counting numbers
@@ -132,7 +129,7 @@ class RegionGraph : public FactorGraph {
          *  subset of the variables in the outer region. The counting numbers for the inner
          *  regions are calculated by calcCountingNumbers() and satisfy the Moebius formula.
          */
-        RegionGraph( const FactorGraph &fg, const std::vector<VarSet> &cl ) : FactorGraph(), G(), ORs(), IRs(), fac2OR() {
+        RegionGraph( const FactorGraph& fg, const std::vector<VarSet>& cl ) : FactorGraph(), _G(), _ORs(), _IRs(), _fac2OR() {
             constructCVM( fg, cl );
 
             // Check counting numbers
@@ -148,24 +145,43 @@ class RegionGraph : public FactorGraph {
     /// \name Queries
     //@{
         /// Returns number of outer regions
-        size_t nrORs() const { return ORs.size(); }
+        size_t nrORs() const { return _ORs.size(); }
         /// Returns number of inner regions
-        size_t nrIRs() const { return IRs.size(); }
+        size_t nrIRs() const { return _IRs.size(); }
 
         /// Returns constant reference to outer region \a alpha
-        const FRegion & OR(size_t alpha) const { return ORs[alpha]; }
+        const FRegion& OR( size_t alpha ) const {
+            DAI_DEBASSERT( alpha < nrORs() );
+            return _ORs[alpha]; 
+        }
         /// Returns reference to outer region \a alpha
-        FRegion & OR(size_t alpha) { return ORs[alpha]; }
+        FRegion& OR( size_t alpha ) {
+            DAI_DEBASSERT( alpha < nrORs() );
+            return _ORs[alpha]; 
+        }
 
         /// Returns constant reference to inner region \a beta
-        const Region & IR(size_t beta) const { return IRs[beta]; }
+        const Region& IR( size_t beta ) const {
+            DAI_DEBASSERT( beta < nrIRs() );
+            return _IRs[beta]; 
+        }
         /// Returns reference to inner region \a beta
-        Region & IR(size_t beta) { return IRs[beta]; }
+        Region& IR( size_t beta ) {
+            DAI_DEBASSERT( beta < nrIRs() );
+            return _IRs[beta];
+        }
+
+        /// Returns the index of the outer region to which the \a I 'th factor corresponds
+        size_t fac2OR( size_t I ) {
+            DAI_DEBASSERT( I < nrFactors() );
+            DAI_DEBASSERT( I < _fac2OR.size() );
+            return _fac2OR[I];
+        }
 
         /// Returns constant reference to the neighbors of outer region \a alpha
-        const Neighbors & nbOR( size_t alpha ) const { return G.nb1(alpha); }
+        const Neighbors& nbOR( size_t alpha ) const { return _G.nb1(alpha); }
         /// Returns constant reference to the neighbors of inner region \a beta
-        const Neighbors & nbIR( size_t beta ) const { return G.nb2(beta); }
+        const Neighbors& nbIR( size_t beta ) const { return _G.nb2(beta); }
 
         /// Check whether the counting numbers are valid
         /** Counting numbers are said to be (variable) valid if for each variable \f$x\f$,
@@ -179,13 +195,13 @@ class RegionGraph : public FactorGraph {
     /// \name Operations
     //@{
         /// Set the content of the \a I 'th factor and make a backup of its old content if \a backup == \c true
-        virtual void setFactor( size_t I, const Factor &newFactor, bool backup = false ) {
+        virtual void setFactor( size_t I, const FactornewFactor, bool backup = false ) {
             FactorGraph::setFactor( I, newFactor, backup );
             RecomputeOR( I );
         }
 
         /// Set the contents of all factors as specified by \a facs and make a backup of the old contents if \a backup == \c true
-        virtual void setFactors( const std::map<size_t, Factor> & facs, bool backup = false ) {
+        virtual void setFactors( const std::map<size_t, Factor>& facs, bool backup = false ) {
             FactorGraph::setFactors( facs, backup );
             VarSet ns;
             for( std::map<size_t, Factor>::const_iterator fac = facs.begin(); fac != facs.end(); fac++ )
@@ -201,7 +217,7 @@ class RegionGraph : public FactorGraph {
         /// Recompute all outer regions involving the variables in \a vs
         /** The factor contents of each outer region involving at least one of the variables in \a vs is set to the product of the factors belonging to that region.
          */
-        void RecomputeORs( const VarSet &vs );
+        void RecomputeORs( const VarSetvs );
 
         /// Recompute all outer regions involving factor \a I
         /** The factor contents of each outer region involving the \a I 'th factor is set to the product of the factors belonging to that region.
@@ -221,15 +237,15 @@ class RegionGraph : public FactorGraph {
     /// \name Input/output
     //@{
         /// Writes a RegionGraph to an output stream
-        friend std::ostream & operator << ( std::ostream & os, const RegionGraph & rg );
+        friend std::ostream& operator << ( std::ostream& os, const RegionGraph& rg );
     //@}
 
     protected:
         /// Helper function for constructors
-        void construct( const FactorGraph &fg, const std::vector<VarSet> &ors, const std::vector<Region> &irs, const std::vector<std::pair<size_t,size_t> > &edges );
+        void construct( const FactorGraph& fg, const std::vector<VarSet>& ors, const std::vector<Region>& irs, const std::vector<std::pair<size_t,size_t> >& edges );
 
         /// Helper function for constructors (CVM style)
-        void constructCVM( const FactorGraph &fg, const std::vector<VarSet> &cl );
+        void constructCVM( const FactorGraph& fg, const std::vector<VarSet>& cl );
 };
 
 
index 506120d..ce5eee7 100644 (file)
@@ -81,7 +81,7 @@ class SmallSet {
     //@{
         /// Inserts \a t into \c *this
         SmallSet& insert( const T& t ) {
-            SmallSet::iterator it = std::lower_bound( _elements.begin(), _elements.end(), t );
+            typename SmallSet<T>::iterator it = std::lower_bound( _elements.begin(), _elements.end(), t );
             if( (it == _elements.end()) || (*it != t) )
                 _elements.insert( it, t );
             return *this;
index 60470fe..2df387b 100644 (file)
@@ -200,7 +200,7 @@ class TreeEP : public JTree {
         /// Helper function for constructors
         void construct( const RootedTree &tree );
         /// Returns \c true if factor \a I is not part of the tree
-        bool offtree( size_t I ) const { return (fac2OR[I] == -1U); }
+        bool offtree( size_t I ) const { return (fac2OR(I) == -1U); }
 };
 
 
index 285fb51..0c83060 100644 (file)
@@ -421,7 +421,7 @@ Real HAK::doDoubleLoop() {
     double tic = toc();
 
     // Save original outer regions
-    vector<FRegion> org_ORs = ORs;
+    vector<FRegion> org_ORs = _ORs;
 
     // Save original inner counting numbers and set negative counting numbers to zero
     vector<Real> org_IR_cs( nrIRs(), 0.0 );
@@ -494,7 +494,7 @@ Real HAK::doDoubleLoop() {
         _maxdiff = maxDiff;
 
     // Restore original outer regions
-    ORs = org_ORs;
+    _ORs = org_ORs;
 
     // Restore original inner counting numbers
     for( size_t beta = 0; beta < nrIRs(); ++beta )
index b132302..08537b1 100644 (file)
@@ -61,18 +61,18 @@ string JTree::printProperties() const {
 }
 
 
-JTree::JTree( const FactorGraph &fg, const PropertySet &opts, bool automatic ) : DAIAlgRG(fg), _mes(), _logZ(), RTree(), Qa(), Qb(), props() {
+JTree::JTree( const FactorGraph &fg, const PropertySet &opts, bool automatic ) : DAIAlgRG(), _mes(), _logZ(), RTree(), Qa(), Qb(), props() {
     setProperties( opts );
 
-    if( !isConnected() )
+    if( !fg.isConnected() )
        DAI_THROW(FACTORGRAPH_NOT_CONNECTED);
 
     if( automatic ) {
         // Create ClusterGraph which contains factors as clusters
         vector<VarSet> cl;
         cl.reserve( fg.nrFactors() );
-        for( size_t I = 0; I < nrFactors(); I++ )
-            cl.push_back( factor(I).vars() );
+        for( size_t I = 0; I < fg.nrFactors(); I++ )
+            cl.push_back( fg.factor(I).vars() );
         ClusterGraph _cg( cl );
 
         if( props.verbose >= 3 )
@@ -106,12 +106,15 @@ JTree::JTree( const FactorGraph &fg, const PropertySet &opts, bool automatic ) :
             cerr << "VarElim result: " << ElimVec << endl;
 
         // Generate the junction tree corresponding to the elimination sequence
-        GenerateJT( ElimVec );
+        GenerateJT( fg, ElimVec );
     }
 }
 
 
-void JTree::construct( const std::vector<VarSet> &cl, bool verify ) {
+void JTree::construct( const FactorGraph &fg, const std::vector<VarSet> &cl, bool verify ) {
+    // Copy the factor graph
+    FactorGraph::operator=( fg );
+
     // Construct a weighted graph (each edge is weighted with the cardinality
     // of the intersection of the nodes, where the nodes are the elements of cl).
     WeightedGraph<int> JuncGraph;
@@ -128,20 +131,20 @@ void JTree::construct( const std::vector<VarSet> &cl, bool verify ) {
     // Construct corresponding region graph
 
     // Create outer regions
-    ORs.clear();
-    ORs.reserve( cl.size() );
+    _ORs.clear();
+    _ORs.reserve( cl.size() );
     for( size_t i = 0; i < cl.size(); i++ )
-        ORs.push_back( FRegion( Factor(cl[i], 1.0), 1.0 ) );
+        _ORs.push_back( FRegion( Factor(cl[i], 1.0), 1.0 ) );
 
     // For each factor, find an outer region that subsumes that factor.
     // Then, multiply the outer region with that factor.
-    fac2OR.clear();
-    fac2OR.resize( nrFactors(), -1U );
+    _fac2OR.clear();
+    _fac2OR.resize( nrFactors(), -1U );
     for( size_t I = 0; I < nrFactors(); I++ ) {
         size_t alpha;
         for( alpha = 0; alpha < nrORs(); alpha++ )
             if( OR(alpha).vars() >> factor(I).vars() ) {
-                fac2OR[I] = alpha;
+                _fac2OR[I] = alpha;
                 break;
             }
         if( verify )
@@ -150,19 +153,19 @@ void JTree::construct( const std::vector<VarSet> &cl, bool verify ) {
     RecomputeORs();
 
     // Create inner regions and edges
-    IRs.clear();
-    IRs.reserve( RTree.size() );
+    _IRs.clear();
+    _IRs.reserve( RTree.size() );
     vector<Edge> edges;
     edges.reserve( 2 * RTree.size() );
     for( size_t i = 0; i < RTree.size(); i++ ) {
         edges.push_back( Edge( RTree[i].first, nrIRs() ) );
         edges.push_back( Edge( RTree[i].second, nrIRs() ) );
         // inner clusters have counting number -1
-        IRs.push_back( Region( cl[RTree[i].first] & cl[RTree[i].second], -1.0 ) );
+        _IRs.push_back( Region( cl[RTree[i].first] & cl[RTree[i].second], -1.0 ) );
     }
 
     // create bipartite graph
-    G.construct( nrORs(), nrIRs(), edges.begin(), edges.end() );
+    _G.construct( nrORs(), nrIRs(), edges.begin(), edges.end() );
 
     // Check counting numbers
 #ifdef DAI_DEBUG
@@ -182,8 +185,8 @@ void JTree::construct( const std::vector<VarSet> &cl, bool verify ) {
 }
 
 
-void JTree::GenerateJT( const std::vector<VarSet> &cl ) {
-    construct( cl, true );
+void JTree::GenerateJT( const FactorGraph &fg, const std::vector<VarSet> &cl ) {
+    construct( fg, cl, true );
 
     // Create messages
     _mes.clear();
index 0755fd8..bd955ad 100644 (file)
@@ -28,23 +28,23 @@ void RegionGraph::construct( const FactorGraph &fg, const std::vector<VarSet> &o
     FactorGraph::operator=( fg );
 
     // Copy inner regions
-    IRs = irs;
+    _IRs = irs;
 
     // Construct outer regions (giving them counting number 1.0)
-    ORs.clear();
-    ORs.reserve( ors.size() );
+    _ORs.clear();
+    _ORs.reserve( ors.size() );
     foreach( const VarSet &alpha, ors )
-        ORs.push_back( FRegion(Factor(alpha, 1.0), 1.0) );
+        _ORs.push_back( FRegion(Factor(alpha, 1.0), 1.0) );
 
     // For each factor, find an outer region that subsumes that factor.
     // Then, multiply the outer region with that factor.
-    fac2OR.clear();
-    fac2OR.reserve( nrFactors() );
+    _fac2OR.clear();
+    _fac2OR.reserve( nrFactors() );
     for( size_t I = 0; I < nrFactors(); I++ ) {
         size_t alpha;
         for( alpha = 0; alpha < nrORs(); alpha++ )
             if( OR(alpha).vars() >> factor(I).vars() ) {
-                fac2OR.push_back( alpha );
+                _fac2OR.push_back( alpha );
                 break;
             }
         DAI_ASSERT( alpha != nrORs() );
@@ -52,7 +52,7 @@ void RegionGraph::construct( const FactorGraph &fg, const std::vector<VarSet> &o
     RecomputeORs();
 
     // Create bipartite graph
-    G.construct( nrORs(), nrIRs(), edges.begin(), edges.end() );
+    _G.construct( nrORs(), nrIRs(), edges.begin(), edges.end() );
 }
 
 
@@ -167,8 +167,8 @@ void RegionGraph::RecomputeORs() {
     for( size_t alpha = 0; alpha < nrORs(); alpha++ )
         OR(alpha).fill( 1.0 );
     for( size_t I = 0; I < nrFactors(); I++ )
-        if( fac2OR[I] != -1U )
-            OR( fac2OR[I] ) *= factor( I );
+        if( fac2OR(I) != -1U )
+            OR( fac2OR(I) ) *= factor( I );
 }
 
 
@@ -177,19 +177,19 @@ void RegionGraph::RecomputeORs( const VarSet &ns ) {
         if( OR(alpha).vars().intersects( ns ) )
             OR(alpha).fill( 1.0 );
     for( size_t I = 0; I < nrFactors(); I++ )
-        if( fac2OR[I] != -1U )
-            if( OR( fac2OR[I] ).vars().intersects( ns ) )
-                OR( fac2OR[I] ) *= factor( I );
+        if( fac2OR(I) != -1U )
+            if( OR( fac2OR(I) ).vars().intersects( ns ) )
+                OR( fac2OR(I) ) *= factor( I );
 }
 
 
 void RegionGraph::RecomputeOR( size_t I ) {
     DAI_ASSERT( I < nrFactors() );
-    if( fac2OR[I] != -1U ) {
-        size_t alpha = fac2OR[I];
+    if( fac2OR(I) != -1U ) {
+        size_t alpha = fac2OR(I);
         OR(alpha).fill( 1.0 );
         for( size_t J = 0; J < nrFactors(); J++ )
-            if( fac2OR[J] == alpha )
+            if( fac2OR(J) == alpha )
                 OR(alpha) *= factor( J );
     }
 }
index 9c52cb7..8fa5304 100644 (file)
@@ -123,7 +123,7 @@ void TreeEP::construct( const RootedTree &tree ) {
 
     // If no outer region can be found subsuming that factor, label the
     // factor as off-tree.
-    JTree::construct( cl, false );
+    JTree::construct( *this, cl, false );
 
     if( props.verbose >= 1 )
         cerr << "TreeEP::construct: The tree has size " << JTree::RTree.size() << endl;
index 5df6991..72fe390 100644 (file)
@@ -708,9 +708,9 @@ BOOST_AUTO_TEST_CASE( MiscOperationsTest ) {
     BOOST_CHECK_EQUAL( y[1], (x[1] + x[3] + x[5]) / x.sum() );
     y = x.marginal( v2 );
     BOOST_CHECK( y.vars() == VarSet( v2 ) );
-    BOOST_CHECK_EQUAL( y[0], (x[0] + x[1]) / x.sum() );
-    BOOST_CHECK_EQUAL( y[1], (x[2] + x[3]) / x.sum() );
-    BOOST_CHECK_EQUAL( y[2], (x[4] + x[5]) / x.sum() );
+    BOOST_CHECK_CLOSE( y[0], (x[0] + x[1]) / x.sum(), tol );
+    BOOST_CHECK_CLOSE( y[1], (x[2] + x[3]) / x.sum(), tol );
+    BOOST_CHECK_CLOSE( y[2], (x[4] + x[5]) / x.sum(), tol );
     y = x.marginal( VarSet() );
     BOOST_CHECK( y.vars() == VarSet() );
     BOOST_CHECK_EQUAL( y[0], 1.0 );
@@ -723,9 +723,9 @@ BOOST_AUTO_TEST_CASE( MiscOperationsTest ) {
     BOOST_CHECK_EQUAL( y[1], (x[1] + x[3] + x[5]) / x.sum() );
     y = x.marginal( v2, true );
     BOOST_CHECK( y.vars() == VarSet( v2 ) );
-    BOOST_CHECK_EQUAL( y[0], (x[0] + x[1]) / x.sum() );
-    BOOST_CHECK_EQUAL( y[1], (x[2] + x[3]) / x.sum() );
-    BOOST_CHECK_EQUAL( y[2], (x[4] + x[5]) / x.sum() );
+    BOOST_CHECK_CLOSE( y[0], (x[0] + x[1]) / x.sum(), tol );
+    BOOST_CHECK_CLOSE( y[1], (x[2] + x[3]) / x.sum(), tol );
+    BOOST_CHECK_CLOSE( y[2], (x[4] + x[5]) / x.sum(), tol );
     y = x.marginal( VarSet(), true );
     BOOST_CHECK( y.vars() == VarSet() );
     BOOST_CHECK_EQUAL( y[0], 1.0 );
diff --git a/tests/unit/regiongraph.cpp b/tests/unit/regiongraph.cpp
new file mode 100644 (file)
index 0000000..fdd55c0
--- /dev/null
@@ -0,0 +1,740 @@
+/*  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/regiongraph.h>
+#include <vector>
+#include <strstream>
+
+
+using namespace dai;
+
+
+const double tol = 1e-8;
+
+
+#define BOOST_TEST_MODULE RegionGraphTest
+
+
+#include <boost/test/unit_test.hpp>
+#include <boost/test/floating_point_comparison.hpp>
+
+
+BOOST_AUTO_TEST_CASE( ConstructorsTest ) {
+    RegionGraph G;
+    BOOST_CHECK_EQUAL( G.vars(), std::vector<Var>() );
+    BOOST_CHECK_EQUAL( G.factors(), std::vector<Factor>() );
+
+    std::vector<Factor> facs;
+    facs.push_back( Factor( VarSet( Var(0, 2), Var(1, 2) ) ) );
+    facs.push_back( Factor( VarSet( Var(0, 2), Var(2, 2) ) ) );
+    facs.push_back( Factor( VarSet( Var(1, 2), Var(2, 2) ) ) );
+    facs.push_back( Factor( VarSet( Var(1, 2) ) ) );
+    std::vector<Var> vars;
+    vars.push_back( Var( 0, 2 ) );
+    vars.push_back( Var( 1, 2 ) );
+    vars.push_back( Var( 2, 2 ) );
+
+    FactorGraph fg(facs);
+    RegionGraph G1( fg, fg.maximalFactorDomains() );
+    BOOST_CHECK_EQUAL( G1.vars(), vars );
+    BOOST_CHECK_EQUAL( G1.factors(), facs );
+
+    RegionGraph *G2 = G1.clone();
+    BOOST_CHECK_EQUAL( G2->vars(), vars );
+    BOOST_CHECK_EQUAL( G2->factors(), facs );
+    delete G2;
+
+    RegionGraph G3 = G1;
+    BOOST_CHECK_EQUAL( G3.vars(), vars );
+    BOOST_CHECK_EQUAL( G3.factors(), facs );
+
+    RegionGraph G4( G1 );
+    BOOST_CHECK_EQUAL( G4.vars(), vars );
+    BOOST_CHECK_EQUAL( G4.factors(), facs );
+}
+
+
+BOOST_AUTO_TEST_CASE( AccMutIterTest ) {
+    std::vector<Factor> facs;
+    facs.push_back( Factor( VarSet( Var(0, 2), Var(1, 2) ) ) );
+    facs.push_back( Factor( VarSet( Var(0, 2), Var(2, 2) ) ) );
+    facs.push_back( Factor( VarSet( Var(1, 2), Var(2, 2) ) ) );
+    facs.push_back( Factor( VarSet( Var(1, 2) ) ) );
+    std::vector<Var> vars;
+    vars.push_back( Var( 0, 2 ) );
+    vars.push_back( Var( 1, 2 ) );
+    vars.push_back( Var( 2, 2 ) );
+
+    FactorGraph fg( facs );
+    RegionGraph G( fg, fg.maximalFactorDomains() );
+    BOOST_CHECK_EQUAL( G.var(0), Var(0, 2) );
+    BOOST_CHECK_EQUAL( G.var(1), Var(1, 2) );
+    BOOST_CHECK_EQUAL( G.var(2), Var(2, 2) );
+    BOOST_CHECK_EQUAL( G.vars(), vars );
+    BOOST_CHECK_EQUAL( G.factor(0), facs[0] );
+    BOOST_CHECK_EQUAL( G.factor(1), facs[1] );
+    BOOST_CHECK_EQUAL( G.factor(2), facs[2] );
+    BOOST_CHECK_EQUAL( G.factor(3), facs[3] );
+    BOOST_CHECK_EQUAL( G.factors(), facs );
+    BOOST_CHECK_EQUAL( G.nbV(0).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbV(0,0), 0 );
+    BOOST_CHECK_EQUAL( G.nbV(0,1), 1 );
+    BOOST_CHECK_EQUAL( G.nbV(1).size(), 3 );
+    BOOST_CHECK_EQUAL( G.nbV(1,0), 0 );
+    BOOST_CHECK_EQUAL( G.nbV(1,1), 2 );
+    BOOST_CHECK_EQUAL( G.nbV(1,2), 3 );
+    BOOST_CHECK_EQUAL( G.nbV(0).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbV(2,0), 1 );
+    BOOST_CHECK_EQUAL( G.nbV(2,1), 2 );
+    BOOST_CHECK_EQUAL( G.nbF(0).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbF(0,0), 0 );
+    BOOST_CHECK_EQUAL( G.nbF(0,1), 1 );
+    BOOST_CHECK_EQUAL( G.nbF(1).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbF(1,0), 0 );
+    BOOST_CHECK_EQUAL( G.nbF(1,1), 2 );
+    BOOST_CHECK_EQUAL( G.nbF(2).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbF(2,0), 1 );
+    BOOST_CHECK_EQUAL( G.nbF(2,1), 2 );
+    BOOST_CHECK_EQUAL( G.nbF(3).size(), 1 );
+    BOOST_CHECK_EQUAL( G.nbF(3,0), 1 );
+
+    RegionGraph::const_iterator cit = G.begin();
+    RegionGraph::iterator it = G.begin();
+    for( size_t I = 0; I < G.nrFactors(); I++, cit++, it++ ) {
+        BOOST_CHECK_EQUAL( *cit, G.factor(I) );
+        BOOST_CHECK_EQUAL( *it, G.factor(I) );
+    }
+    BOOST_CHECK( cit == G.end() );
+    BOOST_CHECK( it == G.end() );
+}
+
+
+BOOST_AUTO_TEST_CASE( QueriesTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
+    VarSet v012 = v01 | v2;
+
+    RegionGraph G0;
+    BOOST_CHECK_EQUAL( G0.nrVars(), 0 );
+    BOOST_CHECK_EQUAL( G0.nrFactors(), 0 );
+    BOOST_CHECK_EQUAL( G0.nrEdges(), 0 );
+    BOOST_CHECK_THROW( G0.findVar( v0 ), Exception );
+    BOOST_CHECK_THROW( G0.findVars( v01 ), Exception );
+    BOOST_CHECK_THROW( G0.findFactor( v01 ), Exception );
+#ifdef DAI_DBEUG
+    BOOST_CHECK_THROW( G0.delta( 0 ), Exception );
+    BOOST_CHECK_THROW( G0.Delta( 0 ), Exception );
+    BOOST_CHECK_THROW( G0.delta( v0 ), Exception );
+    BOOST_CHECK_THROW( G0.Delta( v0 ), Exception );
+#endif
+    BOOST_CHECK( G0.isConnected() );
+    BOOST_CHECK( G0.isTree() );
+    BOOST_CHECK( G0.isBinary() );
+    BOOST_CHECK( G0.isPairwise() );
+    BOOST_CHECK( G0.MarkovGraph() == GraphAL() );
+    BOOST_CHECK( G0.bipGraph() == BipartiteGraph() );
+    BOOST_CHECK_EQUAL( G0.maximalFactorDomains().size(), 0 );
+
+    std::vector<Factor> facs;
+    facs.push_back( Factor( v01 ) );
+    facs.push_back( Factor( v12 ) );
+    facs.push_back( Factor( v1 ) );
+    std::vector<Var> vars;
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+    GraphAL H(3);
+    H.addEdge( 0, 1 );
+    H.addEdge( 1, 2 );
+    BipartiteGraph K(3, 3);
+    K.addEdge( 0, 0 );
+    K.addEdge( 1, 0 );
+    K.addEdge( 1, 1 );
+    K.addEdge( 2, 1 );
+    K.addEdge( 1, 2 );
+
+    FactorGraph fg( facs );
+    RegionGraph G1( fg, fg.maximalFactorDomains() );
+    BOOST_CHECK_EQUAL( G1.nrVars(), 3 );
+    BOOST_CHECK_EQUAL( G1.nrFactors(), 3 );
+    BOOST_CHECK_EQUAL( G1.nrEdges(), 5 );
+    BOOST_CHECK_EQUAL( G1.findVar( v0 ), 0 );
+    BOOST_CHECK_EQUAL( G1.findVar( v1 ), 1 );
+    BOOST_CHECK_EQUAL( G1.findVar( v2 ), 2 );
+    BOOST_CHECK_EQUAL( G1.findVars( v01 ), SmallSet<size_t>( 0, 1 ) );
+    BOOST_CHECK_EQUAL( G1.findVars( v02 ), SmallSet<size_t>( 0, 2 ) );
+    BOOST_CHECK_EQUAL( G1.findVars( v12 ), SmallSet<size_t>( 1, 2 ) );
+    BOOST_CHECK_EQUAL( G1.findFactor( v01 ), 0 );
+    BOOST_CHECK_EQUAL( G1.findFactor( v12 ), 1 );
+    BOOST_CHECK_EQUAL( G1.findFactor( v1 ), 2 );
+    BOOST_CHECK_THROW( G1.findFactor( v02 ), Exception );
+    BOOST_CHECK_EQUAL( G1.delta( 0 ), v1 );
+    BOOST_CHECK_EQUAL( G1.delta( 1 ), v02 );
+    BOOST_CHECK_EQUAL( G1.delta( 2 ), v1 );
+    BOOST_CHECK_EQUAL( G1.Delta( 0 ), v01 );
+    BOOST_CHECK_EQUAL( G1.Delta( 1 ), v012 );
+    BOOST_CHECK_EQUAL( G1.Delta( 2 ), v12 );
+    BOOST_CHECK_EQUAL( G1.delta( v0 ), v1 );
+    BOOST_CHECK_EQUAL( G1.delta( v1 ), v02 );
+    BOOST_CHECK_EQUAL( G1.delta( v2 ), v1 );
+    BOOST_CHECK_EQUAL( G1.delta( v01 ), v2 );
+    BOOST_CHECK_EQUAL( G1.delta( v02 ), v1 );
+    BOOST_CHECK_EQUAL( G1.delta( v12 ), v0 );
+    BOOST_CHECK_EQUAL( G1.delta( v012 ), VarSet() );
+    BOOST_CHECK_EQUAL( G1.Delta( v0 ), v01 );
+    BOOST_CHECK_EQUAL( G1.Delta( v1 ), v012 );
+    BOOST_CHECK_EQUAL( G1.Delta( v2 ), v12 );
+    BOOST_CHECK_EQUAL( G1.Delta( v01 ), v012 );
+    BOOST_CHECK_EQUAL( G1.Delta( v02 ), v012 );
+    BOOST_CHECK_EQUAL( G1.Delta( v12 ), v012 );
+    BOOST_CHECK_EQUAL( G1.Delta( v012 ), v012 );
+    BOOST_CHECK( G1.isConnected() );
+    BOOST_CHECK( G1.isTree() );
+    BOOST_CHECK( G1.isBinary() );
+    BOOST_CHECK( G1.isPairwise() );
+    BOOST_CHECK( G1.MarkovGraph() == H );
+    BOOST_CHECK( G1.bipGraph() == K );
+    BOOST_CHECK_EQUAL( G1.maximalFactorDomains().size(), 2 );
+    BOOST_CHECK_EQUAL( G1.maximalFactorDomains()[0], v01 );
+    BOOST_CHECK_EQUAL( G1.maximalFactorDomains()[1], v12 );
+
+    facs.push_back( Factor( v02 ) );
+    H.addEdge( 0, 2 );
+    K.addNode2();
+    K.addEdge( 0, 3 );
+    K.addEdge( 2, 3 );
+    fg = FactorGraph( facs );
+    RegionGraph G2( fg, fg.maximalFactorDomains() );
+    BOOST_CHECK_EQUAL( G2.nrVars(), 3 );
+    BOOST_CHECK_EQUAL( G2.nrFactors(), 4 );
+    BOOST_CHECK_EQUAL( G2.nrEdges(), 7 );
+    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.findVars( v01 ), SmallSet<size_t>( 0, 1 ) );
+    BOOST_CHECK_EQUAL( G2.findVars( v02 ), SmallSet<size_t>( 0, 2 ) );
+    BOOST_CHECK_EQUAL( G2.findVars( v12 ), SmallSet<size_t>( 1, 2 ) );
+    BOOST_CHECK_EQUAL( G2.findFactor( v01 ), 0 );
+    BOOST_CHECK_EQUAL( G2.findFactor( v12 ), 1 );
+    BOOST_CHECK_EQUAL( G2.findFactor( v1 ), 2 );
+    BOOST_CHECK_EQUAL( G2.findFactor( v02 ), 3 );
+    BOOST_CHECK_EQUAL( G2.delta( 0 ), v12 );
+    BOOST_CHECK_EQUAL( G2.delta( 1 ), v02 );
+    BOOST_CHECK_EQUAL( G2.delta( 2 ), v01 );
+    BOOST_CHECK_EQUAL( G2.Delta( 0 ), v012 );
+    BOOST_CHECK_EQUAL( G2.Delta( 1 ), v012 );
+    BOOST_CHECK_EQUAL( G2.Delta( 2 ), v012 );
+    BOOST_CHECK( G2.isConnected() );
+    BOOST_CHECK( !G2.isTree() );
+    BOOST_CHECK( G2.isBinary() );
+    BOOST_CHECK( G2.isPairwise() );
+    BOOST_CHECK( G2.MarkovGraph() == H );
+    BOOST_CHECK( G2.bipGraph() == K );
+    BOOST_CHECK_EQUAL( G2.maximalFactorDomains().size(), 3 );
+    BOOST_CHECK_EQUAL( G2.maximalFactorDomains()[0], v01 );
+    BOOST_CHECK_EQUAL( G2.maximalFactorDomains()[1], v12 );
+    BOOST_CHECK_EQUAL( G2.maximalFactorDomains()[2], v02 );
+
+    Var v3( 3, 3 );
+    VarSet v03( v0, v3 );
+    VarSet v13( v1, v3 );
+    VarSet v23( v2, v3 );
+    VarSet v013 = v01 | v3;
+    VarSet v023 = v02 | v3;
+    VarSet v123 = v12 | v3;
+    VarSet v0123 = v012 | v3;
+    vars.push_back( v3 );
+    facs.push_back( Factor( v3 ) );
+    H.addNode();
+    K.addNode1();
+    K.addNode2();
+    K.addEdge( 3, 4 );
+    fg = FactorGraph( facs );
+    RegionGraph G3( fg, fg.maximalFactorDomains() );
+    BOOST_CHECK_EQUAL( G3.nrVars(), 4 );
+    BOOST_CHECK_EQUAL( G3.nrFactors(), 5 );
+    BOOST_CHECK_EQUAL( G3.nrEdges(), 8 );
+    BOOST_CHECK_EQUAL( G3.findVar( v0 ), 0 );
+    BOOST_CHECK_EQUAL( G3.findVar( v1 ), 1 );
+    BOOST_CHECK_EQUAL( G3.findVar( v2 ), 2 );
+    BOOST_CHECK_EQUAL( G3.findVar( v3 ), 3 );
+    BOOST_CHECK_EQUAL( G3.findVars( v01 ), SmallSet<size_t>( 0, 1 ) );
+    BOOST_CHECK_EQUAL( G3.findVars( v02 ), SmallSet<size_t>( 0, 2 ) );
+    BOOST_CHECK_EQUAL( G3.findVars( v12 ), SmallSet<size_t>( 1, 2 ) );
+    BOOST_CHECK_EQUAL( G3.findFactor( v01 ), 0 );
+    BOOST_CHECK_EQUAL( G3.findFactor( v12 ), 1 );
+    BOOST_CHECK_EQUAL( G3.findFactor( v1 ), 2 );
+    BOOST_CHECK_EQUAL( G3.findFactor( v02 ), 3 );
+    BOOST_CHECK_EQUAL( G3.findFactor( v3 ), 4 );
+    BOOST_CHECK_THROW( G3.findFactor( v23 ), Exception );
+    BOOST_CHECK_EQUAL( G3.delta( 0 ), v12 );
+    BOOST_CHECK_EQUAL( G3.delta( 1 ), v02 );
+    BOOST_CHECK_EQUAL( G3.delta( 2 ), v01 );
+    BOOST_CHECK_EQUAL( G3.delta( 3 ), VarSet() );
+    BOOST_CHECK_EQUAL( G3.Delta( 0 ), v012 );
+    BOOST_CHECK_EQUAL( G3.Delta( 1 ), v012 );
+    BOOST_CHECK_EQUAL( G3.Delta( 2 ), v012 );
+    BOOST_CHECK_EQUAL( G3.Delta( 3 ), v3 );
+    BOOST_CHECK( !G3.isConnected() );
+    BOOST_CHECK( !G3.isTree() );
+    BOOST_CHECK( !G3.isBinary() );
+    BOOST_CHECK( G3.isPairwise() );
+    BOOST_CHECK( G3.MarkovGraph() == H );
+    BOOST_CHECK( G3.bipGraph() == K );
+    BOOST_CHECK_EQUAL( G3.maximalFactorDomains().size(), 4 );
+    BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[0], v01 );
+    BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[1], v12 );
+    BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[2], v02 );
+    BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[3], v3 );
+
+    facs.push_back( Factor( v123 ) );
+    H.addEdge( 1, 3 );
+    H.addEdge( 2, 3 );
+    K.addNode2();
+    K.addEdge( 1, 5 );
+    K.addEdge( 2, 5 );
+    K.addEdge( 3, 5 );
+    fg = FactorGraph( facs );
+    RegionGraph G4( fg, fg.maximalFactorDomains() );
+    BOOST_CHECK_EQUAL( G4.nrVars(), 4 );
+    BOOST_CHECK_EQUAL( G4.nrFactors(), 6 );
+    BOOST_CHECK_EQUAL( G4.nrEdges(), 11 );
+    BOOST_CHECK_EQUAL( G4.findVar( v0 ), 0 );
+    BOOST_CHECK_EQUAL( G4.findVar( v1 ), 1 );
+    BOOST_CHECK_EQUAL( G4.findVar( v2 ), 2 );
+    BOOST_CHECK_EQUAL( G4.findVar( v3 ), 3 );
+    BOOST_CHECK_EQUAL( G4.findVars( v01 ), SmallSet<size_t>( 0, 1 ) );
+    BOOST_CHECK_EQUAL( G4.findVars( v02 ), SmallSet<size_t>( 0, 2 ) );
+    BOOST_CHECK_EQUAL( G4.findVars( v12 ), SmallSet<size_t>( 1, 2 ) );
+    BOOST_CHECK_EQUAL( G4.findFactor( v01 ), 0 );
+    BOOST_CHECK_EQUAL( G4.findFactor( v12 ), 1 );
+    BOOST_CHECK_EQUAL( G4.findFactor( v1 ), 2 );
+    BOOST_CHECK_EQUAL( G4.findFactor( v02 ), 3 );
+    BOOST_CHECK_EQUAL( G4.findFactor( v3 ), 4 );
+    BOOST_CHECK_EQUAL( G4.findFactor( v123 ), 5 );
+    BOOST_CHECK_THROW( G4.findFactor( v23 ), Exception );
+    BOOST_CHECK_EQUAL( G4.delta( 0 ), v12 );
+    BOOST_CHECK_EQUAL( G4.delta( 1 ), v023 );
+    BOOST_CHECK_EQUAL( G4.delta( 2 ), v013 );
+    BOOST_CHECK_EQUAL( G4.delta( 3 ), v12 );
+    BOOST_CHECK_EQUAL( G4.Delta( 0 ), v012 );
+    BOOST_CHECK_EQUAL( G4.Delta( 1 ), v0123 );
+    BOOST_CHECK_EQUAL( G4.Delta( 2 ), v0123 );
+    BOOST_CHECK_EQUAL( G4.Delta( 3 ), v123 );
+    BOOST_CHECK( G4.isConnected() );
+    BOOST_CHECK( !G4.isTree() );
+    BOOST_CHECK( !G4.isBinary() );
+    BOOST_CHECK( !G4.isPairwise() );
+    BOOST_CHECK( G4.MarkovGraph() == H );
+    BOOST_CHECK( G4.bipGraph() == K );
+    BOOST_CHECK_EQUAL( G4.maximalFactorDomains().size(), 3 );
+    BOOST_CHECK_EQUAL( G4.maximalFactorDomains()[0], v01 );
+    BOOST_CHECK_EQUAL( G4.maximalFactorDomains()[1], v02 );
+    BOOST_CHECK_EQUAL( G4.maximalFactorDomains()[2], v123 );
+}
+
+
+BOOST_AUTO_TEST_CASE( BackupRestoreTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
+    VarSet v012 = v01 | v2;
+
+    std::vector<Factor> facs;
+    facs.push_back( Factor( v01 ) );
+    facs.push_back( Factor( v12 ) );
+    facs.push_back( Factor( v1 ) );
+    std::vector<Var> vars;
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+
+    FactorGraph fg( facs );
+    RegionGraph G( fg, fg.maximalFactorDomains() );
+    RegionGraph Gorg( G );
+
+    BOOST_CHECK_THROW( G.setFactor( 0, Factor( v0 ), false ), Exception );
+    G.setFactor( 0, Factor( v01, 2.0 ), false );
+    BOOST_CHECK_THROW( G.restoreFactor( 0 ), Exception );
+    G.setFactor( 0, Factor( v01, 3.0 ), true );
+    G.restoreFactor( 0 );
+    BOOST_CHECK_EQUAL( G.factor( 0 )[0], 2.0 );
+    G.setFactor( 0, Gorg.factor( 0 ), false );
+    G.backupFactor( 0 );
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    G.setFactor( 0, Factor( v01, 2.0 ), false );
+    BOOST_CHECK_EQUAL( G.factor( 0 )[0], 2.0 );
+    G.restoreFactor( 0 );
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+
+    std::map<size_t, Factor> fs;
+    fs[0] = Factor( v01, 3.0 );
+    fs[2] = Factor( v1, 2.0 );
+    G.setFactors( fs, false );
+    BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
+    BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    G.restoreFactors();
+    BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
+    BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    G = Gorg;
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+    G.setFactors( fs, true );
+    BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
+    BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    G.restoreFactors();
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+    std::set<size_t> fsind;
+    fsind.insert( 0 );
+    fsind.insert( 2 );
+    G.backupFactors( fsind );
+    G.setFactors( fs, false );
+    BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
+    BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    G.restoreFactors();
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+
+    G.backupFactors( v2 );
+    G.setFactor( 1, Factor(v12, 5.0) );
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Factor(v12, 5.0) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+    G.restoreFactors( v2 );
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+
+    G.backupFactors( v1 );
+    fs[1] = Factor( v12, 5.0 );
+    G.setFactors( fs, false );
+    BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
+    BOOST_CHECK_EQUAL( G.factor(1), fs[1] );
+    BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    G.restoreFactors();
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+    G.setFactors( fs, true );
+    BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
+    BOOST_CHECK_EQUAL( G.factor(1), fs[1] );
+    BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    G.restoreFactors( v1 );
+    BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.factor(1), Gorg.factor(1) );
+    BOOST_CHECK_EQUAL( G.factor(2), Gorg.factor(2) );
+}
+
+
+BOOST_AUTO_TEST_CASE( TransformationsTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
+    VarSet v012 = v01 | v2;
+
+    std::vector<Factor> facs;
+    facs.push_back( Factor( v01 ).randomize() );
+    facs.push_back( Factor( v12 ).randomize() );
+    facs.push_back( Factor( v1 ).randomize() );
+    std::vector<Var> vars;
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+
+    FactorGraph fg( facs );
+    RegionGraph G( fg, fg.maximalFactorDomains() );
+
+    FactorGraph Gsmall = G.maximalFactors();
+    BOOST_CHECK_EQUAL( Gsmall.nrVars(), 3 );
+    BOOST_CHECK_EQUAL( Gsmall.nrFactors(), 2 );
+    BOOST_CHECK_EQUAL( Gsmall.factor( 0 ), G.factor( 0 ) * G.factor( 2 ) );
+    BOOST_CHECK_EQUAL( Gsmall.factor( 1 ), G.factor( 1 ) );
+
+    size_t i = 0;
+    for( size_t x = 0; x < 2; x++ ) {
+        FactorGraph Gcl = G.clamped( i, x );
+        BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+        BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+        BOOST_CHECK_EQUAL( Gcl.factor(0), createFactorDelta(vars[i], x) );
+        BOOST_CHECK_EQUAL( Gcl.factor(1), G.factor(0).slice(vars[i], x) * G.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.factor(2), G.factor(1) );
+    }
+    i = 1;
+    for( size_t x = 0; x < 2; x++ ) {
+        FactorGraph Gcl = G.clamped( i, x );
+        BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+        BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+        BOOST_CHECK_EQUAL( Gcl.factor(0), createFactorDelta(vars[i], x) * G.factor(2).slice(vars[i],x) );
+        BOOST_CHECK_EQUAL( Gcl.factor(1), G.factor(0).slice(vars[i], x) );
+        BOOST_CHECK_EQUAL( Gcl.factor(2), G.factor(1).slice(vars[i], x) );
+    }
+    i = 2;
+    for( size_t x = 0; x < 2; x++ ) {
+        FactorGraph Gcl = G.clamped( i, x );
+        BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+        BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+        BOOST_CHECK_EQUAL( Gcl.factor(0), createFactorDelta(vars[i], x) );
+        BOOST_CHECK_EQUAL( Gcl.factor(1), G.factor(0) );
+        BOOST_CHECK_EQUAL( Gcl.factor(2), G.factor(1).slice(vars[i], x) * G.factor(2) );
+    }
+}
+
+
+BOOST_AUTO_TEST_CASE( OperationsTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
+    VarSet v012 = v01 | v2;
+
+    std::vector<Factor> facs;
+    facs.push_back( Factor( v01 ).randomize() );
+    facs.push_back( Factor( v12 ).randomize() );
+    facs.push_back( Factor( v1 ).randomize() );
+    std::vector<Var> vars;
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+
+    FactorGraph fg( facs );
+    RegionGraph G( fg, fg.maximalFactorDomains() );
+
+    // clamp
+    RegionGraph Gcl = G;
+    for( size_t i = 0; i < 3; i++ )
+        for( size_t x = 0; x < 2; x++ ) {
+            Gcl.clamp( i, x, true );
+            Factor delta = createFactorDelta( vars[i], x );
+            BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+            BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+            for( size_t j = 0; j < 3; j++ )
+                if( G.factor(j).vars().contains( vars[i] ) )
+                    BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) * delta );
+                else
+                    BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+
+            Gcl.restoreFactors();
+            for( size_t j = 0; j < 3; j++ )
+                BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+        }
+
+    // clampVar
+    for( size_t i = 0; i < 3; i++ )
+        for( size_t x = 0; x < 2; x++ ) {
+            Gcl.clampVar( i, std::vector<size_t>(1, x), true );
+            Factor delta = createFactorDelta( vars[i], x );
+            BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+            BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+            for( size_t j = 0; j < 3; j++ )
+                if( G.factor(j).vars().contains( vars[i] ) )
+                    BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) * delta );
+                else
+                    BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+
+            Gcl.restoreFactors();
+            for( size_t j = 0; j < 3; j++ )
+                BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+        }
+    for( size_t i = 0; i < 3; i++ )
+        for( size_t x = 0; x < 2; x++ ) {
+            Gcl.clampVar( i, std::vector<size_t>(), true );
+            BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+            BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+            for( size_t j = 0; j < 3; j++ )
+                if( G.factor(j).vars().contains( vars[i] ) )
+                    BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) * 0.0 );
+                else
+                    BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+
+            Gcl.restoreFactors();
+            for( size_t j = 0; j < 3; j++ )
+                BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+        }
+    std::vector<size_t> both;
+    both.push_back( 0 );
+    both.push_back( 1 );
+    for( size_t i = 0; i < 3; i++ )
+        for( size_t x = 0; x < 2; x++ ) {
+            Gcl.clampVar( i, both, true );
+            BOOST_CHECK_EQUAL( Gcl.nrVars(), 3 );
+            BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
+            for( size_t j = 0; j < 3; j++ )
+                BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            Gcl.restoreFactors();
+        }
+
+    // clampFactor
+    for( size_t x = 0; x < 4; x++ ) {
+        Gcl.clampFactor( 0, std::vector<size_t>(1,x), true );
+        Factor mask( v01, 0.0 );
+        mask.set( x, 1.0 );
+        BOOST_CHECK_EQUAL( Gcl.factor(0), G.factor(0) * mask );
+        BOOST_CHECK_EQUAL( Gcl.factor(1), G.factor(1) );
+        BOOST_CHECK_EQUAL( Gcl.factor(2), G.factor(2) );
+        Gcl.restoreFactor( 0 );
+    }
+    for( size_t x = 0; x < 4; x++ ) {
+        Gcl.clampFactor( 1, std::vector<size_t>(1,x), true );
+        Factor mask( v12, 0.0 );
+        mask.set( x, 1.0 );
+        BOOST_CHECK_EQUAL( Gcl.factor(0), G.factor(0) );
+        BOOST_CHECK_EQUAL( Gcl.factor(1), G.factor(1) * mask );
+        BOOST_CHECK_EQUAL( Gcl.factor(2), G.factor(2) );
+        Gcl.restoreFactor( 1 );
+    }
+    for( size_t x = 0; x < 2; x++ ) {
+        Gcl.clampFactor( 2, std::vector<size_t>(1,x), true );
+        Factor mask( v1, 0.0 );
+        mask.set( x, 1.0 );
+        BOOST_CHECK_EQUAL( Gcl.factor(0), G.factor(0) );
+        BOOST_CHECK_EQUAL( Gcl.factor(1), G.factor(1) );
+        BOOST_CHECK_EQUAL( Gcl.factor(2), G.factor(2) * mask );
+        Gcl.restoreFactors();
+    }
+
+    // makeCavity
+    RegionGraph Gcav( G );
+    Gcav.makeCavity( 0, true );
+    BOOST_CHECK_EQUAL( Gcav.factor(0), Factor( v01, 1.0 ) );
+    BOOST_CHECK_EQUAL( Gcav.factor(1), G.factor(1) );
+    BOOST_CHECK_EQUAL( Gcav.factor(2), G.factor(2) );
+    Gcav.restoreFactors();
+    Gcav.makeCavity( 1, true );
+    BOOST_CHECK_EQUAL( Gcav.factor(0), Factor( v01, 1.0 ) );
+    BOOST_CHECK_EQUAL( Gcav.factor(1), Factor( v12, 1.0 ) );
+    BOOST_CHECK_EQUAL( Gcav.factor(2), Factor( v1, 1.0 ) );
+    Gcav.restoreFactors();
+    Gcav.makeCavity( 2, true );
+    BOOST_CHECK_EQUAL( Gcav.factor(0), G.factor(0) );
+    BOOST_CHECK_EQUAL( Gcav.factor(1), Factor( v12, 1.0 ) );
+    BOOST_CHECK_EQUAL( Gcav.factor(2), G.factor(2) );
+    Gcav.restoreFactors();
+}
+
+
+BOOST_AUTO_TEST_CASE( IOTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
+    VarSet v012 = v01 | v2;
+
+    std::vector<Factor> facs;
+    facs.push_back( Factor( v01 ).randomize() );
+    facs.push_back( Factor( v12 ).randomize() );
+    facs.push_back( Factor( v1 ).randomize() );
+    std::vector<Var> vars;
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+
+    FactorGraph fg( facs );
+    RegionGraph G( fg, fg.maximalFactorDomains() );
+
+    G.WriteToFile( "regiongraph_test.fg" );
+    RegionGraph G2;
+    G2.ReadFromFile( "regiongraph_test.fg" );
+
+    BOOST_CHECK( G.vars() == G2.vars() );
+    BOOST_CHECK( G.bipGraph() == G2.bipGraph() );
+    BOOST_CHECK_EQUAL( G.nrFactors(), G2.nrFactors() );
+    for( size_t I = 0; I < G.nrFactors(); I++ ) {
+        BOOST_CHECK( G.factor(I).vars() == G2.factor(I).vars() );
+        for( size_t s = 0; s < G.factor(I).nrStates(); s++ )
+            BOOST_CHECK_CLOSE( G.factor(I)[s], G2.factor(I)[s], tol );
+    }
+
+    std::stringstream ss;
+    std::string s;
+    G.printDot( ss );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "graph G {" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "node[shape=circle,width=0.4,fixedsize=true];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "node[shape=box,width=0.3,height=0.3,fixedsize=true];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tf0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tf1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tf2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv0 -- f0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv1 -- f0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv1 -- f1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv1 -- f2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tv2 -- f1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "}" );
+
+/*    G.factor(0).fill(1.0);
+    G.factor(1).fill(2.0);
+    G.factor(2).fill(3.0);
+    ss << G;
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "3" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "0 1 " );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2 2 " );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "4" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "0          1" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "1          1" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2          1" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "3          1" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "1 2 " );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2 2 " );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "4" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "0          2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "1          2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2          2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "3          2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "1" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "1 " );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2 " );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "2" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "0          3" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "1          3" );
+*/
+/*    ss << G;
+    RegionGraph G3;
+    ss >> G3;
+    BOOST_CHECK( G.vars() == G3.vars() );
+    BOOST_CHECK( G.bipGraph() == G3.bipGraph() );
+    BOOST_CHECK_EQUAL( G.nrFactors(), G3.nrFactors() );
+    for( size_t I = 0; I < G.nrFactors(); I++ ) {
+        BOOST_CHECK( G.factor(I).vars() == G3.factor(I).vars() );
+        for( size_t s = 0; s < G.factor(I).nrStates(); s++ )
+            BOOST_CHECK_CLOSE( G.factor(I)[s], G3.factor(I)[s], tol );
+    }*/
+}