Miscellaneous improvements to regiongraph, factorgraph and bipgraph, and finished...
authorJoris Mooij <joris.mooij@tuebingen.mpg.de>
Thu, 8 Apr 2010 10:38:39 +0000 (12:38 +0200)
committerJoris Mooij <joris.mooij@tuebingen.mpg.de>
Thu, 8 Apr 2010 10:38:39 +0000 (12:38 +0200)
* Improved regiongraph.h/cpp:
  - Renamed RegionGraph::calcCountingNumbers() into
    RegionGraph::calcCVMCountingNumbers() and made it protected, because
    exposing it to the world serves no purpose.
  - Added RegionGraph::DAG() accessor which returns a reference to the
    region graph DAG structure (currently implemented as a BipartiteGraph)
  - Made RegionGraph::RecomputeORs(), RegionGraph::RecomputeORs( const VarSet& )
    and RegionGraph::RecomputeOR( size_t ) protected and renamed them by changing
    the "Recompute" into "recompute", because exposing them to the world serves no
    purpose.
  - RegionGraph::WriteToFile( const char* ), RegionGraph::ReadFromFile( const char * )
    and RegionGraph::printDot( std::ostream& ) incorrectly called their respective
    FactorGraph ancestor methods; this has been corrected by letting them throw a
    NOT_IMPLEMENTED exception.
  - Changed the way a RegionGraph is streamed to an std::ostream.
* Improved factorgraph.h/cpp:
  - Deprecated the iterator interface:
    o FactorGraph::iterator typedef
    o FactorGraph::const_iterator typedef
    o FactorGraph::begin() members
    o FactorGraph::end() members
  - Deprecated FactorGraph::factor(size_t) which offered write access to a factor
    because using this functionality in the descendant RegionGraph is dangerous,
    as the outer regions are not automatically recomputed.
* Improved bipgraph.h/cpp:
  - Added operator<<( std::ostream&, const BipartiteGraph& )

21 files changed:
ChangeLog
Makefile
README
include/dai/bipgraph.h
include/dai/doc.h
include/dai/factorgraph.h
include/dai/graph.h
include/dai/hak.h
include/dai/regiongraph.h
include/dai/var.h
include/dai/weightedgraph.h
src/bbp.cpp
src/bipgraph.cpp
src/factorgraph.cpp
src/jtree.cpp
src/regiongraph.cpp
tests/unit/bipgraph.cpp
tests/unit/clustergraph.cpp
tests/unit/factorgraph.cpp
tests/unit/graph.cpp
tests/unit/regiongraph.cpp

index f1c3c66..e0345fa 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -16,6 +16,20 @@ git HEAD
     RegionGraph::IRs and RegionGraph::fac2OR protected.
   - Removed partial constructor RegionGraph::RegionGraph( const FactorGraph& fg )
   - Added some error checking code
+  - Renamed RegionGraph::calcCountingNumbers() into
+    RegionGraph::calcCVMCountingNumbers() and made it protected, because
+    exposing it to the world serves no purpose.
+  - Added RegionGraph::DAG() accessor which returns a reference to the 
+    region graph DAG structure (currently implemented as a BipartiteGraph)
+  - Made RegionGraph::RecomputeORs(), RegionGraph::RecomputeORs( const VarSet& )
+    and RegionGraph::RecomputeOR( size_t ) protected and renamed them by changing
+    the "Recompute" into "recompute", because exposing them to the world serves no 
+    purpose.
+  - RegionGraph::WriteToFile( const char* ), RegionGraph::ReadFromFile( const char * )
+    and RegionGraph::printDot( std::ostream& ) incorrectly called their respective
+    FactorGraph ancestor methods; this has been corrected by letting them throw a
+    NOT_IMPLEMENTED exception.
+  - Changed the way a RegionGraph is streamed to an std::ostream.
 * Improved clustergraph.h/cpp:
   - Made (previously public) members ClusterGraph::G, ClusterGraph::vars and 
     ClusterGraph::clusters private and added ClusterGraph::bipGraph(), 
@@ -35,6 +49,14 @@ git HEAD
     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.
+  - Deprecated the iterator interface:
+    o FactorGraph::iterator typedef
+    o FactorGraph::const_iterator typedef
+    o FactorGraph::begin() members
+    o FactorGraph::end() members
+  - Deprecated FactorGraph::factor(size_t) which offered write access to a factor
+    because using this functionality in the descendant RegionGraph is dangerous,
+    as the outer regions are not automatically recomputed.
 * Improved factor.h/cpp:
   - Fixed bug in min( const TFactor<T> &, const TFactor<T> & )
   - Fixed bug in max( const TFactor<T> &, const TFactor<T> & )
@@ -99,6 +121,7 @@ git HEAD
     return a SmallSet<size_t> instead of a vector<size_t>
   - Added BipartiteGraph::operator==( const BipartiteGraph& )
   - Added BipartiteGraph( size_t nr1, size_t nr2 ) constructor
+  - Added operator<<( std::ostream&, const BipartiteGraph& )
 * Improved graph.h/cpp:
   - Fixed bug in GraphAL::nrEdges()
   - Fixed bug in GraphAL::addEdge()
@@ -120,7 +143,7 @@ git HEAD
 * Compressed Makefile
 * Added unit tests for var, smallset, varset, graph, bipgraph,
   weightedgraph, enum, util, exceptions, properties, index, prob, 
-  factor, factorgraph, clustergraph
+  factor, factorgraph, clustergraph, regiongraph
 * Added unit testing framework
 * Added initialization of TRWBP weights by sampling spanning trees
 * Cleaned up MR code:
index b84a36e..e83da4e 100644 (file)
--- a/Makefile
+++ b/Makefile
@@ -142,11 +142,6 @@ unittests : tests/unit/var$(EE) tests/unit/smallset$(EE) tests/unit/varset$(EE)
        tests/unit/factorgraph$(EE)
        tests/unit/clustergraph$(EE)
        tests/unit/regiongraph$(EE)
-ifneq ($(OS),WINDOWS)
-       rm factorgraph_test.fg
-else
-       del factorgraph_test.fg
-endif
 
 tests : tests/testdai$(EE) tests/testem/testem$(EE) tests/testbbp$(EE) $(unittests)
 
@@ -287,19 +282,19 @@ TAGS :
 # CLEAN
 ########
 
-ifneq ($(OS),WINDOWS)
 .PHONY : clean
+ifneq ($(OS),WINDOWS)
 clean :
        -rm $(OBJECTS)
        -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) tests/unit/regiongraph$(EE)
+       -rm factorgraph_test.fg
        -rm utils/fg2dot$(EE) utils/createfg$(EE) utils/fginfo$(EE)
        -rm -R doc
        -rm -R lib
 else
-.PHONY : clean
 clean :
        -del *.obj
        -del *.ilk
@@ -341,6 +336,7 @@ clean :
        -del tests\unit\factorgraph$(EE)
        -del tests\unit\clustergraph$(EE)
        -del tests\unit\regiongraph$(EE)
+       -del factorgraph_test.fg
        -del $(LIB)\libdai$(LE)
        -rmdir lib
 endif
diff --git a/README b/README
index 9dd8b2e..efe975b 100644 (file)
--- a/README
+++ b/README
@@ -3,7 +3,7 @@ libDAI  -  A free/open source C++ library for Discrete Approximate Inference
 -------------------------------------------------------------------------------
 
 Version:  git HEAD
-Date:     February 11, 2010 - or later
+Date:     April 1, 2010 - or later
 See also: http://www.libdai.org
 
 -------------------------------------------------------------------------------
index 47aa8e3..0cc2796 100644 (file)
@@ -155,10 +155,11 @@ class BipartiteGraph {
          *  \param nrNodes2 The number of nodes of type 2.
          *  \param begin Points to the first edge.
          *  \param end Points just beyond the last edge.
+         *  \param check Whether to only add an edge if it does not exist already.
          */
         template<typename EdgeInputIterator>
-        BipartiteGraph( size_t nrNodes1, size_t nrNodes2, EdgeInputIterator begin, EdgeInputIterator end ) : _nb1(), _nb2() {
-            construct( nrNodes1, nrNodes2, begin, end );
+        BipartiteGraph( size_t nrNodes1, size_t nrNodes2, EdgeInputIterator begin, EdgeInputIterator end, bool check=true ) : _nb1(), _nb2() {
+            construct( nrNodes1, nrNodes2, begin, end, check );
         }
     //@}
 
@@ -221,9 +222,10 @@ class BipartiteGraph {
          *  \param nrNodes2 The number of nodes of type 2.
          *  \param begin Points to the first edge.
          *  \param end Points just beyond the last edge.
+         *  \param check Whether to only add an edge if it does not exist already.
          */
         template<typename EdgeInputIterator>
-        void construct( size_t nrNodes1, size_t nrNodes2, EdgeInputIterator begin, EdgeInputIterator end );
+        void construct( size_t nrNodes1, size_t nrNodes2, EdgeInputIterator begin, EdgeInputIterator end, bool check=true );
 
         /// Adds a node of type 1 without neighbors and returns the index of the added node.
         size_t addNode1() { _nb1.push_back( Neighbors() ); return _nb1.size() - 1; }
@@ -400,24 +402,25 @@ class BipartiteGraph {
     //@{
         /// Writes this BipartiteGraph to an output stream in GraphViz .dot syntax
         void printDot( std::ostream& os ) const;
+
+        /// Writes this BipartiteGraph to an output stream
+        friend std::ostream& operator<<( std::ostream& os, const BipartiteGraph& g ) {
+            g.printDot( os );
+            return os;
+        }
     //@}
 };
 
 
 template<typename EdgeInputIterator>
-void BipartiteGraph::construct( size_t nrNodes1, size_t nrNodes2, EdgeInputIterator begin, EdgeInputIterator end ) {
+void BipartiteGraph::construct( size_t nrNodes1, size_t nrNodes2, EdgeInputIterator begin, EdgeInputIterator end, bool check ) {
     _nb1.clear();
     _nb1.resize( nrNodes1 );
     _nb2.clear();
     _nb2.resize( nrNodes2 );
 
-    for( EdgeInputIterator e = begin; e != end; e++ ) {
-#ifdef DAI_DEBUG
-        addEdge( e->first, e->second, true );
-#else
-        addEdge( e->first, e->second, false );
-#endif
-    }
+    for( EdgeInputIterator e = begin; e != end; e++ )
+        addEdge( e->first, e->second, check );
 }
 
 
index d6f5ed2..d25bad3 100644 (file)
@@ -44,7 +44,7 @@
 /** \mainpage Reference manual for libDAI - A free/open source C++ library for Discrete Approximate Inference methods
  *  \author Joris Mooij, Frederik Eaton
  *  \version git HEAD
- *  \date February 11, 2010 - or later
+ *  \date April 1, 2010 - or later
  *
  *  <hr size="1">
  *  \section about About libDAI
index 0ee02b3..e7bf19e 100644 (file)
@@ -75,9 +75,13 @@ class FactorGraph {
         typedef BipartiteGraph::Edge      Edge;
 
         /// Iterator over factors
+        /** \deprecated The FactorGraph iterator interface will be removed in future versions of libDAI
+         */
         typedef std::vector<Factor>::iterator iterator;
 
         /// Constant iterator over factors
+        /** \deprecated The FactorGraph iterator interface will be removed in future versions of libDAI
+         */
         typedef std::vector<Factor>::const_iterator const_iterator;
 
     private:
@@ -126,6 +130,9 @@ class FactorGraph {
         const std::vector<Var>& vars() const { return _vars; }
 
         /// Returns reference to \a I 'th factor
+        /** \deprecated Please use the const member dai::FactorGraph::factor( size_t ) instead,
+         *  or dai::FactorGraph::setFactor() for changing a factor.
+         */
         Factor& factor( size_t I ) { 
             DAI_DEBASSERT( I < nrFactors() );
             return _factors[I]; 
@@ -152,12 +159,20 @@ class FactorGraph {
     /// \name Iterator interface
     //@{
         /// Returns iterator pointing to first factor
+        /** \deprecated The FactorGraph iterator interface will be removed in future versions of libDAI
+         */
         iterator begin() { return _factors.begin(); }
         /// Returns constant iterator pointing to first factor
+        /** \deprecated The FactorGraph iterator interface will be removed in future versions of libDAI
+         */
         const_iterator begin() const { return _factors.begin(); }
         /// Returns iterator pointing beyond last factor
+        /** \deprecated The FactorGraph iterator interface will be removed in future versions of libDAI
+         */
         iterator end() { return _factors.end(); }
         /// Returns constant iterator pointing beyond last factor
+        /** \deprecated The FactorGraph iterator interface will be removed in future versions of libDAI
+         */
         const_iterator end() const { return _factors.end(); }
     //@}
 
@@ -346,13 +361,13 @@ class FactorGraph {
          *  \throw CANNOT_READ_FILE if the file cannot be opened
          *  \throw INVALID_FACTORGRAPH_FILE if the file is not valid
          */
-        void ReadFromFile( const char *filename );
+        virtual void ReadFromFile( const char *filename );
 
         /// Writes a factor graph to a file
         /** \see \ref fileformats-factorgraph
          *  \throw CANNOT_WRITE_FILE if the file cannot be written
          */
-        void WriteToFile( const char *filename, size_t precision=15 ) const;
+        virtual void WriteToFile( const char *filename, size_t precision=15 ) const;
 
         /// Writes a factor graph to an output stream
         /** \see \ref fileformats-factorgraph
@@ -366,7 +381,7 @@ class FactorGraph {
         friend std::istream& operator>> (std::istream& is, FactorGraph& fg );
 
         /// Writes a factor graph to a GraphViz .dot file
-        void printDot( std::ostream& os ) const;
+        virtual void printDot( std::ostream& os ) const;
     //@}
 
     private:
index cd6dd4d..b72a2b4 100644 (file)
@@ -133,10 +133,11 @@ class GraphAL {
          *  \param nr The number of nodes.
          *  \param begin Points to the first edge.
          *  \param end Points just beyond the last edge.
+         *  \param check Whether to only add an edge if it does not exist already.
          */
         template<typename EdgeInputIterator>
-        GraphAL( size_t nr, EdgeInputIterator begin, EdgeInputIterator end ) : _nb() {
-            construct( nr, begin, end );
+        GraphAL( size_t nr, EdgeInputIterator begin, EdgeInputIterator end, bool check=true ) : _nb() {
+            construct( nr, begin, end, check );
         }
     //@}
 
@@ -174,9 +175,10 @@ class GraphAL {
          *  \param nr The number of nodes.
          *  \param begin Points to the first edge.
          *  \param end Points just beyond the last edge.
+         *  \param check Whether to only add an edge if it does not exist already.
          */
         template<typename EdgeInputIterator>
-        void construct( size_t nr, EdgeInputIterator begin, EdgeInputIterator end );
+        void construct( size_t nr, EdgeInputIterator begin, EdgeInputIterator end, bool check=true );
 
         /// Adds a node without neighbors and returns the index of the added node.
         size_t addNode() { _nb.push_back( Neighbors() ); return _nb.size() - 1; }
@@ -300,17 +302,12 @@ class GraphAL {
 
 
 template<typename EdgeInputIterator>
-void GraphAL::construct( size_t nr, EdgeInputIterator begin, EdgeInputIterator end ) {
+void GraphAL::construct( size_t nr, EdgeInputIterator begin, EdgeInputIterator end, bool check ) {
     _nb.clear();
     _nb.resize( nr );
 
-    for( EdgeInputIterator e = begin; e != end; e++ ) {
-#ifdef DAI_DEBUG
-        addEdge( e->first, e->second, true );
-#else
-        addEdge( e->first, e->second, false );
-#endif
-    }
+    for( EdgeInputIterator e = begin; e != end; e++ )
+        addEdge( e->first, e->second, check );
 }
 
 
index bfb01b8..891a455 100644 (file)
@@ -12,6 +12,7 @@
 /// \file
 /// \brief Defines class HAK, which implements a variant of Generalized Belief Propagation.
 /// \idea Implement more general region graphs and corresponding Generalized Belief Propagation updates as described in [\ref YFW05].
+/// \todo Implement TRW / FBP using HAK.
 
 
 #ifndef __defined_libdai_hak_h
index 0006fa8..94e4448 100644 (file)
@@ -30,20 +30,20 @@ namespace dai {
 class Region : public VarSet {
     private:
         /// Counting number
-        Real          _c;
+        Real _c;
 
     public:
         /// Default constructor
         Region() : VarSet(), _c(1.0) {}
 
         /// Construct from a set of variables and a counting number
-        Region( const VarSet &x, Real c ) : VarSet(x), _c(c) {}
+        Region( const VarSetx, Real c ) : VarSet(x), _c(c) {}
 
         /// Returns constant reference to counting number
-        const Real & c() const { return _c; }
+        const Real& c() const { return _c; }
 
         /// Returns reference to counting number
-        Real & c() { return _c; }
+        Real& c() { return _c; }
 };
 
 
@@ -58,13 +58,13 @@ class FRegion : public Factor {
         FRegion() : Factor(), _c(1.0) {}
 
         /// Constructs from a factor and a counting number
-        FRegion( const Factor & x, Real c ) : Factor(x), _c(c) {}
+        FRegion( const Factor& x, Real c ) : Factor(x), _c(c) {}
 
         /// Returns constant reference to counting number
-        const Real & c() const { return _c; }
+        const Real& c() const { return _c; }
 
         /// Returns reference to counting number
-        Real & c() { return _c; }
+        Real& c() { return _c; }
 };
 
 
@@ -83,6 +83,13 @@ class FRegion : public Factor {
  *
  *  Each factor in the factor graph belongs to an outer region; normally, the factor contents
  *  of an outer region would be the product of all the factors that belong to that region.
+ *  \idea Generalize the definition of region graphs to the one given in [\ref YFW05], i.e., replace
+ *  the current implementation which uses a BipartiteGraph with one that uses a DAG.
+ *  \idea The outer regions are products of factors; right now, this product is constantly cached:
+ *  changing one factor results in an update of all relevant outer regions. This may not be the most
+ *  efficient approach; an alternative would be to only precompute the factor products at the start
+ *  of an inference algorithm - e.g., in init(). This has the additional advantage that FactorGraph
+ e  can offer write access to its factors.
  */
 class RegionGraph : public FactorGraph {
     protected:
@@ -142,7 +149,7 @@ class RegionGraph : public FactorGraph {
         virtual RegionGraph* clone() const { return new RegionGraph(*this); }
     //@}
 
-    /// \name Queries
+    /// \name Accessors and mutators
     //@{
         /// Returns number of outer regions
         size_t nrORs() const { return _ORs.size(); }
@@ -180,9 +187,20 @@ class RegionGraph : public FactorGraph {
 
         /// Returns constant reference to the neighbors of outer region \a 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); }
 
+        /// Returns DAG structure of the region graph
+        /** \note Currently, the DAG is implemented as a BipartiteGraph; the nodes of
+         *  type 1 are the outer regions, the nodes of type 2 the inner regions, and
+         *  edges correspond with arrows from nodes of type 1 to type 2.
+         */
+        const BipartiteGraph& DAG() const { return _G; }
+    //@}
+
+    /// \name Queries
+    //@{
         /// Check whether the counting numbers are valid
         /** Counting numbers are said to be (variable) valid if for each variable \f$x\f$,
          *    \f[\sum_{\alpha \ni x} c_\alpha + \sum_{\beta \ni x} c_\beta = 1\f]
@@ -197,7 +215,7 @@ class RegionGraph : public FactorGraph {
         /// 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 ) {
             FactorGraph::setFactor( I, newFactor, backup );
-            RecomputeOR( I );
+            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
@@ -206,38 +224,63 @@ class RegionGraph : public FactorGraph {
             VarSet ns;
             for( std::map<size_t, Factor>::const_iterator fac = facs.begin(); fac != facs.end(); fac++ )
                 ns |= fac->second.vars();
-            RecomputeORs( ns );
+            recomputeORs( ns );
         }
 
+        /// Calculates counting numbers of inner regions based upon counting numbers of outer regions
+        /** The counting numbers of the inner regions are set using the Moebius inversion formula:
+         *    \f[ c_\beta := 1 - \sum_{\gamma \in \mathrm{an}(\beta)} c_\gamma \f]
+         *  where \f$\mathrm{an}(\beta)\f$ are the ancestors of inner region \f$\beta\f$ according to
+         *  the partial ordering induced by the subset relation (i.e., a region is a child of another
+         *  region if its variables are a subset of the variables of its parent region).
+         *  \deprecated This functionality has been protected.
+         */
+        void calcCountingNumbers() { calcCVMCountingNumbers(); }
+
         /// Recompute all outer regions
         /** The factor contents of each outer region is set to the product of the factors belonging to that region.
+         *  \deprecated This functionality has been protected.
          */
-        void RecomputeORs();
+        void RecomputeORs() { recomputeORs(); }
 
         /// 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.
+         *  \deprecated This functionality has been protected.
          */
-        void RecomputeORs( const VarSet& vs );
+        void RecomputeORs( const VarSet& vs ) { recomputeORs( vs ); }
 
         /// 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.
+         *  \deprecated This functionality has been protected.
          */
-        void RecomputeOR( size_t I );
-
-        /// Calculates counting numbers of inner regions based upon counting numbers of outer regions
-        /** The counting numbers of the inner regions are set using the Moebius inversion formula:
-         *    \f[ c_\beta := 1 - \sum_{\gamma \in \mathrm{an}(\beta)} c_\gamma \f]
-         *  where \f$\mathrm{an}(\beta)\f$ are the ancestors of inner region \f$\beta\f$ according to
-         *  the partial ordering induced by the subset relation (i.e., a region is a child of another
-         *  region if its variables are a subset of the variables of its parent region).
-         */
-        void calcCountingNumbers();
+        void RecomputeOR( size_t I ) { recomputeOR( I ); }
     //@}
 
     /// \name Input/output
     //@{
+        /// Reads a region graph from a file
+        /** \note Not implemented yet
+         */
+        virtual void ReadFromFile( const char* /*filename*/ ) {
+            DAI_THROW(NOT_IMPLEMENTED);
+        }
+
+        /// Writes a factor graph to a file
+        /** \note Not implemented yet
+         */
+        virtual void WriteToFile( const char* /*filename*/, size_t /*precision*/=15 ) const {
+            DAI_THROW(NOT_IMPLEMENTED);
+        }
+
         /// 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 );
+
+        /// Writes a region graph to a GraphViz .dot file
+        /** \note Not implemented yet
+         */
+        virtual void printDot( std::ostream& /*os*/ ) const {
+            DAI_THROW(NOT_IMPLEMENTED);
+        }
     //@}
 
     protected:
@@ -246,6 +289,31 @@ class RegionGraph : public FactorGraph {
 
         /// Helper function for constructors (CVM style)
         void constructCVM( const FactorGraph& fg, const std::vector<VarSet>& cl );
+
+        /// Recompute all outer regions
+        /** The factor contents of each outer region is set to the product of the factors belonging to that region.
+         */
+        void recomputeORs();
+
+        /// 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 );
+
+        /// 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.
+         */
+        void recomputeOR( size_t I );
+
+        /// Calculates counting numbers of inner regions based upon counting numbers of outer regions
+        /** The counting numbers of the inner regions are set using the Moebius inversion formula:
+         *    \f[ c_\beta := 1 - \sum_{\gamma \in \mathrm{an}(\beta)} c_\gamma \f]
+         *  where \f$\mathrm{an}(\beta)\f$ are the ancestors of inner region \f$\beta\f$ according to
+         *  the partial ordering induced by the subset relation (i.e., a region is a child of another
+         *  region if its variables are a subset of the variables of its parent region).
+         */
+        void calcCVMCountingNumbers();
+
 };
 
 
index 74dea69..bfc936a 100644 (file)
@@ -57,30 +57,30 @@ class Var {
         size_t& label() { return _label; }
 
         /// Returns the number of states
-        size_t states () const { return _states; }
+        size_t states() const { return _states; }
         /// Returns reference to number of states
-        size_t& states () { return _states; }
+        size_t& states() { return _states; }
 
         /// Smaller-than operator (only compares labels)
-        bool operator < ( const Var& n ) const { 
+        bool operator< ( const Var& n ) const { 
 #ifdef DAI_DEBUG
             if( _label == n._label )
                 DAI_ASSERT( _states == n._states );
 #endif
-            return( _label <  n._label );
+            return( _label < n._label );
         }
 
         /// Larger-than operator (only compares labels)
-        bool operator > ( const Var& n ) const { 
+        bool operator> ( const Var& n ) const { 
 #ifdef DAI_DEBUG
             if( _label == n._label )
                 DAI_ASSERT( _states == n._states );
 #endif
-            return( _label >  n._label ); 
+            return( _label > n._label ); 
         }
 
         /// Smaller-than-or-equal-to operator (only compares labels)
-        bool operator <= ( const Var& n ) const {
+        bool operator<= ( const Var& n ) const {
 #ifdef DAI_DEBUG
             if( _label == n._label )
                 DAI_ASSERT( _states == n._states );
@@ -89,7 +89,7 @@ class Var {
         }
 
         /// Larger-than-or-equal-to operator (only compares labels)
-        bool operator >= ( const Var& n ) const {
+        bool operator>= ( const Var& n ) const {
 #ifdef DAI_DEBUG
             if( _label == n._label )
                 DAI_ASSERT( _states == n._states );
@@ -98,7 +98,7 @@ class Var {
         }
 
         /// Not-equal-to operator (only compares labels)
-        bool operator != ( const Var& n ) const {
+        bool operator!= ( const Var& n ) const {
 #ifdef DAI_DEBUG
             if( _label == n._label )
                 DAI_ASSERT( _states == n._states );
@@ -107,7 +107,7 @@ class Var {
         }
 
         /// Equal-to operator (only compares labels)
-        bool operator == ( const Var& n ) const {
+        bool operator== ( const Var& n ) const {
 #ifdef DAI_DEBUG
             if( _label == n._label )
                 DAI_ASSERT( _states == n._states );
index fd58358..2333850 100644 (file)
@@ -11,7 +11,7 @@
 
 /** \file
  *  \brief Defines some utility functions for (weighted) undirected graphs, trees and rooted trees.
- *  \todo Improve general support for graphs and trees.
+ *  \todo Improve general support for graphs and trees: implement a DAG and a Tree.
  */
 
 
index 5a4049d..9bf9741 100644 (file)
@@ -1014,7 +1014,8 @@ Real numericBBPTest( const InfAlg &bp, const std::vector<size_t> *state, const P
                 psi_1_prb.set( xi, psi_1_prb[xi] + h );
 //                 psi_1_prb.normalize();
                 size_t I = bp_prb->fg().nbV(i)[0]; // use first factor in list of neighbors of i
-                bp_prb->fg().factor(I) *= Factor( bp_prb->fg().var(i), psi_1_prb );
+                Factor tmp = bp_prb->fg().factor(I) * Factor( bp_prb->fg().var(i), psi_1_prb );
+                bp_prb->fg().setFactor( I, tmp );
 
                 // call 'init' on the perturbed variables
                 bp_prb->init( bp_prb->fg().var(i) );
index 4389b45..e7cf3c4 100644 (file)
@@ -289,7 +289,7 @@ bool BipartiteGraph::isTree() const {
 
 
 void BipartiteGraph::printDot( std::ostream& os ) const {
-    os << "graph G {" << endl;
+    os << "graph BipartiteGraph {" << endl;
     os << "node[shape=circle,width=0.4,fixedsize=true];" << endl;
     for( size_t n1 = 0; n1 < nrNodes1(); n1++ )
         os << "\tx" << n1 << ";" << endl;
index cc25844..15b3adf 100644 (file)
@@ -253,7 +253,7 @@ void FactorGraph::WriteToFile( const char *filename, size_t precision ) const {
 
 
 void FactorGraph::printDot( std::ostream &os ) const {
-    os << "graph G {" << endl;
+    os << "graph FactorGraph {" << endl;
     os << "node[shape=circle,width=0.4,fixedsize=true];" << endl;
     for( size_t i = 0; i < nrVars(); i++ )
         os << "\tv" << var(i).label() << ";" << endl;
index 08537b1..050d03f 100644 (file)
@@ -150,7 +150,7 @@ void JTree::construct( const FactorGraph &fg, const std::vector<VarSet> &cl, boo
         if( verify )
             DAI_ASSERT( alpha != nrORs() );
     }
-    RecomputeORs();
+    recomputeORs();
 
     // Create inner regions and edges
     _IRs.clear();
index bd955ad..05c4c16 100644 (file)
@@ -49,7 +49,7 @@ void RegionGraph::construct( const FactorGraph &fg, const std::vector<VarSet> &o
             }
         DAI_ASSERT( alpha != nrORs() );
     }
-    RecomputeORs();
+    recomputeORs();
 
     // Create bipartite graph
     _G.construct( nrORs(), nrIRs(), edges.begin(), edges.end() );
@@ -100,11 +100,11 @@ void RegionGraph::constructCVM( const FactorGraph &fg, const std::vector<VarSet>
     construct( fg, cg.clusters(), irs, edges );
 
     // Calculate counting numbers
-    calcCountingNumbers();
+    calcCVMCountingNumbers();
 }
 
 
-void RegionGraph::calcCountingNumbers() {
+void RegionGraph::calcCVMCountingNumbers() {
     // Calculates counting numbers of inner regions based upon counting numbers of outer regions
 
     vector<vector<size_t> > ancestors(nrIRs());
@@ -163,7 +163,7 @@ bool RegionGraph::checkCountingNumbers() const {
 }
 
 
-void RegionGraph::RecomputeORs() {
+void RegionGraph::recomputeORs() {
     for( size_t alpha = 0; alpha < nrORs(); alpha++ )
         OR(alpha).fill( 1.0 );
     for( size_t I = 0; I < nrFactors(); I++ )
@@ -172,7 +172,7 @@ void RegionGraph::RecomputeORs() {
 }
 
 
-void RegionGraph::RecomputeORs( const VarSet &ns ) {
+void RegionGraph::recomputeORs( const VarSet &ns ) {
     for( size_t alpha = 0; alpha < nrORs(); alpha++ )
         if( OR(alpha).vars().intersects( ns ) )
             OR(alpha).fill( 1.0 );
@@ -183,7 +183,7 @@ void RegionGraph::RecomputeORs( const VarSet &ns ) {
 }
 
 
-void RegionGraph::RecomputeOR( size_t I ) {
+void RegionGraph::recomputeOR( size_t I ) {
     DAI_ASSERT( I < nrFactors() );
     if( fac2OR(I) != -1U ) {
         size_t alpha = fac2OR(I);
@@ -197,20 +197,18 @@ void RegionGraph::RecomputeOR( size_t I ) {
 
 /// Send RegionGraph to output stream
 ostream & operator << (ostream & os, const RegionGraph & rg) {
-    os << "Outer regions" << endl;
+    os << "digraph RegionGraph {" << endl;
+    os << "node[shape=box];" << endl;
     for( size_t alpha = 0; alpha < rg.nrORs(); alpha++ )
-        os << alpha << ": " << rg.OR(alpha).vars() << ": c = " << rg.OR(alpha).c() << endl;
-
-    os << "Inner regions" << endl;
+        os << "\ta" << alpha << " [label=\"a" << alpha << ": " << rg.OR(alpha).vars() << ", c=" << rg.OR(alpha).c() << "\"];" << endl;
+    os << "node[shape=ellipse];" << endl;
     for( size_t beta = 0; beta < rg.nrIRs(); beta++ )
-        os << beta << ": " << (VarSet)rg.IR(beta) << ": c = " << rg.IR(beta).c() << endl;
-
-    os << "Edges" << endl;
+        os << "\tb" << beta << " [label=\"b" << beta << ": " << (VarSet)rg.IR(beta) << ", c=" << rg.IR(beta).c() << "\"];" << endl;
     for( size_t alpha = 0; alpha < rg.nrORs(); alpha++ )
         foreach( const RegionGraph::Neighbor &beta, rg.nbOR(alpha) )
-            os << alpha << "->" << beta << endl;
-
-    return(os);
+            os << "\ta" << alpha << " -> b" << beta << ";" << endl;
+    os << "}" << endl;
+    return os;
 }
 
 
index 2cd2fc6..c72d271 100644 (file)
@@ -387,33 +387,35 @@ BOOST_AUTO_TEST_CASE( StreamTest ) {
     BipartiteGraph G( 2, 3, edges.begin(), edges.end() );
 
     std::stringstream ss;
+    std::string s;
+
     G.printDot( ss );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "graph BipartiteGraph {" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "node[shape=circle,width=0.4,fixedsize=true];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1;" );
+    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, "\ty0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ty1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ty2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0 -- y0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0 -- y1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1 -- y1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1 -- y2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "}" );
 
-    std::string s;
-    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, "\tx0;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx1;" );
-    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, "\ty0;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\ty1;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\ty2;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx0 -- y0;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx0 -- y1;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx1 -- y1;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx1 -- y2;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "}" );
+    ss << G;
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "graph BipartiteGraph {" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "node[shape=circle,width=0.4,fixedsize=true];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1;" );
+    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, "\ty0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ty1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ty2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0 -- y0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0 -- y1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1 -- y1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1 -- y2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "}" );
 }
index dcdecf1..2f84617 100644 (file)
@@ -35,8 +35,10 @@ BOOST_AUTO_TEST_CASE( ConstructorsTest ) {
     BOOST_CHECK( G.bipGraph() == BipartiteGraph() );
     BOOST_CHECK_EQUAL( G.nrVars(), 0 );
     BOOST_CHECK_EQUAL( G.nrClusters(), 0 );
+#ifdef DAI_DEBUG
     BOOST_CHECK_THROW( G.var( 0 ), Exception );
     BOOST_CHECK_THROW( G.cluster( 0 ), Exception );
+#endif
     BOOST_CHECK_THROW( G.findVar( Var( 0, 2 ) ), Exception );
 
     Var v0( 0, 2 );
index 497ca71..618ff83 100644 (file)
@@ -68,7 +68,7 @@ BOOST_AUTO_TEST_CASE( ConstructorsTest ) {
 }
 
 
-BOOST_AUTO_TEST_CASE( AccMutIterTest ) {
+BOOST_AUTO_TEST_CASE( AccMutTest ) {
     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) ) ) );
@@ -110,15 +110,6 @@ BOOST_AUTO_TEST_CASE( AccMutIterTest ) {
     BOOST_CHECK_EQUAL( G.nbF(2,1), 2 );
     BOOST_CHECK_EQUAL( G.nbF(3).size(), 1 );
     BOOST_CHECK_EQUAL( G.nbF(3,0), 1 );
-
-    FactorGraph::const_iterator cit = G.begin();
-    FactorGraph::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() );
 }
 
 
@@ -674,7 +665,7 @@ BOOST_AUTO_TEST_CASE( IOTest ) {
     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, "graph FactorGraph {" );
     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;" );
@@ -690,9 +681,9 @@ BOOST_AUTO_TEST_CASE( IOTest ) {
     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);
+    G.setFactor( 0, Factor( G.factor(0).vars(), 1.0 ) );
+    G.setFactor( 1, Factor( G.factor(1).vars(), 2.0 ) );
+    G.setFactor( 2, Factor( G.factor(2).vars(), 3.0 ) );
     ss << G;
     std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "3" );
     std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "" );
index 4d8471f..4e1aff8 100644 (file)
@@ -295,9 +295,9 @@ BOOST_AUTO_TEST_CASE( QueriesAndCreationTest ) {
         }
 
     // createGraphGrid3D
-    for( size_t N1 = 0; N1 < 10; N1++ )
-        for( size_t N2 = 0; N2 < 10; N2++ )
-            for( size_t N3 = 0; N3 < 10; N3++ ) {
+    for( size_t N1 = 0; N1 < 8; N1++ )
+        for( size_t N2 = 0; N2 < 8; N2++ )
+            for( size_t N3 = 0; N3 < 8; N3++ ) {
                 GraphAL G = createGraphGrid3D( N1, N2, N3, false );
                 BOOST_CHECK_EQUAL( G.nrNodes(), N1 * N2 * N3 );
                 BOOST_CHECK_EQUAL( G.nrEdges(), (N1 > 0 && N2 > 0 && N3 > 0) ? 3 * (N1-1) * (N2-1) * (N3-1) + 2 * (N1-1) * (N2-1) + 2 * (N1-1) * (N3-1) + 2 *  (N2-1) * (N3-1) + (N1-1) + (N2-1) + (N3-1) : 0 );
@@ -394,8 +394,8 @@ BOOST_AUTO_TEST_CASE( QueriesAndCreationTest ) {
     }
 
     // createGraphRegular
-    for( size_t N = 0; N < 100; N++ ) {
-        for( size_t d = 0; d < N && d <= 20; d++ ) {
+    for( size_t N = 0; N < 50; N++ ) {
+        for( size_t d = 0; d < N && d <= 15; d++ ) {
             if( (N * d) % 2 == 0 ) {
                 GraphAL G = createGraphRegular( N, d );
                 BOOST_CHECK_EQUAL( G.nrNodes(), N );
@@ -435,26 +435,15 @@ BOOST_AUTO_TEST_CASE( StreamTest ) {
     G.printDot( ss );
 
     std::string s;
-    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, "\tx0;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx1;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx2;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx3;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx0 -- x1;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx0 -- x2;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx1 -- x3;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "\tx2 -- x3;" );
-    std::getline( ss, s );
-    BOOST_CHECK_EQUAL( s, "}" );
+    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, "\tx0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx3;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0 -- x1;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx0 -- x2;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx1 -- x3;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tx2 -- x3;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "}" );
 }
index fdd55c0..ef193c2 100644 (file)
@@ -30,62 +30,192 @@ const double tol = 1e-8;
 
 
 BOOST_AUTO_TEST_CASE( ConstructorsTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
     RegionGraph G;
     BOOST_CHECK_EQUAL( G.vars(), std::vector<Var>() );
     BOOST_CHECK_EQUAL( G.factors(), std::vector<Factor>() );
+    BOOST_CHECK( G.bipGraph() == BipartiteGraph() );
+    BOOST_CHECK( G.DAG() == BipartiteGraph() );
+    BOOST_CHECK_EQUAL( G.nrORs(), 0 );
+    BOOST_CHECK_EQUAL( G.nrIRs(), 0 );
 
     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) ) ) );
+    facs.push_back( Factor( v01 ) );
+    facs.push_back( Factor( v02 ) );
+    facs.push_back( Factor( v12 ) );
+    facs.push_back( Factor( v1 ) );
     std::vector<Var> vars;
-    vars.push_back( Var( 0, 2 ) );
-    vars.push_back( Var( 1, 2 ) );
-    vars.push_back( Var( 2, 2 ) );
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+    BipartiteGraph dag( 3, 3 );
+    dag.addEdge( 0, 0 );
+    dag.addEdge( 0, 1 );
+    dag.addEdge( 1, 0 );
+    dag.addEdge( 1, 2 );
+    dag.addEdge( 2, 1 );
+    dag.addEdge( 2, 2 );
+    BipartiteGraph bipgraph( dag );
+    bipgraph.addNode2();
+    bipgraph.addEdge( 1, 3 );
 
     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_CHECK( G1.bipGraph() == bipgraph );
+    BOOST_CHECK( G1.DAG() == dag );
+    BOOST_CHECK_EQUAL( G1.nrORs(), 3 );
+    BOOST_CHECK_EQUAL( G1.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G1.OR(0), facs[0] * facs[3] );
+    BOOST_CHECK_EQUAL( G1.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G1.OR(1), facs[1] );
+    BOOST_CHECK_EQUAL( G1.OR(2).c(), 1 );
+    BOOST_CHECK_EQUAL( G1.OR(2), facs[2] );
+    BOOST_CHECK_EQUAL( G1.nrIRs(), 3 );
+    BOOST_CHECK_EQUAL( G1.IR(0).c(), -1 );
+    BOOST_CHECK_EQUAL( G1.IR(0), v0 );
+    BOOST_CHECK_EQUAL( G1.IR(1).c(), -1 );
+    BOOST_CHECK_EQUAL( G1.IR(1), v1 );
+    BOOST_CHECK_EQUAL( G1.IR(2).c(), -1 );
+    BOOST_CHECK_EQUAL( G1.IR(2), v2 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(1), 1 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(2), 2 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(3), 0 );
+    BOOST_CHECK( G1.checkCountingNumbers() );
+
+    std::vector<VarSet> ors;
+    std::vector<Region> irs;
+    ors.push_back( v01 | v2 );
+    irs.push_back( Region( v01, 0.0 ) );
+    irs.push_back( Region( v2, 0.0 ) );
+    typedef std::pair<size_t,size_t> edge;
+    std::vector<edge> edges;
+    edges.push_back( edge( 0, 0 ) );
+    edges.push_back( edge( 0, 1 ) );
+    BipartiteGraph dag2( 1, 2, edges.begin(), edges.end() );
+    RegionGraph G2( fg, ors, irs, edges );
+    BOOST_CHECK_EQUAL( G2.vars(), vars );
+    BOOST_CHECK_EQUAL( G2.factors(), facs );
+    BOOST_CHECK( G2.bipGraph() == bipgraph );
+    BOOST_CHECK( G2.DAG() == dag2 );
+    BOOST_CHECK_EQUAL( G2.nrORs(), 1 );
+    BOOST_CHECK_EQUAL( G2.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G2.OR(0), facs[0] * facs[1] * facs[2] * facs[3] );
+    BOOST_CHECK_EQUAL( G2.nrIRs(), 2 );
+    BOOST_CHECK_EQUAL( G2.IR(0), irs[0] );
+    BOOST_CHECK_EQUAL( G2.IR(1), irs[1] );
+    BOOST_CHECK_EQUAL( G2.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(1), 0 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(2), 0 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(3), 0 );
+    BOOST_CHECK( G2.checkCountingNumbers() );
+
+    RegionGraph *G3 = G1.clone();
+    BOOST_CHECK_EQUAL( G3->vars(), G1.vars() );
+    BOOST_CHECK_EQUAL( G3->factors(), G1.factors() );
+    BOOST_CHECK( G3->bipGraph() == bipgraph );
+    BOOST_CHECK( G3->DAG() == dag );
+    BOOST_CHECK_EQUAL( G3->nrORs(), G1.nrORs() );
+    BOOST_CHECK_EQUAL( G3->OR(0), G1.OR(0) );
+    BOOST_CHECK_EQUAL( G3->OR(1), G1.OR(1) );
+    BOOST_CHECK_EQUAL( G3->OR(2), G1.OR(2) );
+    BOOST_CHECK_EQUAL( G3->nrIRs(), G1.nrIRs() );
+    BOOST_CHECK_EQUAL( G3->IR(0), G1.IR(0) );
+    BOOST_CHECK_EQUAL( G3->IR(1), G1.IR(1) );
+    BOOST_CHECK_EQUAL( G3->IR(2), G1.IR(2) );
+    BOOST_CHECK_EQUAL( G3->fac2OR(0), G1.fac2OR(0) );
+    BOOST_CHECK_EQUAL( G3->fac2OR(1), G1.fac2OR(1) );
+    BOOST_CHECK_EQUAL( G3->fac2OR(2), G1.fac2OR(2) );
+    BOOST_CHECK_EQUAL( G3->fac2OR(3), G1.fac2OR(3) );
+    BOOST_CHECK( G3->checkCountingNumbers() );
+    delete G3;
+
+    RegionGraph G4 = G1;
+    BOOST_CHECK_EQUAL( G4.vars(), G1.vars() );
+    BOOST_CHECK_EQUAL( G4.factors(), G1.factors() );
+    BOOST_CHECK( G4.bipGraph() == bipgraph );
+    BOOST_CHECK( G4.DAG() == dag );
+    BOOST_CHECK_EQUAL( G4.nrORs(), G1.nrORs() );
+    BOOST_CHECK_EQUAL( G4.OR(0), G1.OR(0) );
+    BOOST_CHECK_EQUAL( G4.OR(1), G1.OR(1) );
+    BOOST_CHECK_EQUAL( G4.OR(2), G1.OR(2) );
+    BOOST_CHECK_EQUAL( G4.nrIRs(), G1.nrIRs() );
+    BOOST_CHECK_EQUAL( G4.IR(0), G1.IR(0) );
+    BOOST_CHECK_EQUAL( G4.IR(1), G1.IR(1) );
+    BOOST_CHECK_EQUAL( G4.IR(2), G1.IR(2) );
+    BOOST_CHECK_EQUAL( G4.fac2OR(0), G1.fac2OR(0) );
+    BOOST_CHECK_EQUAL( G4.fac2OR(1), G1.fac2OR(1) );
+    BOOST_CHECK_EQUAL( G4.fac2OR(2), G1.fac2OR(2) );
+    BOOST_CHECK_EQUAL( G4.fac2OR(3), G1.fac2OR(3) );
+    BOOST_CHECK( G4.checkCountingNumbers() );
+
+    RegionGraph G5( G1 );
+    BOOST_CHECK_EQUAL( G5.vars(), G1.vars() );
+    BOOST_CHECK_EQUAL( G5.factors(), G1.factors() );
+    BOOST_CHECK( G5.bipGraph() == bipgraph );
+    BOOST_CHECK( G5.DAG() == dag );
+    BOOST_CHECK_EQUAL( G5.nrORs(), G1.nrORs() );
+    BOOST_CHECK_EQUAL( G5.OR(0), G1.OR(0) );
+    BOOST_CHECK_EQUAL( G5.OR(1), G1.OR(1) );
+    BOOST_CHECK_EQUAL( G5.OR(2), G1.OR(2) );
+    BOOST_CHECK_EQUAL( G5.nrIRs(), G1.nrIRs() );
+    BOOST_CHECK_EQUAL( G5.IR(0), G1.IR(0) );
+    BOOST_CHECK_EQUAL( G5.IR(1), G1.IR(1) );
+    BOOST_CHECK_EQUAL( G5.IR(2), G1.IR(2) );
+    BOOST_CHECK_EQUAL( G5.fac2OR(0), G1.fac2OR(0) );
+    BOOST_CHECK_EQUAL( G5.fac2OR(1), G1.fac2OR(1) );
+    BOOST_CHECK_EQUAL( G5.fac2OR(2), G1.fac2OR(2) );
+    BOOST_CHECK_EQUAL( G5.fac2OR(3), G1.fac2OR(3) );
+    BOOST_CHECK( G5.checkCountingNumbers() );
 }
 
 
-BOOST_AUTO_TEST_CASE( AccMutIterTest ) {
+BOOST_AUTO_TEST_CASE( AccMutTest ) {
+    Var v0( 0, 2 );
+    Var v1( 1, 2 );
+    Var v2( 2, 2 );
+    VarSet v01( v0, v1 );
+    VarSet v02( v0, v2 );
+    VarSet v12( v1, v2 );
     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) ) ) );
+    facs.push_back( Factor( v01 ) );
+    facs.push_back( Factor( v02 ) );
+    facs.push_back( Factor( v12 ) );
+    facs.push_back( Factor( v1 ) );
     std::vector<Var> vars;
-    vars.push_back( Var( 0, 2 ) );
-    vars.push_back( Var( 1, 2 ) );
-    vars.push_back( Var( 2, 2 ) );
+    vars.push_back( v0 );
+    vars.push_back( v1 );
+    vars.push_back( v2 );
+    BipartiteGraph dag( 3, 3 );
+    dag.addEdge( 0, 0 );
+    dag.addEdge( 0, 1 );
+    dag.addEdge( 1, 0 );
+    dag.addEdge( 1, 2 );
+    dag.addEdge( 2, 1 );
+    dag.addEdge( 2, 2 );
+    BipartiteGraph bipgraph( dag );
+    bipgraph.addNode2();
+    bipgraph.addEdge( 1, 3 );
 
     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.var(0), v0 );
+    BOOST_CHECK_EQUAL( G.var(1), v1 );
+    BOOST_CHECK_EQUAL( G.var(2), v2 );
     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( G.bipGraph() == bipgraph );
     BOOST_CHECK_EQUAL( G.nbV(0).size(), 2 );
     BOOST_CHECK_EQUAL( G.nbV(0,0), 0 );
     BOOST_CHECK_EQUAL( G.nbV(0,1), 1 );
@@ -107,15 +237,46 @@ BOOST_AUTO_TEST_CASE( AccMutIterTest ) {
     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_CHECK( G.DAG() == dag );
+    BOOST_CHECK_EQUAL( G.nrORs(), 3 );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0).vars(), v01 );
+    BOOST_CHECK_EQUAL( G.OR(0), facs[0] * facs[3] );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1).vars(), v02 );
+    BOOST_CHECK_EQUAL( G.OR(1), facs[1] );
+    BOOST_CHECK_EQUAL( G.OR(2).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(2).vars(), v12 );
+    BOOST_CHECK_EQUAL( G.OR(2), facs[2] );
+    BOOST_CHECK_EQUAL( G.nrIRs(), 3 );
+    BOOST_CHECK_EQUAL( G.IR(0).c(), -1 );
+    BOOST_CHECK_EQUAL( G.IR(0), v0 );
+    BOOST_CHECK_EQUAL( G.IR(1).c(), -1 );
+    BOOST_CHECK_EQUAL( G.IR(1), v1 );
+    BOOST_CHECK_EQUAL( G.IR(2).c(), -1 );
+    BOOST_CHECK_EQUAL( G.IR(2), v2 );
+    BOOST_CHECK_EQUAL( G.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G.fac2OR(1), 1 );
+    BOOST_CHECK_EQUAL( G.fac2OR(2), 2 );
+    BOOST_CHECK_EQUAL( G.fac2OR(3), 0 );
+    BOOST_CHECK_EQUAL( G.nbOR(0).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbOR(0)[0], 0 );
+    BOOST_CHECK_EQUAL( G.nbOR(0)[1], 1 );
+    BOOST_CHECK_EQUAL( G.nbOR(1).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbOR(1)[0], 0 );
+    BOOST_CHECK_EQUAL( G.nbOR(1)[1], 2 );
+    BOOST_CHECK_EQUAL( G.nbOR(2).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbOR(2)[0], 1 );
+    BOOST_CHECK_EQUAL( G.nbOR(2)[1], 2 );
+    BOOST_CHECK_EQUAL( G.nbIR(0).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbIR(0)[0], 0 );
+    BOOST_CHECK_EQUAL( G.nbIR(0)[1], 1 );
+    BOOST_CHECK_EQUAL( G.nbIR(1).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbIR(1)[0], 0 );
+    BOOST_CHECK_EQUAL( G.nbIR(1)[1], 2 );
+    BOOST_CHECK_EQUAL( G.nbIR(2).size(), 2 );
+    BOOST_CHECK_EQUAL( G.nbIR(2)[0], 1 );
+    BOOST_CHECK_EQUAL( G.nbIR(2)[1], 2 );
 }
 
 
@@ -140,6 +301,7 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     BOOST_CHECK_THROW( G0.Delta( 0 ), Exception );
     BOOST_CHECK_THROW( G0.delta( v0 ), Exception );
     BOOST_CHECK_THROW( G0.Delta( v0 ), Exception );
+    BOOST_CHECK_THROW( G0.fac2OR( 0 ), Exception );
 #endif
     BOOST_CHECK( G0.isConnected() );
     BOOST_CHECK( G0.isTree() );
@@ -148,6 +310,10 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     BOOST_CHECK( G0.MarkovGraph() == GraphAL() );
     BOOST_CHECK( G0.bipGraph() == BipartiteGraph() );
     BOOST_CHECK_EQUAL( G0.maximalFactorDomains().size(), 0 );
+    BOOST_CHECK( G0.DAG() == BipartiteGraph() );
+    BOOST_CHECK( G0.checkCountingNumbers() );
+    BOOST_CHECK_EQUAL( G0.nrORs(), 0 );
+    BOOST_CHECK_EQUAL( G0.nrIRs(), 0 );
 
     std::vector<Factor> facs;
     facs.push_back( Factor( v01 ) );
@@ -166,6 +332,9 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     K.addEdge( 1, 1 );
     K.addEdge( 2, 1 );
     K.addEdge( 1, 2 );
+    BipartiteGraph dag(2, 1);
+    dag.addEdge( 0, 0 );
+    dag.addEdge( 1, 0 );
 
     FactorGraph fg( facs );
     RegionGraph G1( fg, fg.maximalFactorDomains() );
@@ -211,12 +380,32 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     BOOST_CHECK_EQUAL( G1.maximalFactorDomains().size(), 2 );
     BOOST_CHECK_EQUAL( G1.maximalFactorDomains()[0], v01 );
     BOOST_CHECK_EQUAL( G1.maximalFactorDomains()[1], v12 );
+    BOOST_CHECK( G1.DAG() == dag );
+    BOOST_CHECK( G1.checkCountingNumbers() );
+    BOOST_CHECK_EQUAL( G1.nrORs(), 2 );
+    BOOST_CHECK_EQUAL( G1.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G1.OR(0), facs[0] * facs[2] );
+    BOOST_CHECK_EQUAL( G1.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G1.OR(1), facs[1] );
+    BOOST_CHECK_EQUAL( G1.nrIRs(), 1 );
+    BOOST_CHECK_EQUAL( G1.IR(0).c(), -1 );
+    BOOST_CHECK_EQUAL( G1.IR(0), v1 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(1), 1 );
+    BOOST_CHECK_EQUAL( G1.fac2OR(2), 0 );
 
     facs.push_back( Factor( v02 ) );
     H.addEdge( 0, 2 );
     K.addNode2();
     K.addEdge( 0, 3 );
     K.addEdge( 2, 3 );
+    dag = BipartiteGraph( 3, 3 );
+    dag.addEdge( 0, 0 );
+    dag.addEdge( 0, 1 );
+    dag.addEdge( 1, 1 );
+    dag.addEdge( 1, 2 );
+    dag.addEdge( 2, 0 );
+    dag.addEdge( 2, 2 );
     fg = FactorGraph( facs );
     RegionGraph G2( fg, fg.maximalFactorDomains() );
     BOOST_CHECK_EQUAL( G2.nrVars(), 3 );
@@ -248,6 +437,26 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     BOOST_CHECK_EQUAL( G2.maximalFactorDomains()[0], v01 );
     BOOST_CHECK_EQUAL( G2.maximalFactorDomains()[1], v12 );
     BOOST_CHECK_EQUAL( G2.maximalFactorDomains()[2], v02 );
+    BOOST_CHECK( G2.DAG() == dag );
+    BOOST_CHECK( G2.checkCountingNumbers() );
+    BOOST_CHECK_EQUAL( G2.nrORs(), 3 );
+    BOOST_CHECK_EQUAL( G2.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G2.OR(0), facs[0] * facs[2] );
+    BOOST_CHECK_EQUAL( G2.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G2.OR(1), facs[1] );
+    BOOST_CHECK_EQUAL( G2.OR(2).c(), 1 );
+    BOOST_CHECK_EQUAL( G2.OR(2), facs[3] );
+    BOOST_CHECK_EQUAL( G2.nrIRs(), 3 );
+    BOOST_CHECK_EQUAL( G2.IR(0).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G2.IR(0), v0 );
+    BOOST_CHECK_EQUAL( G2.IR(1).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G2.IR(1), v1 );
+    BOOST_CHECK_EQUAL( G2.IR(2).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G2.IR(2), v2 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(1), 1 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(2), 0 );
+    BOOST_CHECK_EQUAL( G2.fac2OR(3), 2 );
 
     Var v3( 3, 3 );
     VarSet v03( v0, v3 );
@@ -263,6 +472,7 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     K.addNode1();
     K.addNode2();
     K.addEdge( 3, 4 );
+    dag.addNode1();
     fg = FactorGraph( facs );
     RegionGraph G3( fg, fg.maximalFactorDomains() );
     BOOST_CHECK_EQUAL( G3.nrVars(), 4 );
@@ -300,6 +510,29 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[1], v12 );
     BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[2], v02 );
     BOOST_CHECK_EQUAL( G3.maximalFactorDomains()[3], v3 );
+    BOOST_CHECK( G3.DAG() == dag );
+    BOOST_CHECK( G3.checkCountingNumbers() );
+    BOOST_CHECK_EQUAL( G3.nrORs(), 4 );
+    BOOST_CHECK_EQUAL( G3.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G3.OR(0), facs[0] * facs[2] );
+    BOOST_CHECK_EQUAL( G3.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G3.OR(1), facs[1] );
+    BOOST_CHECK_EQUAL( G3.OR(2).c(), 1 );
+    BOOST_CHECK_EQUAL( G3.OR(2), facs[3] );
+    BOOST_CHECK_EQUAL( G3.OR(3).c(), 1 );
+    BOOST_CHECK_EQUAL( G3.OR(3), facs[4] );
+    BOOST_CHECK_EQUAL( G3.nrIRs(), 3 );
+    BOOST_CHECK_EQUAL( G3.IR(0).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G3.IR(0), v0 );
+    BOOST_CHECK_EQUAL( G3.IR(1).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G3.IR(1), v1 );
+    BOOST_CHECK_EQUAL( G3.IR(2).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G3.IR(2), v2 );
+    BOOST_CHECK_EQUAL( G3.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G3.fac2OR(1), 1 );
+    BOOST_CHECK_EQUAL( G3.fac2OR(2), 0 );
+    BOOST_CHECK_EQUAL( G3.fac2OR(3), 2 );
+    BOOST_CHECK_EQUAL( G3.fac2OR(4), 3 );
 
     facs.push_back( Factor( v123 ) );
     H.addEdge( 1, 3 );
@@ -308,6 +541,13 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     K.addEdge( 1, 5 );
     K.addEdge( 2, 5 );
     K.addEdge( 3, 5 );
+    dag = BipartiteGraph( 3, 3 );
+    dag.addEdge( 0, 0 );
+    dag.addEdge( 0, 1 );
+    dag.addEdge( 1, 0 );
+    dag.addEdge( 1, 2 );
+    dag.addEdge( 2, 1 );
+    dag.addEdge( 2, 2 );
     fg = FactorGraph( facs );
     RegionGraph G4( fg, fg.maximalFactorDomains() );
     BOOST_CHECK_EQUAL( G4.nrVars(), 4 );
@@ -345,6 +585,28 @@ BOOST_AUTO_TEST_CASE( QueriesTest ) {
     BOOST_CHECK_EQUAL( G4.maximalFactorDomains()[0], v01 );
     BOOST_CHECK_EQUAL( G4.maximalFactorDomains()[1], v02 );
     BOOST_CHECK_EQUAL( G4.maximalFactorDomains()[2], v123 );
+    BOOST_CHECK( G4.DAG() == dag );
+    BOOST_CHECK( G4.checkCountingNumbers() );
+    BOOST_CHECK_EQUAL( G4.nrORs(), 3 );
+    BOOST_CHECK_EQUAL( G4.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G4.OR(0), facs[0] * facs[2] );
+    BOOST_CHECK_EQUAL( G4.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G4.OR(1), facs[3] );
+    BOOST_CHECK_EQUAL( G4.OR(2).c(), 1 );
+    BOOST_CHECK_EQUAL( G4.OR(2), facs[1] * facs[4] * facs[5] );
+    BOOST_CHECK_EQUAL( G4.nrIRs(), 3 );
+    BOOST_CHECK_EQUAL( G4.IR(0).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G4.IR(0), v0 );
+    BOOST_CHECK_EQUAL( G4.IR(1).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G4.IR(1), v1 );
+    BOOST_CHECK_EQUAL( G4.IR(2).c(), -1.0 );
+    BOOST_CHECK_EQUAL( G4.IR(2), v2 );
+    BOOST_CHECK_EQUAL( G4.fac2OR(0), 0 );
+    BOOST_CHECK_EQUAL( G4.fac2OR(1), 2 );
+    BOOST_CHECK_EQUAL( G4.fac2OR(2), 0 );
+    BOOST_CHECK_EQUAL( G4.fac2OR(3), 1 );
+    BOOST_CHECK_EQUAL( G4.fac2OR(4), 2 );
+    BOOST_CHECK_EQUAL( G4.fac2OR(5), 2 );
 }
 
 
@@ -369,20 +631,48 @@ BOOST_AUTO_TEST_CASE( BackupRestoreTest ) {
     FactorGraph fg( facs );
     RegionGraph G( fg, fg.maximalFactorDomains() );
     RegionGraph Gorg( G );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
 
     BOOST_CHECK_THROW( G.setFactor( 0, Factor( v0 ), false ), Exception );
     G.setFactor( 0, Factor( v01, 2.0 ), false );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     BOOST_CHECK_THROW( G.restoreFactor( 0 ), Exception );
     G.setFactor( 0, Factor( v01, 3.0 ), true );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     G.restoreFactor( 0 );
     BOOST_CHECK_EQUAL( G.factor( 0 )[0], 2.0 );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     G.setFactor( 0, Gorg.factor( 0 ), false );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     G.backupFactor( 0 );
     BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
     G.setFactor( 0, Factor( v01, 2.0 ), false );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     BOOST_CHECK_EQUAL( G.factor( 0 )[0], 2.0 );
     G.restoreFactor( 0 );
     BOOST_CHECK_EQUAL( G.factor(0), Gorg.factor(0) );
+    BOOST_CHECK_EQUAL( G.OR(0).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1).c(), 1 );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
 
     std::map<size_t, Factor> fs;
     fs[0] = Factor( v01, 3.0 );
@@ -391,22 +681,32 @@ BOOST_AUTO_TEST_CASE( BackupRestoreTest ) {
     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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     std::set<size_t> fsind;
     fsind.insert( 0 );
     fsind.insert( 2 );
@@ -415,20 +715,28 @@ BOOST_AUTO_TEST_CASE( BackupRestoreTest ) {
     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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
 
     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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
 
     G.backupFactors( v1 );
     fs[1] = Factor( v12, 5.0 );
@@ -436,18 +744,26 @@ BOOST_AUTO_TEST_CASE( BackupRestoreTest ) {
     BOOST_CHECK_EQUAL( G.factor(0), fs[0] );
     BOOST_CHECK_EQUAL( G.factor(1), fs[1] );
     BOOST_CHECK_EQUAL( G.factor(2), fs[2] );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.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) );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     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] );
+    BOOST_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
     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_CHECK_EQUAL( G.OR(0), G.factor(0) * G.factor(2) );
+    BOOST_CHECK_EQUAL( G.OR(1), G.factor(1) );
 }
 
 
@@ -542,10 +858,14 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
                     BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) * delta );
                 else
                     BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
 
             Gcl.restoreFactors();
             for( size_t j = 0; j < 3; j++ )
                 BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         }
 
     // clampVar
@@ -560,10 +880,14 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
                     BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) * delta );
                 else
                     BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
 
             Gcl.restoreFactors();
             for( size_t j = 0; j < 3; j++ )
                 BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         }
     for( size_t i = 0; i < 3; i++ )
         for( size_t x = 0; x < 2; x++ ) {
@@ -575,10 +899,14 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
                     BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) * 0.0 );
                 else
                     BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
 
             Gcl.restoreFactors();
             for( size_t j = 0; j < 3; j++ )
                 BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         }
     std::vector<size_t> both;
     both.push_back( 0 );
@@ -590,7 +918,11 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
             BOOST_CHECK_EQUAL( Gcl.nrFactors(), 3 );
             for( size_t j = 0; j < 3; j++ )
                 BOOST_CHECK_EQUAL( Gcl.factor(j), G.factor(j) );
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
             Gcl.restoreFactors();
+            BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+            BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         }
 
     // clampFactor
@@ -601,7 +933,11 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
         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) );
+        BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         Gcl.restoreFactor( 0 );
+        BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
     }
     for( size_t x = 0; x < 4; x++ ) {
         Gcl.clampFactor( 1, std::vector<size_t>(1,x), true );
@@ -610,7 +946,11 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
         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) );
+        BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         Gcl.restoreFactor( 1 );
+        BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
     }
     for( size_t x = 0; x < 2; x++ ) {
         Gcl.clampFactor( 2, std::vector<size_t>(1,x), true );
@@ -619,7 +959,11 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
         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 );
+        BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
         Gcl.restoreFactors();
+        BOOST_CHECK_EQUAL( Gcl.OR(0), Gcl.factor(0) * Gcl.factor(2) );
+        BOOST_CHECK_EQUAL( Gcl.OR(1), Gcl.factor(1) );
     }
 
     // makeCavity
@@ -628,17 +972,29 @@ BOOST_AUTO_TEST_CASE( OperationsTest ) {
     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) );
+    BOOST_CHECK_EQUAL( Gcav.OR(0), Gcav.factor(0) * Gcav.factor(2) );
+    BOOST_CHECK_EQUAL( Gcav.OR(1), Gcav.factor(1) );
     Gcav.restoreFactors();
+    BOOST_CHECK_EQUAL( Gcav.OR(0), Gcav.factor(0) * Gcav.factor(2) );
+    BOOST_CHECK_EQUAL( Gcav.OR(1), Gcav.factor(1) );
     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 ) );
+    BOOST_CHECK_EQUAL( Gcav.OR(0), Gcav.factor(0) * Gcav.factor(2) );
+    BOOST_CHECK_EQUAL( Gcav.OR(1), Gcav.factor(1) );
     Gcav.restoreFactors();
+    BOOST_CHECK_EQUAL( Gcav.OR(0), Gcav.factor(0) * Gcav.factor(2) );
+    BOOST_CHECK_EQUAL( Gcav.OR(1), Gcav.factor(1) );
     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) );
+    BOOST_CHECK_EQUAL( Gcav.OR(0), Gcav.factor(0) * Gcav.factor(2) );
+    BOOST_CHECK_EQUAL( Gcav.OR(1), Gcav.factor(1) );
     Gcav.restoreFactors();
+    BOOST_CHECK_EQUAL( Gcav.OR(0), Gcav.factor(0) * Gcav.factor(2) );
+    BOOST_CHECK_EQUAL( Gcav.OR(1), Gcav.factor(1) );
 }
 
 
@@ -652,9 +1008,9 @@ BOOST_AUTO_TEST_CASE( IOTest ) {
     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() );
+    facs.push_back( Factor( v01, 1.0 ) );
+    facs.push_back( Factor( v12, 2.0 ) );
+    facs.push_back( Factor( v1, 3.0 ) );
     std::vector<Var> vars;
     vars.push_back( v0 );
     vars.push_back( v1 );
@@ -662,79 +1018,20 @@ BOOST_AUTO_TEST_CASE( IOTest ) {
 
     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 );
-    }
+    BOOST_CHECK_THROW( G.WriteToFile( "regiongraph_test.fg" ), Exception );
+    BOOST_CHECK_THROW( G.ReadFromFile( "regiongraph_test.fg" ), Exception );
+    BOOST_CHECK_THROW( G.printDot( std::cout ), Exception );
 
     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 );
-    }*/
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "digraph RegionGraph {" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "node[shape=box];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ta0 [label=\"a0: {x0, x1}, c=1\"];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ta1 [label=\"a1: {x1, x2}, c=1\"];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "node[shape=ellipse];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\tb0 [label=\"b0: {x1}, c=-1\"];" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ta0 -> b0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "\ta1 -> b0;" );
+    std::getline( ss, s ); BOOST_CHECK_EQUAL( s, "}" );
 }