Merged regiongraph.* and daialg.* from SVN head,
authorJoris Mooij <jorism@osun.tuebingen.mpg.de>
Fri, 26 Sep 2008 09:44:29 +0000 (11:44 +0200)
committerJoris Mooij <jorism@osun.tuebingen.mpg.de>
Fri, 26 Sep 2008 09:44:29 +0000 (11:44 +0200)
and did some small fixes here and there

15 files changed:
STATUS
include/dai/bp.h
include/dai/daialg.h
include/dai/exactinf.h
include/dai/factorgraph.h
include/dai/hak.h
include/dai/jtree.h
include/dai/lc.h
include/dai/mf.h
include/dai/mr.h
include/dai/regiongraph.h
include/dai/treeep.h
src/daialg.cpp
src/factorgraph.cpp
src/regiongraph.cpp

diff --git a/STATUS b/STATUS
index db2f534..a532b03 100644 (file)
--- a/STATUS
+++ b/STATUS
@@ -10,6 +10,9 @@ class ExtFactorGraph : public FactorGraph {
        // blabla
 }
 A disadvantage of this approach may be worse cachability.
+A problem is if there are nog properties for nodes (type 1), nodes (type 2)
+or edges. Maybe this can be solved using specializations, or using variadac
+template arguments or something similar?
 - BipartiteGraph::isConnected should be optimized.
 - http://www.boost.org/development/requirements.html#Design_and_Programming
 - Would it be a good idea to cache second-order neighborhoods (delta's) in BipGraph?
@@ -18,6 +21,9 @@ A disadvantage of this approach may be worse cachability.
   No, a better idea is to avoid calls to findVar() as much as possible.
 - Can the FactorGraph constructors be optimized further?
 - Remove x2x?
+- Add iterations (like maxdiff, but counts iterations) to some algorithms
+- A DAIAlg<T> should not inherit from a FactorGraph/RegionGraph, but should store a
+reference to it
 
 TODO FOR RELEASE:
 - Test Visual C++ compatibility
@@ -56,15 +62,15 @@ factor.h
 prob.h
 factorgraph.h
 factorgraph.cpp
+regiongraph.h
+regiongraph.cpp
+daialg.h
+daialg.cpp
 
 FILES IN SVN HEAD THAT ARE STILL RELEVANT:
 ChangeLog
 README
 TODO
-regiongraph.h
-regiongraph.cpp
-daialg.h
-daialg.cpp
 
 bp.h
 bp.cpp
index 48c2957..324e644 100644 (file)
@@ -60,8 +60,10 @@ class BP : public DAIAlgFG {
         BP() : DAIAlgFG(), edges(), props(), maxdiff(0.0) {};
         /// Copy constructor
         BP( const BP & x ) : DAIAlgFG(x), edges(x.edges), props(x.props), maxdiff(x.maxdiff) {};
-        /// Clone *this
+        /// Clone *this (virtual copy constructor)
         BP* clone() const { return new BP(*this); }
+        /// Create (virtual constructor)
+        virtual BP* create() const { return new BP; }
         /// Construct from FactorGraph fg and PropertySet opts
         BP( const FactorGraph & fg, const PropertySet &opts ) : DAIAlgFG(fg), edges(), props(), maxdiff(0.0) {
             setProperties( opts );
@@ -92,6 +94,8 @@ class BP : public DAIAlgFG {
         std::string identify() const;
         void create();
         void init();
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &ns );
         double run();
 
         void findMaxResidual( size_t &i, size_t &_I );
@@ -103,7 +107,6 @@ class BP : public DAIAlgFG {
         std::vector<Factor> beliefs() const;
         Real logZ() const;
 
-        void init( const VarSet &ns );
         void restoreFactors( const VarSet &ns ) { FactorGraph::restoreFactors(ns); init(ns); }
 
         /// Set Props according to the PropertySet opts, where the values can be stored as std::strings or as the type of the corresponding Props member
index 5d3226e..ba3edee 100644 (file)
@@ -41,6 +41,9 @@ class InfAlg {
         /// Clone (virtual copy constructor)
         virtual InfAlg* clone() const = 0;
 
+        /// Create (virtual constructor)
+        virtual InfAlg* create() const = 0;
+        
         /// Virtual desctructor
         // (this is needed because this class contains virtual functions)
         virtual ~InfAlg() {}
@@ -63,6 +66,9 @@ class InfAlg {
         /// Clear messages and beliefs
         virtual void init() = 0;
 
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &ns ) = 0;
+
         /// The actual approximate inference algorithm
         virtual double run() = 0;
 
@@ -77,10 +83,10 @@ class InfAlg {
         virtual void restoreFactors( const VarSet &ns ) = 0;
 
         /// Clamp variable n to value i (i.e. multiply with a Kronecker delta \f$\delta_{x_n, i}\f$)
-        virtual void clamp( const Var & n, size_t i ) = 0;
+        virtual void clamp( const Var & n, size_t i, bool backup = false ) = 0;
 
         /// Set all factors interacting with var(i) to 1
-        virtual void makeCavity( size_t i ) = 0;
+        virtual void makeCavity( size_t i, bool backup = false ) = 0;
 
         /// Get reference to underlying FactorGraph
         virtual FactorGraph &fg() = 0;
@@ -105,6 +111,15 @@ class DAIAlg : public InfAlg, public T {
         /// Copy constructor
         DAIAlg( const DAIAlg & x ) : InfAlg(x), T(x) {}
 
+        /// Assignment operator
+        DAIAlg & operator=( const DAIAlg &x ) {
+            if( this != &x ) {
+                InfAlg::operator=(x);
+                T::operator=(x);
+            }
+            return *this;
+        }
+
         /// Save factor I (using T::backupFactor)
         void backupFactor( size_t I ) { T::backupFactor( I ); }
         /// Save Factors involving ns (using T::backupFactors)
@@ -116,10 +131,10 @@ class DAIAlg : public InfAlg, public T {
         void restoreFactors( const VarSet &ns ) { T::restoreFactors( ns ); }
 
         /// Clamp variable n to value i (i.e. multiply with a Kronecker delta \f$\delta_{x_n, i}\f$) (using T::clamp)
-        void clamp( const Var & n, size_t i ) { T::clamp( n, i ); }
+        void clamp( const Var & n, size_t i, bool backup = false ) { T::clamp( n, i, backup ); }
 
         /// Set all factors interacting with var(i) to 1 (using T::makeCavity)
-        void makeCavity( size_t i ) { T::makeCavity( i ); }
+        void makeCavity( size_t i, bool backup = false ) { T::makeCavity( i, backup ); }
 
         /// Get reference to underlying FactorGraph
         FactorGraph &fg() { return (FactorGraph &)(*this); }
index 204c28f..7dd0bcb 100644 (file)
@@ -73,11 +73,11 @@ class ExactInf : public DAIAlgFG {
             return *this;
         }
 
-/*        /// Create (virtual constructor)
+        /// Create (virtual constructor)
         virtual ExactInf* create() const {
             return new ExactInf();
         }
-*/
+
         /// Return maximum difference between single node 
         /// beliefs for two consecutive iterations
         virtual double maxDiff() const {
@@ -108,7 +108,9 @@ class ExactInf : public DAIAlgFG {
         virtual void init();
 
         /// Clear messages and beliefs corresponding to the nodes in ns
-        virtual void init( const VarSet &/*ns*/ ) {}
+        virtual void init( const VarSet &/*ns*/ ) {
+            DAI_THROW(NOT_IMPLEMENTED);
+        }
 
         /// The actual approximate inference algorithm
         virtual double run();
index 2a99be6..3cc0a9f 100644 (file)
@@ -36,19 +36,19 @@ class FactorGraph {
     public:
         BipartiteGraph         G;
         std::vector<Var>       vars;
-        std::vector<Factor>    factors;
         typedef BipartiteGraph::Neighbor  Neighbor;
         typedef BipartiteGraph::Neighbors Neighbors;
         typedef BipartiteGraph::Edge      Edge;
 
     private:
+        std::vector<Factor>      _factors;
         std::map<size_t,Factor>  _backup;
 
     public:
         /// Default constructor
-        FactorGraph() : G(), vars(), factors(), _backup() {}
+        FactorGraph() : G(), vars(), _factors(), _backup() {}
         /// Copy constructor
-        FactorGraph(const FactorGraph & x) : G(x.G), vars(x.vars), factors(x.factors), _backup(x._backup) {}
+        FactorGraph(const FactorGraph & x) : G(x.G), vars(x.vars), _factors(x._factors), _backup(x._backup) {}
         /// Construct FactorGraph from vector of Factors
         FactorGraph(const std::vector<Factor> &P);
         // Construct a FactorGraph from given factor and variable iterators
@@ -60,8 +60,8 @@ class FactorGraph {
             if( this != &x ) {
                 G          = x.G;
                 vars       = x.vars;
-                factors    = x.factors;
-                _backup = x._backup;
+                _factors   = x._factors;
+                _backup    = x._backup;
             }
             return *this;
         }
@@ -81,14 +81,17 @@ class FactorGraph {
         Var & var(size_t i) { return vars[i]; }
         /// Get const reference to i'th variable
         const Var & var(size_t i) const { return vars[i]; }
-        Factor & factor(size_t I) { return factors[I]; }
         /// Get const reference to I'th factor
-        const Factor & factor(size_t I) const { return factors[I]; }
+        Factor & factor(size_t I) { return _factors[I]; }
+        /// Get const reference to I'th factor
+        const Factor & factor(size_t I) const { return _factors[I]; }
+        /// Get const reference to all factors
+        const std::vector<Factor> & factors() const { return _factors; }
 
         /// Get number of variables
         size_t nrVars() const { return vars.size(); }
         /// Get number of factors
-        size_t nrFactors() const { return factors.size(); }
+        size_t nrFactors() const { return _factors.size(); }
         size_t nrEdges() const { return G.nrEdges(); }
 
         /// Provides read access to neighbors of variable
@@ -149,10 +152,10 @@ class FactorGraph {
 
         /// Set the content of the I'th factor and make a backup of its old content if backup == true
         virtual void setFactor( size_t I, const Factor &newFactor, bool backup = false ) {
-            assert( newFactor.vars() == factors[I].vars() ); 
+            assert( newFactor.vars() == factor(I).vars() ); 
             if( backup )
                 backupFactor( I );
-            factors[I] = newFactor; 
+            _factors[I] = newFactor; 
         }
 
         /// Set the contents of all factors as specified by facs and make a backup of the old contents if backup == true
@@ -213,9 +216,9 @@ template<typename FactorInputIterator, typename VarInputIterator>
 FactorGraph::FactorGraph(FactorInputIterator fact_begin, FactorInputIterator fact_end, VarInputIterator var_begin, VarInputIterator var_end, size_t nr_fact_hint, size_t nr_var_hint ) : G(), _backup() {
     // add factors
     size_t nrEdges = 0;
-    factors.reserve( nr_fact_hint );
+    _factors.reserve( nr_fact_hint );
     for( FactorInputIterator p2 = fact_begin; p2 != fact_end; ++p2 ) {
-        factors.push_back( *p2 );
+        _factors.push_back( *p2 );
         nrEdges += p2->vars().size();
     }
 
index 0c5fd5e..112c881 100644 (file)
@@ -63,6 +63,9 @@ class HAK : public DAIAlgRG {
         /// Clone function
         HAK* clone() const { return new HAK(*this); }
         
+        /// Create (virtual constructor)
+        virtual HAK* create() const { return new HAK(); }
+
         /// Construct from RegionGraph
         HAK(const RegionGraph & rg, const PropertySet &opts);
 
@@ -94,13 +97,14 @@ class HAK : public DAIAlgRG {
         double doDoubleLoop();
         double run();
         void init();
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &ns );
         std::string identify() const;
         Factor belief( const Var &n ) const;
         Factor belief( const VarSet &ns ) const;
         std::vector<Factor> beliefs() const;
         Real logZ () const;
 
-        void init( const VarSet &ns );
         void restoreFactors( const VarSet &ns ) { RegionGraph::restoreFactors( ns ); init( ns ); }
         void setProperties( const PropertySet &opts );
         PropertySet getProperties() const;
index 17d8cf2..20c9dd5 100644 (file)
@@ -57,6 +57,8 @@ class JTree : public DAIAlgRG {
         JTree() : DAIAlgRG(), _RTree(), _Qa(), _Qb(), _mes(), _logZ(), props() {}
         JTree( const JTree& x ) : DAIAlgRG(x), _RTree(x._RTree), _Qa(x._Qa), _Qb(x._Qb), _mes(x._mes), _logZ(x._logZ), props(x.props) {}
         JTree* clone() const { return new JTree(*this); }
+        /// Create (virtual constructor)
+        virtual JTree* create() const { return new JTree(); }
         JTree & operator=( const JTree& x ) {
             if( this != &x ) {
                 DAIAlgRG::operator=(x);
@@ -78,6 +80,8 @@ class JTree : public DAIAlgRG {
         static const char *Name;
         std::string identify() const;
         void init() {}
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &/*ns*/ ) {}
         void runHUGIN();
         void runShaferShenoy();
         double run();
@@ -86,7 +90,6 @@ class JTree : public DAIAlgRG {
         std::vector<Factor> beliefs() const;
         Real logZ() const;
 
-        void init( const VarSet &/*ns*/ ) {}
         void restoreFactors( const VarSet &ns ) { RegionGraph::restoreFactors( ns ); init( ns ); }
 
         size_t findEfficientTree( const VarSet& ns, DEdgeVec &Tree, size_t PreviousRoot=(size_t)-1 ) const;
index 0648ce1..86034f7 100644 (file)
@@ -66,6 +66,8 @@ class LC : public DAIAlgFG {
         LC(const LC & x) : DAIAlgFG(x), _pancakes(x._pancakes), _cavitydists(x._cavitydists), _phis(x._phis), _beliefs(x._beliefs), props(x.props), maxdiff(x.maxdiff) {}
         /// Clone function
         LC* clone() const { return new LC(*this); }
+        /// Create (virtual constructor)
+        virtual LC* create() const { return new LC(); }
         /// Construct LC object from a FactorGraph and parameters
         LC( const FactorGraph & fg, const PropertySet &opts );
         /// Assignment operator
@@ -88,6 +90,10 @@ class LC : public DAIAlgFG {
         long SetCavityDists( std::vector<Factor> &Q );
 
         void init();
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &/*ns*/ ) {
+            DAI_THROW(NOT_IMPLEMENTED);
+        }
         Factor NewPancake (size_t i, size_t _I, bool & hasNaNs);
         double run();
 
index 2b6113d..6318886 100644 (file)
@@ -50,6 +50,8 @@ class MF : public DAIAlgFG {
         // copy constructor
         MF( const MF& x ) : DAIAlgFG(x), _beliefs(x._beliefs), props(x.props), maxdiff(x.maxdiff) {}
         MF* clone() const { return new MF(*this); }
+        /// Create (virtual constructor)
+        virtual MF* create() const { return new MF(); }
         // construct MF object from FactorGraph
         MF( const FactorGraph & fg, const PropertySet &opts ) : DAIAlgFG(fg), _beliefs(), props(), maxdiff(0.0) {
             setProperties( opts );
@@ -70,6 +72,8 @@ class MF : public DAIAlgFG {
         std::string identify() const;
         void create();
         void init();
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &ns );
         double run();
         Factor beliefV (size_t i) const;
         Factor belief (const Var &n) const;
@@ -77,7 +81,6 @@ class MF : public DAIAlgFG {
         std::vector<Factor> beliefs() const;
         Real logZ() const;
 
-        void init( const VarSet &ns );
         void restoreFactors( const VarSet &ns ) { FactorGraph::restoreFactors(ns); init(ns); }
         void setProperties( const PropertySet &opts );
         PropertySet getProperties() const;
index 2c09295..2b55e96 100644 (file)
@@ -68,6 +68,7 @@ class MR : public DAIAlgFG {
         double maxdiff;
 
     public:
+        MR() {}
         MR( const FactorGraph & fg, const PropertySet &opts );
         void init(size_t Nin, double *_w, double *_th);
         void makekindex();
@@ -88,6 +89,10 @@ class MR : public DAIAlgFG {
             return 0.0; 
         }
         void init() {}
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &/*ns*/ ) {
+            DAI_THROW(NOT_IMPLEMENTED);
+        }
         static const char *Name;
         std::string identify() const;
         double _tJ(size_t i, sub_nb A);
@@ -103,6 +108,8 @@ class MR : public DAIAlgFG {
 
         double sign(double a) { return (a >= 0) ? 1.0 : -1.0; }
         MR* clone() const { return new MR(*this); }
+        /// Create (virtual constructor)
+        virtual MR* create() const { return new MR(); }
 
         void setProperties( const PropertySet &opts );
         PropertySet getProperties() const;
index 24cf3d2..d83bf04 100644 (file)
@@ -136,6 +136,32 @@ class RegionGraph : public FactorGraph {
             return *this;
         }
 
+        /// Create (virtual default constructor)
+        virtual RegionGraph* create() const {
+            return new RegionGraph();
+        }
+
+        /// Clone (virtual copy constructor)
+        virtual RegionGraph* clone() const {
+            return new RegionGraph(*this);
+        }
+
+        /// Set the content of the I'th factor and make a backup of its old content if backup == true
+        virtual void setFactor( size_t I, const Factor &newFactor, bool backup = false ) {
+            FactorGraph::setFactor( I, newFactor, backup ); 
+            RecomputeOR( I ); 
+        }
+
+        /// Set the contents of all factors as specified by facs and make a backup of the old contents if backup == true
+        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++ )
+                ns |= fac->second.vars();
+            RecomputeORs( ns ); 
+        }
+
+
         /// Provides read access to outer region
         const FRegion & OR(size_t alpha) const { return ORs[alpha]; }
         /// Provides access to outer region
@@ -177,18 +203,6 @@ class RegionGraph : public FactorGraph {
         /// Recompute all outer regions involving factor I
         void RecomputeOR( size_t I );
 
-        /// We have to overload FactorGraph::clamp because the corresponding outer regions have to be recomputed
-        void clamp( const Var &n, size_t i ) { FactorGraph::clamp( n, i ); RecomputeORs( n ); }
-
-        /// We have to overload FactorGraph::makeCavity because the corresponding outer regions have to be recomputed
-        void makeCavity( size_t i ) { FactorGraph::makeCavity( i ); RecomputeORs( var(i) ); }
-
-        /// We have to overload FactorGraph::restoreFactors because the corresponding outer regions have to be recomputed
-        void restoreFactors( const VarSet &ns ) { FactorGraph::restoreFactors( ns ); RecomputeORs( ns ); }
-
-        /// We have to overload FactorGraph::restoreFactor because the corresponding outer regions have to be recomputed
-        void restoreFactor( size_t I ) { FactorGraph::restoreFactor( I ); RecomputeOR( I ); }
-
         /// Send RegionGraph to output stream
         friend std::ostream & operator << ( std::ostream & os, const RegionGraph & rg );
 };
index 684fae3..9120085 100644 (file)
@@ -104,6 +104,8 @@ class TreeEP : public JTree {
                     _Q[I].I() = &factor(I);
         }
         TreeEP* clone() const { return new TreeEP(*this); }
+        /// Create (virtual constructor)
+        virtual TreeEP* create() const { return new TreeEP(); }
         TreeEP & operator=( const TreeEP& x ) {
             if( this != &x ) {
                 JTree::operator=(x);
@@ -122,12 +124,13 @@ class TreeEP : public JTree {
         static const char *Name;
         std::string identify() const;
         void init();
+        /// Clear messages and beliefs corresponding to the nodes in ns
+        virtual void init( const VarSet &/*ns*/ ) { init(); }
         double run();
         Real logZ() const;
 
         bool offtree( size_t I ) const { return (fac2OR[I] == -1U); }
 
-        void init( const VarSet &/*ns*/ ) { init(); }
         void restoreFactors( const VarSet &ns ) { RegionGraph::restoreFactors( ns ); init( ns ); }
 
         void setProperties( const PropertySet &opts );
index b8adfb9..0fabbfb 100644 (file)
@@ -29,8 +29,6 @@ namespace dai {
 using namespace std;
 
 
-/// Calculate the marginal of obj on ns by clamping 
-/// all variables in ns and calculating logZ for each joined state
 Factor calcMarginal( const InfAlg & obj, const VarSet & ns, bool reInit ) {
     Factor Pns (ns);
     
@@ -50,6 +48,8 @@ Factor calcMarginal( const InfAlg & obj, const VarSet & ns, bool reInit ) {
         // run DAIAlg, calc logZ, store in Pns
         if( reInit )
             clamped->init();
+        else
+            clamped->init(ns);
         clamped->run();
 
         Real Z;
@@ -69,7 +69,7 @@ Factor calcMarginal( const InfAlg & obj, const VarSet & ns, bool reInit ) {
 
     delete clamped;
 
-    return( Pns.normalized(Prob::NORMPROB) );
+    return( Pns.normalized() );
 }
 
 
@@ -86,9 +86,9 @@ vector<Factor> calcPairBeliefs( const InfAlg & obj, const VarSet& ns, bool reIni
     for( size_t j = 0; j < N; j++ )
         for( size_t k = 0; k < N; k++ )
             if( j == k )
-                pairbeliefs.push_back(Factor());
+                pairbeliefs.push_back( Factor() );
             else
-                pairbeliefs.push_back(Factor(vns[j] | vns[k]));
+                pairbeliefs.push_back( Factor( vns[j] | vns[k] ) );
 
     InfAlg *clamped = obj.clone();
     if( !reInit )
@@ -98,12 +98,11 @@ vector<Factor> calcPairBeliefs( const InfAlg & obj, const VarSet& ns, bool reIni
     for( size_t j = 0; j < N; j++ ) {
         // clamp Var j to its possible values
         for( size_t j_val = 0; j_val < vns[j].states(); j_val++ ) {
-            // save unclamped factors connected to ns
-            clamped->backupFactors( ns );
-            
-            clamped->clamp( vns[j], j_val );
+            clamped->clamp( vns[j], j_val, true );
             if( reInit )
                 clamped->init();
+            else
+                clamped->init(ns);
             clamped->run();
 
             //if( j == 0 )
@@ -155,7 +154,7 @@ Factor calcMarginal2ndO( const InfAlg & obj, const VarSet& ns, bool reInit ) {
     for( size_t ij = 0; ij < pairbeliefs.size(); ij++ )
         Pns *= pairbeliefs[ij];
     
-    return( Pns.normalized(Prob::NORMPROB) );
+    return( Pns.normalized() );
 }
 
 
@@ -184,6 +183,8 @@ vector<Factor> calcPairBeliefsNew( const InfAlg & obj, const VarSet& ns, bool re
                     clamped->clamp( *nk, k_val );
                     if( reInit )
                         clamped->init();
+                    else
+                        clamped->init(ns);
                     clamped->run();
 
                     double Z_xj = 1.0;
index 922e9e6..7834467 100644 (file)
@@ -41,10 +41,10 @@ using namespace std;
 FactorGraph::FactorGraph( const std::vector<Factor> &P ) : G(), _backup() {
     // add factors, obtain variables
     set<Var> _vars;
-    factors.reserve( P.size() );
+    _factors.reserve( P.size() );
     size_t nrEdges = 0;
     for( vector<Factor>::const_iterator p2 = P.begin(); p2 != P.end(); p2++ ) {
-        factors.push_back( *p2 );
+        _factors.push_back( *p2 );
         copy( p2->vars().begin(), p2->vars().end(), inserter( _vars, _vars.begin() ) );
         nrEdges += p2->vars().size();
     }
@@ -114,7 +114,7 @@ istream& operator >> (istream& is, FactorGraph& fg) {
     long verbose = 0;
 
     try {
-        vector<Factor> factors;
+        vector<Factor> facs;
         size_t nr_Factors;
         string line;
         
@@ -175,7 +175,7 @@ istream& operator >> (istream& is, FactorGraph& fg) {
                     vardims[labels[mi]] = dims[mi];
                 I_vars |= Var(labels[mi], dims[mi]);
             }
-            factors.push_back( Factor( I_vars, 0.0 ) );
+            facs.push_back( Factor( I_vars, 0.0 ) );
             
             // calculate permutation sigma (internally, members are sorted)
             vector<size_t> sigma(nr_members,0);
@@ -210,14 +210,14 @@ istream& operator >> (istream& is, FactorGraph& fg) {
 
                 // store value, but permute indices first according
                 // to internal representation
-                factors.back()[permindex.convert_linear_index( li  )] = val;
+                facs.back()[permindex.convert_linear_index( li  )] = val;
             }
         }
 
         if( verbose >= 3 )
-            cout << "factors:" << factors << endl;
+            cout << "factors:" << facs << endl;
 
-        fg = FactorGraph(factors);
+        fg = FactorGraph(facs);
     } catch (char *e) {
         cout << e << endl;
     }
index 8078279..de21ac2 100644 (file)
@@ -21,6 +21,7 @@
 
 #include <algorithm>
 #include <cmath>
+#include <boost/dynamic_bitset.hpp>
 #include <dai/regiongraph.h>
 #include <dai/factorgraph.h>
 #include <dai/clustergraph.h>
@@ -46,7 +47,6 @@ RegionGraph::RegionGraph( const FactorGraph &fg, const std::vector<Region> &ors,
         for( alpha = 0; alpha < nrORs(); alpha++ )
             if( OR(alpha).vars() >> factor(I).vars() ) {
                 fac2OR.push_back( alpha );
-//              OR(alpha) *= factor(I);
                 break;
             }
         assert( alpha != nrORs() );
@@ -82,7 +82,6 @@ RegionGraph::RegionGraph( const FactorGraph &fg, const std::vector<VarSet> &cl )
         for( alpha = 0; alpha < nrORs(); alpha++ )
             if( OR(alpha).vars() >> factor(I).vars() ) {
                 fac2OR.push_back( alpha );
-//              OR(alpha) *= factor(I);
                 break;
             }
         assert( alpha != nrORs() );
@@ -93,9 +92,9 @@ RegionGraph::RegionGraph( const FactorGraph &fg, const std::vector<VarSet> &cl )
     set<VarSet> betas;
     for( size_t alpha = 0; alpha < cg.clusters.size(); alpha++ )
         for( size_t alpha2 = alpha; (++alpha2) != cg.clusters.size(); ) {
-            VarSet intersect = cg.clusters[alpha] & cg.clusters[alpha2];
-            if( intersect.size() > 0 )
-                betas.insert( intersect );
+            VarSet intersection = cg.clusters[alpha] & cg.clusters[alpha2];
+            if( intersection.size() > 0 )
+                betas.insert( intersection );
         }
 
     // Create inner regions - subsequent passes
@@ -104,9 +103,9 @@ RegionGraph::RegionGraph( const FactorGraph &fg, const std::vector<VarSet> &cl )
         new_betas.clear();
         for( set<VarSet>::const_iterator gamma = betas.begin(); gamma != betas.end(); gamma++ )
             for( set<VarSet>::const_iterator gamma2 = gamma; (++gamma2) != betas.end(); ) {
-                VarSet intersect = (*gamma) & (*gamma2);
-                if( (intersect.size() > 0) && (betas.count(intersect) == 0) )
-                    new_betas.insert( intersect );
+                VarSet intersection = (*gamma) & (*gamma2);
+                if( (intersection.size() > 0) && (betas.count(intersection) == 0) )
+                    new_betas.insert( intersection );
             }
         betas.insert(new_betas.begin(), new_betas.end());
     } while( new_betas.size() );
@@ -142,7 +141,7 @@ void RegionGraph::Calc_Counting_Numbers() {
     // Calculates counting numbers of inner regions based upon counting numbers of outer regions
     
     vector<vector<size_t> > ancestors(nrIRs());
-    vector<bool> assigned(nrIRs(), false);
+    boost::dynamic_bitset<> assigned(nrIRs());
     for( size_t beta = 0; beta < nrIRs(); beta++ ) {
         IR(beta).c() = 0.0;
         for( size_t beta2 = 0; beta2 < nrIRs(); beta2++ )
@@ -166,7 +165,7 @@ void RegionGraph::Calc_Counting_Numbers() {
                     for( vector<size_t>::const_iterator beta2 = ancestors[beta].begin(); beta2 != ancestors[beta].end(); beta2++ )
                         c -= IR(*beta2).c();
                     IR(beta).c() = c;
-                    assigned[beta] = true;
+                    assigned.set(beta, true);
                     new_counting = true;
                 }
             }
@@ -232,11 +231,16 @@ void RegionGraph::RecomputeOR( size_t I ) {
 ostream & operator << (ostream & os, const RegionGraph & rg) {
     os << "Outer regions" << endl;
     for( size_t alpha = 0; alpha < rg.nrORs(); alpha++ )
-        os << rg.OR(alpha).vars() << ": c = " << rg.OR(alpha).c() << endl;
+        os << alpha << ": " << rg.OR(alpha).vars() << ": c = " << rg.OR(alpha).c() << endl;
 
     os << "Inner regions" << endl;
     for( size_t beta = 0; beta < rg.nrIRs(); beta++ )
-        os << (VarSet)rg.IR(beta) << ": c = " << rg.IR(beta).c() << endl;
+        os << beta << ": " << (VarSet)rg.IR(beta) << ": c = " << rg.IR(beta).c() << endl;
+
+    os << "Edges" << endl;
+    for( size_t alpha = 0; alpha < rg.nrORs(); alpha++ )
+        foreach( const RegionGraph::Neighbor &beta, rg.nbOR(alpha) )
+            os << alpha << "->" << beta << endl;   
 
     return(os);
 }