Strengthened convergence criteria of various algorithms
authorJoris Mooij <joris.mooij@tuebingen.mpg.de>
Tue, 12 Jan 2010 11:12:09 +0000 (12:12 +0100)
committerJoris Mooij <joris.mooij@tuebingen.mpg.de>
Tue, 12 Jan 2010 11:12:09 +0000 (12:12 +0100)
22 files changed:
ChangeLog
include/dai/bp.h
include/dai/cbp.h
include/dai/daialg.h
include/dai/doc.h
include/dai/exactinf.h
include/dai/gibbs.h
include/dai/hak.h
include/dai/jtree.h
include/dai/lc.h
include/dai/mf.h
include/dai/mr.h
src/bp.cpp
src/exactinf.cpp
src/gibbs.cpp
src/hak.cpp
src/jtree.cpp
src/lc.cpp
src/mf.cpp
src/mr.cpp
src/treeep.cpp
tests/testfast.out

index 9d2d51c..80e60ee 100644 (file)
--- a/ChangeLog
+++ b/ChangeLog
@@ -1,3 +1,4 @@
+* Strengthened convergence criteria of various algorithms
 * Implemented FBP::logZ()
 * Fixed bug in HAK for trivial region graphs (with only one outer region
   and no inner regions), reported by Alejandro Lage.
index 3bae74d..c7bd160 100644 (file)
@@ -180,8 +180,8 @@ class BP : public DAIAlgFG {
     //@{
         virtual BP* clone() const { return new BP(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &n ) const;
-        virtual Factor belief( const VarSet &ns ) const;
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet &vs ) const;
         virtual Factor beliefV( size_t i ) const;
         virtual Factor beliefF( size_t I ) const;
         virtual std::vector<Factor> beliefs() const;
index f19c2b8..bf1d04a 100644 (file)
@@ -75,8 +75,8 @@ class CBP : public DAIAlgFG {
     //@{
         virtual CBP* clone() const { return new CBP(*this); }
         virtual std::string identify() const { return std::string(Name) + props.toString(); }
-        virtual Factor belief (const Var &n) const { return _beliefsV[findVar(n)]; }
-        virtual Factor belief (const VarSet &) const { DAI_THROW(NOT_IMPLEMENTED); }
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet & ) const { DAI_THROW(NOT_IMPLEMENTED); }
         virtual Factor beliefV( size_t i ) const { return _beliefsV[i]; }
         virtual Factor beliefF( size_t I ) const { return _beliefsF[I]; }
         virtual std::vector<Factor> beliefs() const { return concat(_beliefsV, _beliefsF); }
index 41da1dc..0242c79 100644 (file)
@@ -79,7 +79,7 @@ class InfAlg {
         /// Returns the (approximate) marginal probability distribution of a variable.
         /** \note Before this method is called, run() should have been called.
          */
-        virtual Factor belief( const Var &v ) const = 0;
+        virtual Factor belief( const Var &v ) const { return belief( VarSet(v) ); }
 
         /// Returns the (approximate) marginal probability distribution of a set of variables.
         /** \note Before this method is called, run() should have been called.
index d7b6c05..57b115d 100644 (file)
@@ -55,7 +55,7 @@
 /** \mainpage Reference manual for libDAI - A free/open source C++ library for Discrete Approximate Inference methods
  *  \author Joris Mooij
  *  \version git HEAD
- *  \date November 16, 2009 - or later
+ *  \date January 12, 2010 - or later
  *
  *  <hr size="1">
  *  \section about About libDAI
index 5a1f0bf..58e51dd 100644 (file)
@@ -70,8 +70,8 @@ class ExactInf : public DAIAlgFG {
     //@{
         virtual ExactInf* clone() const { return new ExactInf(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &n ) const { return beliefV( findVar( n ) ); }
-        virtual Factor belief( const VarSet &ns ) const;
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet &vs ) const;
         virtual Factor beliefV( size_t i ) const { return _beliefsV[i]; }
         virtual Factor beliefF( size_t I ) const { return _beliefsF[I]; }
         virtual std::vector<Factor> beliefs() const;
index fbc5f62..ceb4f63 100644 (file)
@@ -75,8 +75,8 @@ class Gibbs : public DAIAlgFG {
     //@{
         virtual Gibbs* clone() const { return new Gibbs(*this); }
         virtual std::string identify() const { return std::string(Name) + printProperties(); }
-        virtual Factor belief( const Var &n ) const;
-        virtual Factor belief( const VarSet &ns ) const;
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet &vs ) const;
         virtual Factor beliefV( size_t i ) const;
         virtual Factor beliefF( size_t I ) const;
         virtual std::vector<Factor> beliefs() const;
index fee0386..9069cdf 100644 (file)
@@ -107,7 +107,6 @@ class HAK : public DAIAlgRG {
     //@{
         virtual HAK* clone() const { return new HAK(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &v ) const;
         virtual Factor belief( const VarSet &vs ) const;
         virtual std::vector<Factor> beliefs() const;
         virtual Real logZ() const;
index 35bbab0..3bed8b1 100644 (file)
@@ -108,8 +108,7 @@ class JTree : public DAIAlgRG {
     //@{
         virtual JTree* clone() const { return new JTree(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &n ) const;
-        virtual Factor belief( const VarSet &ns ) const;
+        virtual Factor belief( const VarSet &vs ) const;
         virtual std::vector<Factor> beliefs() const;
         virtual Real logZ() const;
         virtual void init() {}
index ec8842d..9c97d0b 100644 (file)
@@ -108,8 +108,8 @@ class LC : public DAIAlgFG {
     //@{
         virtual LC* clone() const { return new LC(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &n ) const { return( _beliefs[findVar(n)] ); }
-        virtual Factor belief( const VarSet &/*ns*/ ) const { DAI_THROW(NOT_IMPLEMENTED); return Factor(); }
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet &/*vs*/ ) const { DAI_THROW(NOT_IMPLEMENTED); return Factor(); }
         virtual Factor beliefV( size_t i ) const { return _beliefs[i]; }
         virtual std::vector<Factor> beliefs() const { return _beliefs; }
         virtual Real logZ() const { DAI_THROW(NOT_IMPLEMENTED); return 0.0; }
index 9218579..02208df 100644 (file)
@@ -80,8 +80,8 @@ class MF : public DAIAlgFG {
     //@{
         virtual MF* clone() const { return new MF(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &n ) const;
-        virtual Factor belief( const VarSet &ns ) const;
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet &vs ) const;
         virtual Factor beliefV( size_t i ) const;
         virtual std::vector<Factor> beliefs() const;
         virtual Real logZ() const;
index 3c25b8f..64482d4 100644 (file)
@@ -115,8 +115,9 @@ class MR : public DAIAlgFG {
     //@{
         virtual MR* clone() const { return new MR(*this); }
         virtual std::string identify() const;
-        virtual Factor belief( const Var &n ) const;
-        virtual Factor belief( const VarSet &/*ns*/ ) const { DAI_THROW(NOT_IMPLEMENTED); return Factor(); }
+        virtual Factor belief( const Var &v ) const { return beliefV( findVar( v ) ); }
+        virtual Factor belief( const VarSet &/*vs*/ ) const { DAI_THROW(NOT_IMPLEMENTED); return Factor(); }
+        virtual Factor beliefV( size_t i ) const;
         virtual std::vector<Factor> beliefs() const;
         virtual Real logZ() const { DAI_THROW(NOT_IMPLEMENTED); return 0.0; }
         virtual void init() {}
index fda4c39..2d6fd12 100644 (file)
@@ -232,13 +232,15 @@ Real BP::run() {
         cerr << endl;
 
     double tic = toc();
-    vector<Real> diffs( nrVars(), INFINITY );
     Real maxDiff = INFINITY;
 
-    vector<Factor> old_beliefs;
-    old_beliefs.reserve( nrVars() );
+    vector<Factor> oldBeliefsV, oldBeliefsF;
+    oldBeliefsV.reserve( nrVars() );
     for( size_t i = 0; i < nrVars(); ++i )
-        old_beliefs.push_back( beliefV(i) );
+        oldBeliefsV.push_back( beliefV(i) );
+    oldBeliefsF.reserve( nrFactors() );
+    for( size_t I = 0; I < nrFactors(); ++I )
+        oldBeliefsF.push_back( beliefF(I) );
 
     size_t nredges = nrEdges();
     vector<Edge> update_seq;
@@ -299,12 +301,17 @@ Real BP::run() {
         }
 
         // calculate new beliefs and compare with old ones
+        maxDiff = -INFINITY;
         for( size_t i = 0; i < nrVars(); ++i ) {
-            Factor nb( beliefV(i) );
-            diffs[i] = dist( nb, old_beliefs[i], Prob::DISTLINF );
-            old_beliefs[i] = nb;
+            Factor b( beliefV(i) );
+            maxDiff = std::max( maxDiff, dist( b, oldBeliefsV[i], Prob::DISTLINF ) );
+            oldBeliefsV[i] = b;
+        }
+        for( size_t I = 0; I < nrFactors(); ++I ) {
+            Factor b( beliefF(I) );
+            maxDiff = std::max( maxDiff, dist( b, oldBeliefsF[I], Prob::DISTLINF ) );
+            oldBeliefsF[I] = b;
         }
-        maxDiff = max( diffs );
 
         if( props.verbose >= 3 )
             cerr << Name << "::run:  maxdiff " << maxDiff << " after " << _iters+1 << " passes" << endl;
@@ -411,11 +418,6 @@ Factor BP::beliefF( size_t I ) const {
 }
 
 
-Factor BP::belief( const Var &n ) const {
-    return( beliefV( findVar( n ) ) );
-}
-
-
 vector<Factor> BP::beliefs() const {
     vector<Factor> result;
     for( size_t i = 0; i < nrVars(); ++i )
@@ -427,8 +429,10 @@ vector<Factor> BP::beliefs() const {
 
 
 Factor BP::belief( const VarSet &ns ) const {
-    if( ns.size() == 1 )
-        return belief( *(ns.begin()) );
+    if( ns.size() == 0 )
+        return Factor();
+    else if( ns.size() == 1 )
+        return beliefV( findVar( *(ns.begin() ) ) );
     else {
         size_t I;
         for( I = 0; I < nrFactors(); I++ )
index ffa2e92..8c7d1c7 100644 (file)
@@ -108,7 +108,7 @@ Factor ExactInf::belief( const VarSet &ns ) const {
     if( ns.size() == 0 )
         return Factor();
     else if( ns.size() == 1 ) {
-        return belief( *(ns.begin()) );
+        return beliefV( findVar( *(ns.begin()) ) );
     } else {
         size_t I;
         for( I = 0; I < nrFactors(); I++ )
index 854fbb3..a8c8645 100644 (file)
@@ -205,11 +205,6 @@ Factor Gibbs::beliefF( size_t I ) const {
 }
 
 
-Factor Gibbs::belief( const Var &n ) const {
-    return( beliefV( findVar( n ) ) );
-}
-
-
 vector<Factor> Gibbs::beliefs() const {
     vector<Factor> result;
     for( size_t i = 0; i < nrVars(); ++i )
@@ -221,8 +216,10 @@ vector<Factor> Gibbs::beliefs() const {
 
 
 Factor Gibbs::belief( const VarSet &ns ) const {
-    if( ns.size() == 1 )
-        return belief( *(ns.begin()) );
+    if( ns.size() == 0 )
+        return Factor();
+    else if( ns.size() == 1 )
+        return beliefV( findVar( *(ns.begin()) ) );
     else {
         size_t I;
         for( I = 0; I < nrFactors(); I++ )
index 309e051..6ac23ce 100644 (file)
@@ -264,17 +264,18 @@ Real HAK::doGBP() {
         DAI_ASSERT( nbIR(beta).size() + IR(beta).c() != 0.0 );
 
     // Keep old beliefs to check convergence
-    vector<Factor> old_beliefs;
-    old_beliefs.reserve( nrVars() );
+    vector<Factor> oldBeliefsV;
+    oldBeliefsV.reserve( nrVars() );
     for( size_t i = 0; i < nrVars(); i++ )
-        old_beliefs.push_back( belief( var(i) ) );
-
-    // Differences in single node beliefs
-    vector<Real> diffs( nrVars(), INFINITY );
-    Real maxDiff = INFINITY;
+        oldBeliefsV.push_back( beliefV(i) );
+    vector<Factor> oldBeliefsF;
+    oldBeliefsF.reserve( nrFactors() );
+    for( size_t I = 0; I < nrFactors(); I++ )
+        oldBeliefsF.push_back( beliefF(I) );
 
     // do several passes over the network until maximum number of iterations has
     // been reached or until the maximum belief difference is smaller than tolerance
+    Real maxDiff = INFINITY;
     for( _iters = 0; _iters < props.maxiter && maxDiff > props.tol; _iters++ ) {
         for( size_t beta = 0; beta < nrIRs(); beta++ ) {
             foreach( const Neighbor &alpha, nbIR(beta) ) {
@@ -348,12 +349,17 @@ Real HAK::doGBP() {
         }
 
         // Calculate new single variable beliefs and compare with old ones
-        for( size_t i = 0; i < nrVars(); i++ ) {
-            Factor new_belief = belief( var( i ) );
-            diffs[i] = dist( new_belief, old_beliefs[i], Prob::DISTLINF );
-            old_beliefs[i] = new_belief;
+        maxDiff = -INFINITY;
+        for( size_t i = 0; i < nrVars(); ++i ) {
+            Factor b = beliefV(i);
+            maxDiff = std::max( maxDiff, dist( b, oldBeliefsV[i], Prob::DISTLINF ) );
+            oldBeliefsV[i] = b;
+        }
+        for( size_t I = 0; I < nrFactors(); ++I ) {
+            Factor b = beliefF(I);
+            maxDiff = std::max( maxDiff, dist( b, oldBeliefsF[I], Prob::DISTLINF ) );
+            oldBeliefsF[I] = b;
         }
-        maxDiff = max( diffs );
 
         if( props.verbose >= 3 )
             cerr << Name << "::doGBP:  maxdiff " << maxDiff << " after " << _iters+1 << " passes" << endl;
@@ -398,14 +404,14 @@ Real HAK::doDoubleLoop() {
     }
 
     // Keep old beliefs to check convergence
-    vector<Factor> old_beliefs;
-    old_beliefs.reserve( nrVars() );
+    vector<Factor> oldBeliefsV;
+    oldBeliefsV.reserve( nrVars() );
     for( size_t i = 0; i < nrVars(); i++ )
-        old_beliefs.push_back( belief( var(i) ) );
-
-    // Differences in single node beliefs
-    vector<Real> diffs( nrVars(), INFINITY );
-    Real maxDiff = INFINITY;
+        oldBeliefsV.push_back( beliefV(i) );
+    vector<Factor> oldBeliefsF;
+    oldBeliefsF.reserve( nrFactors() );
+    for( size_t I = 0; I < nrFactors(); I++ )
+        oldBeliefsF.push_back( beliefF(I) );
 
     size_t outer_maxiter   = props.maxiter;
     Real   outer_tol       = props.tol;
@@ -418,6 +424,7 @@ Real HAK::doDoubleLoop() {
 
     size_t outer_iter = 0;
     size_t total_iter = 0;
+    Real maxDiff = INFINITY;
     for( outer_iter = 0; outer_iter < outer_maxiter && maxDiff > outer_tol; outer_iter++ ) {
         // Calculate new outer regions
         for( size_t alpha = 0; alpha < nrORs(); alpha++ ) {
@@ -431,12 +438,17 @@ Real HAK::doDoubleLoop() {
             return 1.0;
 
         // Calculate new single variable beliefs and compare with old ones
+        maxDiff = -INFINITY;
         for( size_t i = 0; i < nrVars(); ++i ) {
-            Factor new_belief = belief( var( i ) );
-            diffs[i] = dist( new_belief, old_beliefs[i], Prob::DISTLINF );
-            old_beliefs[i] = new_belief;
+            Factor b = beliefV(i);
+            maxDiff = std::max( maxDiff, dist( b, oldBeliefsV[i], Prob::DISTLINF ) );
+            oldBeliefsV[i] = b;
+        }
+        for( size_t I = 0; I < nrFactors(); ++I ) {
+            Factor b = beliefF(I);
+            maxDiff = std::max( maxDiff, dist( b, oldBeliefsF[I], Prob::DISTLINF ) );
+            oldBeliefsF[I] = b;
         }
-        maxDiff = max( diffs );
 
         total_iter += Iterations();
 
@@ -503,11 +515,6 @@ Factor HAK::belief( const VarSet &ns ) const {
 }
 
 
-Factor HAK::belief( const Var &n ) const {
-    return belief( (VarSet)n );
-}
-
-
 vector<Factor> HAK::beliefs() const {
     vector<Factor> result;
     for( size_t beta = 0; beta < nrIRs(); beta++ )
index 1c05b5e..4874f6f 100644 (file)
@@ -203,11 +203,6 @@ vector<Factor> JTree::beliefs() const {
 }
 
 
-Factor JTree::belief( const Var &v ) const {
-    return belief( (VarSet)v );
-}
-
-
 void JTree::runHUGIN() {
     for( size_t alpha = 0; alpha < nrORs(); alpha++ )
         Qa[alpha] = OR(alpha);
index c192c28..47c5da4 100644 (file)
@@ -246,8 +246,6 @@ Real LC::run() {
         cerr << endl;
 
     double tic = toc();
-    vector<Real> diffs( nrVars(), INFINITY );
-    Real maxDiff = INFINITY;
 
     Real md = InitCavityDists( props.cavainame, props.cavaiopts );
     if( md > _maxdiff )
@@ -267,9 +265,9 @@ Real LC::run() {
         CalcBelief(i);
     }
 
-    vector<Factor> old_beliefs;
-    for(size_t i=0; i < nrVars(); i++ )
-        old_beliefs.push_back(beliefV(i));
+    vector<Factor> oldBeliefsV;
+    for( size_t i = 0; i < nrVars(); i++ )
+        oldBeliefsV.push_back( beliefV(i) );
 
     bool hasNaNs = false;
     for( size_t i=0; i < nrVars(); i++ )
@@ -291,7 +289,8 @@ Real LC::run() {
 
     // do several passes over the network until maximum number of iterations has
     // been reached or until the maximum belief difference is smaller than tolerance
-    for( _iters=0; _iters < props.maxiter && maxDiff > props.tol; _iters++ ) {
+    Real maxDiff = INFINITY;
+    for( _iters = 0; _iters < props.maxiter && maxDiff > props.tol; _iters++ ) {
         // Sequential updates
         if( props.updates == Properties::UpdateType::SEQRND )
             random_shuffle( update_seq.begin(), update_seq.end() );
@@ -306,11 +305,11 @@ Real LC::run() {
         }
 
         // compare new beliefs with old ones
-        for(size_t i=0; i < nrVars(); i++ ) {
-            diffs[i] = dist( beliefV(i), old_beliefs[i], Prob::DISTLINF );
-            old_beliefs[i] = beliefV(i);
+        maxDiff = -INFINITY;
+        for( size_t i = 0; i < nrVars(); i++ ) {
+            maxDiff = std::max( maxDiff, dist( beliefV(i), oldBeliefsV[i], Prob::DISTLINF ) );
+            oldBeliefsV[i] = beliefV(i);
         }
-        maxDiff = max( diffs );
 
         if( props.verbose >= 3 )
             cerr << Name << "::run:  maxdiff " << maxDiff << " after " << _iters+1 << " passes" << endl;
index 213da2b..b7119d2 100644 (file)
@@ -107,8 +107,6 @@ Real MF::run() {
         cerr << "Starting " << identify() << "...";
 
     double tic = toc();
-    vector<Real> diffs( nrVars(), INFINITY );
-    Real maxDiff = INFINITY;
 
     vector<size_t> update_seq;
     update_seq.reserve( nrVars() );
@@ -117,9 +115,11 @@ Real MF::run() {
 
     // do several passes over the network until maximum number of iterations has
     // been reached or until the maximum belief difference is smaller than tolerance
+    Real maxDiff = INFINITY;
     for( _iters = 0; _iters < props.maxiter && maxDiff > props.tol; _iters++ ) {
         random_shuffle( update_seq.begin(), update_seq.end() );
 
+        maxDiff = -INFINITY;
         foreach( const size_t &i, update_seq ) {
             Factor nb = calcNewBelief( i );
 
@@ -130,11 +130,10 @@ Real MF::run() {
 
             if( props.damping != 0.0 )
                 nb = (nb^(1.0 - props.damping)) * (_beliefs[i]^props.damping);
-            diffs[i] = dist( nb, _beliefs[i], Prob::DISTLINF );
 
+            maxDiff = std::max( maxDiff, dist( nb, _beliefs[i], Prob::DISTLINF ) );
             _beliefs[i] = nb;
         }
-        maxDiff = max( diffs );
 
         if( props.verbose >= 3 )
             cerr << Name << "::run:  maxdiff " << maxDiff << " after " << _iters+1 << " passes" << endl;
@@ -160,16 +159,15 @@ Real MF::run() {
 
 
 Factor MF::beliefV( size_t i ) const {
-    Factor piet;
-    piet = _beliefs[i];
-    piet.normalize();
-    return(piet);
+    return _beliefs[i].normalized();
 }
 
 
 Factor MF::belief (const VarSet &ns) const {
-    if( ns.size() == 1 )
-        return belief( *(ns.begin()) );
+    if( ns.size() == 0 )
+        return Factor();
+    else if( ns.size() == 1 )
+        return beliefV( findVar( *(ns.begin()) ) );
     else {
         DAI_THROW(BELIEF_NOT_AVAILABLE);
         return Factor();
@@ -177,11 +175,6 @@ Factor MF::belief (const VarSet &ns) const {
 }
 
 
-Factor MF::belief (const Var &n) const {
-    return( beliefV( findVar( n ) ) );
-}
-
-
 vector<Factor> MF::beliefs() const {
     vector<Factor> result;
     for( size_t i = 0; i < nrVars(); i++ )
@@ -193,9 +186,9 @@ vector<Factor> MF::beliefs() const {
 Real MF::logZ() const {
     Real s = 0.0;
 
-    for(size_t i=0; i < nrVars(); i++ )
+    for( size_t i = 0; i < nrVars(); i++ )
         s -= beliefV(i).entropy();
-    for(size_t I=0; I < nrFactors(); I++ ) {
+    for( size_t I = 0; I < nrFactors(); I++ ) {
         Factor henk;
         foreach( const Neighbor &j, nbF(I) )  // for all j in I
             henk *= _beliefs[j];
index d57d030..1a74352 100644 (file)
@@ -547,15 +547,13 @@ void MR::makekindex() {
 }
 
 
-Factor MR::belief( const Var &n ) const {
+Factor MR::beliefV( size_t i ) const {
     if( supported ) {
-        size_t i = findVar( n );
-
         Prob x(2);
         x[0] = 0.5 - Mag[i] / 2.0;
         x[1] = 0.5 + Mag[i] / 2.0;
 
-        return Factor( n, x );
+        return Factor( var(i), x );
     } else
         return Factor();
 }
index db6e9b8..1685219 100644 (file)
@@ -369,21 +369,20 @@ void TreeEP::init() {
 Real TreeEP::run() {
     if( props.verbose >= 1 )
         cerr << "Starting " << identify() << "...";
-    if( props.verbose >= 3)
+    if( props.verbose >= 3 )
         cerr << endl;
 
     double tic = toc();
-    vector<Real> diffs( nrVars(), INFINITY );
-    Real maxDiff = INFINITY;
 
-    vector<Factor> old_beliefs;
-    old_beliefs.reserve( nrVars() );
+    vector<Factor> oldBeliefsV;
+    oldBeliefsV.reserve( nrVars() );
     for( size_t i = 0; i < nrVars(); i++ )
-        old_beliefs.push_back(belief(var(i)));
+        oldBeliefsV.push_back( beliefV(i) );
 
     // do several passes over the network until maximum number of iterations has
     // been reached or until the maximum belief difference is smaller than tolerance
-    for( _iters=0; _iters < props.maxiter && maxDiff > props.tol; _iters++ ) {
+    Real maxDiff = INFINITY;
+    for( _iters = 0; _iters < props.maxiter && maxDiff > props.tol; _iters++ ) {
         for( size_t I = 0; I < nrFactors(); I++ )
             if( offtree(I) ) {
                 _Q[I].InvertAndMultiply( Qa, Qb );
@@ -392,12 +391,12 @@ Real TreeEP::run() {
             }
 
         // calculate new beliefs and compare with old ones
+        maxDiff = -INFINITY;
         for( size_t i = 0; i < nrVars(); i++ ) {
-            Factor nb( belief(var(i)) );
-            diffs[i] = dist( nb, old_beliefs[i], Prob::DISTLINF );
-            old_beliefs[i] = nb;
+            Factor nb( beliefV(i) );
+            maxDiff = std::max( maxDiff, dist( nb, oldBeliefsV[i], Prob::DISTLINF ) );
+            oldBeliefsV[i] = nb;
         }
-        maxDiff = max( diffs );
 
         if( props.verbose >= 3 )
             cerr << Name << "::run:  maxdiff " << maxDiff << " after " << _iters+1 << " passes" << endl;
index e21329b..30c081f 100644 (file)
@@ -289,11 +289,11 @@ GBP_LOOP3                                 9.483e-02       3.078e-02       +1.737e-02      1.000e-09
 # ({x13}, (5.266e-01, 4.734e-01))
 # ({x14}, (6.033e-01, 3.967e-01))
 # ({x15}, (1.558e-01, 8.442e-01))
-GBP_LOOP4                                      7.893e-01       3.569e-01       -4.651e-01      1.000e-09       
-# ({x0}, (4.731e-106, 1.000e+00))
+GBP_LOOP4                                      7.893e-01       3.728e-01       -5.824e-01      1.000e-09       
+# ({x0}, (0.000e+00, 1.000e+00))
 # ({x1}, (1.000e+00, 0.000e+00))
 # ({x2}, (0.000e+00, 1.000e+00))
-# ({x3}, (1.000e+00, 2.962e-28))
+# ({x3}, (0.000e+00, 1.000e+00))
 # ({x4}, (1.000e+00, 0.000e+00))
 # ({x5}, (1.000e+00, 0.000e+00))
 # ({x6}, (1.000e+00, 0.000e+00))
@@ -302,15 +302,15 @@ GBP_LOOP4                                 7.893e-01       3.569e-01       -4.651e-01      1.000e-09
 # ({x9}, (0.000e+00, 1.000e+00))
 # ({x10}, (1.000e+00, 0.000e+00))
 # ({x11}, (0.000e+00, 1.000e+00))
-# ({x12}, (6.472e-20, 1.000e+00))
+# ({x12}, (1.000e+00, 0.000e+00))
 # ({x13}, (1.000e+00, 0.000e+00))
 # ({x14}, (0.000e+00, 1.000e+00))
 # ({x15}, (0.000e+00, 1.000e+00))
-GBP_LOOP5                                      7.893e-01       3.569e-01       -4.651e-01      1.000e-09       
-# ({x0}, (4.731e-106, 1.000e+00))
+GBP_LOOP5                                      7.893e-01       3.728e-01       -5.824e-01      1.000e-09       
+# ({x0}, (0.000e+00, 1.000e+00))
 # ({x1}, (1.000e+00, 0.000e+00))
 # ({x2}, (0.000e+00, 1.000e+00))
-# ({x3}, (1.000e+00, 2.962e-28))
+# ({x3}, (0.000e+00, 1.000e+00))
 # ({x4}, (1.000e+00, 0.000e+00))
 # ({x5}, (1.000e+00, 0.000e+00))
 # ({x6}, (1.000e+00, 0.000e+00))
@@ -319,7 +319,7 @@ GBP_LOOP5                                   7.893e-01       3.569e-01       -4.651e-01      1.000e-09
 # ({x9}, (0.000e+00, 1.000e+00))
 # ({x10}, (1.000e+00, 0.000e+00))
 # ({x11}, (0.000e+00, 1.000e+00))
-# ({x12}, (6.472e-20, 1.000e+00))
+# ({x12}, (1.000e+00, 0.000e+00))
 # ({x13}, (1.000e+00, 0.000e+00))
 # ({x14}, (0.000e+00, 1.000e+00))
 # ({x15}, (0.000e+00, 1.000e+00))