Implemented various heuristics for choosing a variable elimination sequence in JTree
[libdai.git] / src / hak.cpp
index 309e051..81698dc 100644 (file)
@@ -155,9 +155,12 @@ HAK::HAK(const FactorGraph & fg, const PropertySet &opts) : DAIAlgRG(), _Qa(), _
     vector<VarSet> cl;
     if( props.clusters == Properties::ClustersType::MIN ) {
         cl = fg.Cliques();
+        constructCVM( fg, cl );
     } else if( props.clusters == Properties::ClustersType::DELTA ) {
+        cl.reserve( fg.nrVars() );
         for( size_t i = 0; i < fg.nrVars(); i++ )
-            cl.push_back(fg.Delta(i));
+            cl.push_back( fg.Delta(i) );
+        constructCVM( fg, cl );
     } else if( props.clusters == Properties::ClustersType::LOOP ) {
         cl = fg.Cliques();
         set<VarSet> scl;
@@ -173,11 +176,36 @@ HAK::HAK(const FactorGraph & fg, const PropertySet &opts) : DAIAlgRG(), _Qa(), _
             for( vector<VarSet>::const_iterator cli = cl.begin(); cli != cl.end(); cli++ )
                 cerr << *cli << endl;
         }
+        constructCVM( fg, cl );
+    } else if( props.clusters == Properties::ClustersType::BETHE ) {
+        // build outer regions (the cliques)
+        cl = fg.Cliques();
+        size_t nrEdges = 0;
+        for( size_t c = 0; c < cl.size(); c++ )
+            nrEdges += cl[c].size();
+
+        // build inner regions (single variables)
+        vector<Region> irs;
+        irs.reserve( fg.nrVars() );
+        for( size_t i = 0; i < fg.nrVars(); i++ )
+            irs.push_back( Region( fg.var(i), 1.0 ) );
+
+        // build edges (an outer and inner region are connected if the outer region contains the inner one)
+        // and calculate counting number for inner regions
+        vector<std::pair<size_t, size_t> > edges;
+        edges.reserve( nrEdges );
+        for( size_t c = 0; c < cl.size(); c++ )
+            for( size_t i = 0; i < irs.size(); i++ )
+                if( cl[c] >> irs[i] ) {
+                    edges.push_back( make_pair( c, i ) );
+                    irs[i].c() -= 1.0;
+                }
+
+        // build region graph
+        RegionGraph::construct( fg, cl, irs, edges );
     } else
         DAI_THROW(UNKNOWN_ENUM_VALUE);
 
-    RegionGraph rg(fg,cl);
-    RegionGraph::operator=(rg);
     construct();
 
     if( props.verbose >= 3 )
@@ -264,17 +292,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 +377,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 +432,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 +452,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 +466,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 +543,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++ )