Merged SVN head ...
[libdai.git] / src / lc.cpp
index d8da46c..2c796fe 100644 (file)
@@ -27,7 +27,6 @@
 #include <dai/diffs.h>
 #include <dai/util.h>
 #include <dai/alldai.h>
-#include <dai/x2x.h>
 
 
 namespace dai {
@@ -57,6 +56,10 @@ void LC::setProperties( const PropertySet &opts ) {
         props.cavaiopts = opts.getStringAs<PropertySet>("cavaiopts");
     if( opts.hasKey("reinit") )
         props.reinit = opts.getStringAs<bool>("reinit");
+    if( opts.hasKey("damping") )
+        props.damping = opts.getStringAs<double>("damping");
+    else
+        props.damping = 0.0;
 }
 
 
@@ -70,19 +73,36 @@ PropertySet LC::getProperties() const {
     opts.Set( "cavainame", props.cavainame );
     opts.Set( "cavaiopts", props.cavaiopts );
     opts.Set( "reinit", props.reinit );
+    opts.Set( "damping", props.damping );
     return opts;
 }
 
 
-LC::LC( const FactorGraph & fg, const PropertySet &opts ) : DAIAlgFG(fg), _pancakes(), _cavitydists(), _phis(), _beliefs(), props(), maxdiff(0.0) {
+string LC::printProperties() const {
+    stringstream s( stringstream::out );
+    s << "[";
+    s << "tol=" << props.tol << ",";
+    s << "maxiter=" << props.maxiter << ",";
+    s << "verbose=" << props.verbose << ",";
+    s << "cavity=" << props.cavity << ",";
+    s << "updates=" << props.updates << ",";
+    s << "cavainame=" << props.cavainame << ",";
+    s << "cavaiopts=" << props.cavaiopts << ",";
+    s << "reinit=" << props.reinit << ",";
+    s << "damping=" << props.damping << "]";
+    return s.str();
+}
+
+
+LC::LC( const FactorGraph & fg, const PropertySet &opts ) : DAIAlgFG(fg), _pancakes(), _cavitydists(), _phis(), _beliefs(), _maxdiff(0.0), _iters(0), props() {
     setProperties( opts );
 
     // create pancakes
-    _pancakes.resize(nrVars());
+    _pancakes.resize( nrVars() );
    
     // create cavitydists
     for( size_t i=0; i < nrVars(); i++ )
-        _cavitydists.push_back(Factor(delta(i)));
+        _cavitydists.push_back(Factor( delta(i) ));
 
     // create phis
     _phis.reserve( nrVars() );
@@ -94,15 +114,14 @@ LC::LC( const FactorGraph & fg, const PropertySet &opts ) : DAIAlgFG(fg), _panca
     }
 
     // create beliefs
+    _beliefs.reserve( nrVars() );
     for( size_t i=0; i < nrVars(); i++ )
         _beliefs.push_back(Factor(var(i)));
 }
 
 
 string LC::identify() const { 
-    stringstream result (stringstream::out);
-    result << Name << getProperties();
-    return result.str();
+    return string(Name) + printProperties();
 }
 
 
@@ -132,34 +151,11 @@ double LC::CalcCavityDist (size_t i, const std::string &name, const PropertySet
             vector<Factor> pairbeliefs = calcPairBeliefsNew( *cav, cav->fg().delta(i), props.reinit );
             for( size_t ij = 0; ij < pairbeliefs.size(); ij++ )
                 Bi *= pairbeliefs[ij];
-        } else if( props.cavity == Properties::CavityType::PAIRINT ) {
-            Bi = calcMarginal( *cav, cav->fg().delta(i), props.reinit );
-            
-            // Set interactions of order > 2 to zero
-            size_t N = delta(i).size();
-            Real *p = &(*Bi.p().p().begin());
-            x2x::p2logp (N, p);
-            x2x::logp2w (N, p);
-            x2x::fill (N, p, 2, 0.0);
-            x2x::w2logp (N, p);
-//            x2x::logpnorm (N, p);
-            x2x::logp2p (N, p);
-        } else if( props.cavity == Properties::CavityType::PAIRCUM ) {
-            Bi = calcMarginal( *cav, cav->fg().delta(i), props.reinit );
-            
-            // Set cumulants of order > 2 to zero
-            size_t N = delta(i).size();
-            Real *p = &(*Bi.p().p().begin());
-            x2x::p2m (N, p);
-            x2x::m2c (N, p, N);
-            x2x::fill (N, p, 2, 0.0);
-            x2x::c2m (N, p, N);
-            x2x::m2p (N, p);
         }
         maxdiff = cav->maxDiff();
         delete cav;
     }
-    Bi.normalize( Prob::NORMPROB );
+    Bi.normalize();
     _cavitydists[i] = Bi;
 
     return maxdiff;
@@ -170,7 +166,7 @@ double LC::InitCavityDists( const std::string &name, const PropertySet &opts ) {
     double tic = toc();
 
     if( props.verbose >= 1 ) {
-        cout << "LC::InitCavityDists:  ";
+        cout << Name << "::InitCavityDists:  ";
         if( props.cavity == Properties::CavityType::UNIFORM )
             cout << "Using uniform initial cavity distributions" << endl;
         else if( props.cavity == Properties::CavityType::FULL )
@@ -190,7 +186,7 @@ double LC::InitCavityDists( const std::string &name, const PropertySet &opts ) {
     init();
 
     if( props.verbose >= 1 ) {
-        cout << "used " << toc() - tic << " clocks." << endl;
+        cout << Name << "::InitCavityDists used " << toc() - tic << " seconds." << endl;
     }
 
     return maxdiff;
@@ -199,7 +195,7 @@ double LC::InitCavityDists( const std::string &name, const PropertySet &opts ) {
 
 long LC::SetCavityDists( std::vector<Factor> &Q ) {
     if( props.verbose >= 1 ) 
-        cout << "LC::SetCavityDists:  Setting initial cavity distributions" << endl;
+        cout << Name << "::SetCavityDists:  Setting initial cavity distributions" << endl;
     if( Q.size() != nrVars() )
         return -1;
     for( size_t i = 0; i < nrVars(); i++ ) {
@@ -229,7 +225,7 @@ void LC::init() {
               _pancakes[i] *= _phis[i][I.iter];
         }
         
-        _pancakes[i].normalize( Prob::NORMPROB );
+        _pancakes[i].normalize();
 
         CalcBelief(i);
     }
@@ -245,19 +241,21 @@ Factor LC::NewPancake (size_t i, size_t _I, bool & hasNaNs) {
     Factor A_I;
     for( VarSet::const_iterator k = Ivars.begin(); k != Ivars.end(); k++ )
         if( var(i) != *k )
-            A_I *= (_pancakes[findVar(*k)] * factor(I).inverse()).part_sum( Ivars / var(i) );
+            A_I *= (_pancakes[findVar(*k)] * factor(I).inverse()).partSum( Ivars / var(i) );
     if( Ivars.size() > 1 )
         A_I ^= (1.0 / (Ivars.size() - 1));
-    Factor A_Ii = (_pancakes[i] * factor(I).inverse() * _phis[i][_I].inverse()).part_sum( Ivars / var(i) );
+    Factor A_Ii = (_pancakes[i] * factor(I).inverse() * _phis[i][_I].inverse()).partSum( Ivars / var(i) );
     Factor quot = A_I.divided_by(A_Ii);
+    if( props.damping != 0.0 )
+        quot = (quot^(1.0 - props.damping)) * (_phis[i][_I]^props.damping);
 
-    piet *= quot.divided_by( _phis[i][_I] ).normalized( Prob::NORMPROB );
-    _phis[i][_I] = quot.normalized( Prob::NORMPROB );
+    piet *= quot.divided_by( _phis[i][_I] ).normalized();
+    _phis[i][_I] = quot.normalized();
 
-    piet.normalize( Prob::NORMPROB );
+    piet.normalize();
 
     if( piet.hasNaNs() ) {
-        cout << "LC::NewPancake(" << i << ", " << _I << "):  has NaNs!" << endl;
+        cout << Name << "::NewPancake(" << i << ", " << _I << "):  has NaNs!" << endl;
         hasNaNs = true;
     }
 
@@ -275,8 +273,8 @@ double LC::run() {
     Diffs diffs(nrVars(), 1.0);
 
     double md = InitCavityDists( props.cavainame, props.cavaiopts );
-    if( md > maxdiff )
-        maxdiff = md;
+    if( md > _maxdiff )
+        _maxdiff = md;
 
     vector<Factor> old_beliefs;
     for(size_t i=0; i < nrVars(); i++ )
@@ -289,8 +287,8 @@ double LC::run() {
             break;
         }
     if( hasNaNs ) {
-        cout << "LC::run:  initial _pancakes has NaNs!" << endl;
-        return -1.0;
+        cout << Name << "::run:  initial _pancakes has NaNs!" << endl;
+        return 1.0;
     }
 
     size_t nredges = nrEdges();
@@ -300,11 +298,9 @@ double LC::run() {
         foreach( const Neighbor &I, nbV(i) )
             update_seq.push_back( Edge( i, I.iter ) );
 
-    size_t iter = 0;
-
     // 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( iter=0; iter < props.maxiter && diffs.maxDiff() > props.tol; iter++ ) {
+    for( _iters=0; _iters < props.maxiter && diffs.maxDiff() > props.tol; _iters++ ) {
         // Sequential updates
         if( props.updates == Properties::UpdateType::SEQRND )
             random_shuffle( update_seq.begin(), update_seq.end() );
@@ -314,7 +310,7 @@ double LC::run() {
             size_t _I = update_seq[t].second;
             _pancakes[i] = NewPancake( i, _I, hasNaNs);
             if( hasNaNs )
-                return -1.0;
+                return 1.0;
             CalcBelief( i );
         }
 
@@ -325,21 +321,21 @@ double LC::run() {
         }
 
         if( props.verbose >= 3 )
-            cout << "LC::run:  maxdiff " << diffs.maxDiff() << " after " << iter+1 << " passes" << endl;
+            cout << Name << "::run:  maxdiff " << diffs.maxDiff() << " after " << _iters+1 << " passes" << endl;
     }
 
-    if( diffs.maxDiff() > maxdiff )
-        maxdiff = diffs.maxDiff();
+    if( diffs.maxDiff() > _maxdiff )
+        _maxdiff = diffs.maxDiff();
 
     if( props.verbose >= 1 ) {
         if( diffs.maxDiff() > props.tol ) {
             if( props.verbose == 1 )
                 cout << endl;
-                cout << "LC::run:  WARNING: not converged within " << props.maxiter << " passes (" << toc() - tic << " clocks)...final maxdiff:" << diffs.maxDiff() << endl;
+                cout << Name << "::run:  WARNING: not converged within " << props.maxiter << " passes (" << toc() - tic << " seconds)...final maxdiff:" << diffs.maxDiff() << endl;
         } else {
             if( props.verbose >= 2 )
-                cout << "LC::run:  ";
-                cout << "converged in " << iter << " passes (" << toc() - tic << " clocks)." << endl;
+                cout << Name << "::run:  ";
+                cout << "converged in " << _iters << " passes (" << toc() - tic << " seconds)." << endl;
         }
     }