Replaced complex numbers by real numbers
authorJoris Mooij <jorism@osun.tuebingen.mpg.de>
Fri, 19 Sep 2008 11:23:20 +0000 (13:23 +0200)
committerJoris Mooij <jorism@osun.tuebingen.mpg.de>
Fri, 19 Sep 2008 11:23:20 +0000 (13:23 +0200)
18 files changed:
include/dai/bp.h
include/dai/daialg.h
include/dai/factor.h
include/dai/hak.h
include/dai/jtree.h
include/dai/lc.h
include/dai/mf.h
include/dai/mr.h
include/dai/prob.h
include/dai/treeep.h
matlab/dai.cpp
src/bp.cpp
src/daialg.cpp
src/hak.cpp
src/jtree.cpp
src/mf.cpp
src/treeep.cpp
tests/test.cpp

index 9bb1b51..c972949 100644 (file)
@@ -91,7 +91,7 @@ class BP : public DAIAlgFG {
         Factor belief (const Var &n) const;
         Factor belief (const VarSet &n) const;
         std::vector<Factor> beliefs() const;
-        Complex logZ() const;
+        Real logZ() const;
 
         void init( const VarSet &ns );
         void undoProbs( const VarSet &ns ) { FactorGraph::undoProbs(ns); init(ns); }
index 32c1ece..049d8a4 100644 (file)
@@ -132,7 +132,7 @@ class InfAlg {
         virtual std::vector<Factor> beliefs() const = 0;
 
         /// Get log partition sum
-        virtual Complex logZ() const = 0;
+        virtual Real logZ() const = 0;
 
         /// Clear messages and beliefs
         virtual void init() = 0;
index 41b5eaf..6f2e5b0 100644 (file)
@@ -36,16 +36,15 @@ namespace dai {
 
 template<typename T> class      TFactor;
 typedef TFactor<Real>           Factor;
-typedef TFactor<Complex>        CFactor;
 
 
 // predefine friends
 template<typename T> Real           dist( const TFactor<T> & x, const TFactor<T> & y, Prob::DistType dt );
-template<typename T> Complex        KL_dist( const TFactor<T> & p, const TFactor<T> & q );
+template<typename T> Real           KL_dist( const TFactor<T> & p, const TFactor<T> & q );
 template<typename T> std::ostream&  operator<< (std::ostream& os, const TFactor<T>& P);
 
         
-// T should be castable from and to double and to complex
+// T should be castable from and to double
 template <typename T> class TFactor {
     protected:
         VarSet      _vs;
@@ -206,13 +205,6 @@ template <typename T> class TFactor {
             return l0; 
         }
 
-        CFactor clog0() const {
-            CFactor l0; 
-            l0._vs = _vs; 
-            l0._p = _p.clog0(); 
-            return l0; 
-        }
-
         T normalize( typename Prob::NormType norm ) { return _p.normalize( norm ); }
         TFactor<T> normalized( typename Prob::NormType norm ) const { 
             TFactor<T> result;
@@ -247,7 +239,7 @@ template <typename T> class TFactor {
         T totalSum() const { return _p.totalSum(); }
         T maxAbs() const { return _p.maxAbs(); }
         T maxVal() const { return _p.maxVal(); }
-        Complex entropy() const { return _p.entropy(); }
+        Real entropy() const { return _p.entropy(); }
         T strength( const Var &i, const Var &j ) const;
 
         friend Real dist( const TFactor<T> & x, const TFactor<T> & y, Prob::DistType dt ) {
@@ -260,7 +252,7 @@ template <typename T> class TFactor {
                 return dist( x._p, y._p, dt );
             }
         }
-        friend Complex KL_dist <> (const TFactor<T> & p, const TFactor<T> & q);
+        friend Real KL_dist <> (const TFactor<T> & p, const TFactor<T> & q);
         template<class U> friend std::ostream& operator<< (std::ostream& os, const TFactor<U>& P);
 };
 
@@ -302,7 +294,7 @@ template<typename T> TFactor<T> TFactor<T>::operator* (const TFactor<T>& Q) cons
 }
 
 
-template<typename T> Complex KL_dist(const TFactor<T> & P, const TFactor<T> & Q) {
+template<typename T> Real KL_dist(const TFactor<T> & P, const TFactor<T> & Q) {
     if( P._vs.empty() || Q._vs.empty() )
         return -1;
     else {
index f1adfc4..9420298 100644 (file)
@@ -88,7 +88,7 @@ class HAK : public DAIAlgRG {
         Factor belief( const Var &n ) const;
         Factor belief( const VarSet &ns ) const;
         std::vector<Factor> beliefs() const;
-        Complex logZ () const;
+        Real logZ () const;
 
         void init( const VarSet &ns );
         void undoProbs( const VarSet &ns ) { RegionGraph::undoProbs( ns ); init( ns ); }
index e6bf12f..3244dd9 100644 (file)
@@ -79,7 +79,7 @@ class JTree : public DAIAlgRG {
         Factor belief( const Var &n ) const;
         Factor belief( const VarSet &ns ) const;
         std::vector<Factor> beliefs() const;
-        Complex logZ() const;
+        Real logZ() const;
 
         void init( const VarSet &/*ns*/ ) {}
         void undoProbs( const VarSet &ns ) { RegionGraph::undoProbs( ns ); init( ns ); }
index fbfaf9a..4af42dc 100644 (file)
@@ -87,7 +87,7 @@ class LC : public DAIAlgFG {
             return Factor(); 
         }
         std::vector<Factor> beliefs() const { return _beliefs; }
-        Complex logZ() const { 
+        Real logZ() const { 
             DAI_THROW(NOT_IMPLEMENTED);
             return 0.0; 
         }
index 7f35546..56e0272 100644 (file)
@@ -64,7 +64,7 @@ class MF : public DAIAlgFG {
         Factor belief (const Var &n) const;
         Factor belief (const VarSet &ns) const;
         std::vector<Factor> beliefs() const;
-        Complex logZ() const;
+        Real logZ() const;
 
         void init( const VarSet &ns );
         void undoProbs( const VarSet &ns ) { FactorGraph::undoProbs(ns); init(ns); }
index ef7b6b0..2dfb2b9 100644 (file)
@@ -77,7 +77,7 @@ class MR : public DAIAlgFG {
             return Factor(); 
         }
         std::vector<Factor> beliefs() const;
-        Complex logZ() const { 
+        Real logZ() const { 
             DAI_THROW(NOT_IMPLEMENTED);
             return 0.0; 
         }
index d11cda1..6037339 100644 (file)
@@ -23,7 +23,6 @@
 #define __defined_libdai_prob_h
 
 
-#include <complex>
 #include <cmath>
 #include <vector>
 #include <ostream>
@@ -38,15 +37,13 @@ namespace dai {
 
 
 typedef double                  Real;
-typedef std::complex<double>    Complex;
 
 template<typename T> class      TProb;
 typedef TProb<Real>             Prob;
-typedef TProb<Complex>          CProb;
 
 
 /// TProb<T> implements a probability vector of type T.
-/// T should be castable from and to double and to complex.
+/// T should be castable from and to double.
 template <typename T> class TProb {
     protected:
         /// The entries
@@ -54,7 +51,7 @@ template <typename T> class TProb {
 
     private:
         /// Calculate x times log(x), or 0 if x == 0
-        Complex xlogx( Real x ) const { return( x == 0.0 ? 0.0 : Complex(x) * std::log(Complex(x))); }
+        Real xlogx( Real x ) const { return( x == 0.0 ? 0.0 : x * std::log(x)); }
 
     public:
         /// NORMPROB means that the sum of all entries should be 1
@@ -331,15 +328,6 @@ template <typename T> class TProb {
             return l0;
         }
 
-        /// Pointwise (complex) log (or 0 if == 0)
-/*      CProb clog0() const {
-            CProb l0;
-            l0._p.reserve( size() );
-            for( size_t i = 0; i < size(); i++ )
-                l0._p.push_back( (_p[i] == 0.0) ? 0.0 : std::log( Complex( _p[i] ) ) );
-            return l0;
-        }*/
-
         /// Return distance of p and q
         friend Real dist( const TProb<T> & p, const TProb<T> & q, DistType dt ) {
 #ifdef DAI_DEBUG
@@ -369,16 +357,16 @@ template <typename T> class TProb {
             return result;
         }
 
-        /// Return (complex) Kullback-Leibler distance with q
-        friend Complex KL_dist( const TProb<T> & p, const TProb<T> & q ) {
+        /// Return Kullback-Leibler distance with q
+        friend Real KL_dist( const TProb<T> & p, const TProb<T> & q ) {
 #ifdef DAI_DEBUG
             assert( p.size() == q.size() );
 #endif
-            Complex result = 0.0;
+            Real result = 0.0;
             for( size_t i = 0; i < p.size(); i++ ) {
                 if( (Real) p[i] != 0.0 ) {
-                    Complex p_i = p[i];
-                    Complex q_i = q[i];
+                    Real p_i = p[i];
+                    Real q_i = q[i];
                     result += p_i * (std::log(p_i) - std::log(q_i));
                 }
             }
@@ -444,9 +432,9 @@ template <typename T> class TProb {
             return (std::find_if( _p.begin(), _p.end(), std::bind2nd( std::less_equal<Real>(), 0.0 ) ) != _p.end());
         }
 
-        /// Returns (complex) entropy
-        Complex entropy() const {
-            Complex S = 0.0;
+        /// Returns entropy
+        Real entropy() const {
+            Real S = 0.0;
             for( size_t i = 0; i < size(); i++ )
                 S -= xlogx(_p[i]);
             return S;
index 4d3e417..edef1ce 100644 (file)
@@ -111,7 +111,7 @@ class TreeEP : public JTree {
         std::string identify() const;
         void init();
         double run();
-        Complex logZ() const;
+        Real logZ() const;
 
         bool offtree( size_t I ) const { return (fac2OR[I] == -1U); }
 
index c4e766b..7ad11b8 100644 (file)
@@ -103,7 +103,7 @@ void mexFunction( int nlhs, mxArray *plhs[], int nrhs, const mxArray*prhs[] )
 
 
     // Save logZ
-    double logZ = real( obj->logZ() );
+    double logZ = obj->logZ();
 
     // Save maxdiff
     double maxdiff = obj->MaxDiff();
index 79d751f..580a4c6 100644 (file)
@@ -390,10 +390,10 @@ Factor BP::beliefF (size_t I) const {
 }
 
 
-Complex BP::logZ() const {
-    Complex sum = 0.0;
+Real BP::logZ() const {
+    Real sum = 0.0;
     for(size_t i = 0; i < nrVars(); ++i )
-        sum += Complex(1.0 - nbV(i).size()) * beliefV(i).entropy();
+        sum += (1.0 - nbV(i).size()) * beliefV(i).entropy();
     for( size_t I = 0; I < nrFactors(); ++I )
         sum -= KL_dist( beliefF(I), factor(I) );
     return sum;
index a688f36..062de76 100644 (file)
@@ -38,7 +38,7 @@ Factor calcMarginal( const InfAlg & obj, const VarSet & ns, bool reInit ) {
     if( !reInit )
         clamped->init();
 
-    Complex logZ0;
+    Real logZ0 = 0.0;
     for( State s(ns); s.valid(); s++ ) {
         // save unclamped factors connected to ns
         clamped->saveProbs( ns );
@@ -54,18 +54,16 @@ Factor calcMarginal( const InfAlg & obj, const VarSet & ns, bool reInit ) {
             clamped->init();
         clamped->run();
 
-        Complex Z;
+        Real Z;
         if( s == 0 ) {
             logZ0 = clamped->logZ();
             Z = 1.0;
         } else {
             // subtract logZ0 to avoid very large numbers
             Z = exp(clamped->logZ() - logZ0);
-            if( fabs(imag(Z)) > 1e-5 )
-                cout << "Marginal:: WARNING: complex Z (" << Z << ")" << endl;
         }
 
-        Pns[s] = real(Z);
+        Pns[s] = Z;
         
         // restore clamped factors
         clamped->undoProbs( ns );
@@ -98,7 +96,7 @@ vector<Factor> calcPairBeliefs( const InfAlg & obj, const VarSet& ns, bool reIni
     if( !reInit )
         clamped->init();
 
-    Complex logZ0;
+    Real logZ0 = 0.0;
     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++ ) {
@@ -120,10 +118,7 @@ vector<Factor> calcPairBeliefs( const InfAlg & obj, const VarSet& ns, bool reIni
                 logZ0 = clamped->logZ();
             } else {
                 // subtract logZ0 to avoid very large numbers
-                Complex Z = exp(clamped->logZ() - logZ0);
-                if( fabs(imag(Z)) > 1e-5 )
-                    cout << "calcPairBelief::  Warning: complex Z: " << Z << endl;
-                Z_xj = real(Z);
+                Z_xj = exp(clamped->logZ() - logZ0);
             }
 
             for( size_t k = 0; k < N; k++ ) 
@@ -177,7 +172,7 @@ vector<Factor> calcPairBeliefsNew( const InfAlg & obj, const VarSet& ns, bool re
     if( !reInit )
         clamped->init();
 
-    Complex logZ0;
+    Real logZ0 = 0.0;
     VarSet::const_iterator nj = ns.begin();
     for( long j = 0; j < (long)ns.size() - 1; j++, nj++ ) {
         size_t k = 0;
@@ -201,10 +196,7 @@ vector<Factor> calcPairBeliefsNew( const InfAlg & obj, const VarSet& ns, bool re
                         logZ0 = clamped->logZ();
                     } else {
                         // subtract logZ0 to avoid very large numbers
-                        Complex Z = exp(clamped->logZ() - logZ0);
-                        if( fabs(imag(Z)) > 1e-5 )
-                            cout << "calcPairBelief::  Warning: complex Z: " << Z << endl;
-                        Z_xj = real(Z);
+                        Z_xj = exp(clamped->logZ() - logZ0);
                     }
 
                     // we assume that j.label() < k.label()
index 5a5a19b..5570dee 100644 (file)
@@ -415,12 +415,12 @@ vector<Factor> HAK::beliefs() const {
 }
 
 
-Complex HAK::logZ() const {
-    Complex sum = 0.0;
+Real HAK::logZ() const {
+    Real sum = 0.0;
     for( size_t beta = 0; beta < nrIRs(); beta++ )
-        sum += Complex(IR(beta).c()) * Qb(beta).entropy();
+        sum += IR(beta).c() * Qb(beta).entropy();
     for( size_t alpha = 0; alpha < nrORs(); alpha++ ) {
-        sum += Complex(OR(alpha).c()) * Qa(alpha).entropy();
+        sum += OR(alpha).c() * Qa(alpha).entropy();
         sum += (OR(alpha).log0() * Qa(alpha)).totalSum();
     }
     return sum;
index c7652e1..ce619f7 100644 (file)
@@ -292,12 +292,12 @@ double JTree::run() {
 }
 
 
-Complex JTree::logZ() const {
-    Complex sum = 0.0;
+Real JTree::logZ() const {
+    Real sum = 0.0;
     for( size_t beta = 0; beta < nrIRs(); beta++ )
-        sum += Complex(IR(beta).c()) * _Qb[beta].entropy();
+        sum += IR(beta).c() * _Qb[beta].entropy();
     for( size_t alpha = 0; alpha < nrORs(); alpha++ ) {
-        sum += Complex(OR(alpha).c()) * _Qa[alpha].entropy();
+        sum += OR(alpha).c() * _Qa[alpha].entropy();
         sum += (OR(alpha).log0() * _Qa[alpha]).totalSum();
     }
     return sum;
index 5587a3e..921692b 100644 (file)
@@ -168,8 +168,8 @@ vector<Factor> MF::beliefs() const {
 }
 
 
-Complex MF::logZ() const {
-    Complex sum = 0.0;
+Real MF::logZ() const {
+    Real sum = 0.0;
     
     for(size_t i=0; i < nrVars(); i++ )
         sum -= beliefV(i).entropy();
@@ -181,7 +181,7 @@ Complex MF::logZ() const {
         Factor piet;
         piet = factor(I).log0();
         piet *= henk;
-        sum -= Complex( piet.totalSum() );
+        sum -= piet.totalSum();
     }
 
     return -sum;
index 1fee924..2a60738 100644 (file)
@@ -210,7 +210,7 @@ TreeEP::TreeEP( const FactorGraph &fg, const Properties &opts ) : JTree(fg, opts
                         if( piet.vars() >> (v_i | *j) ) {
                             piet = piet.marginal( v_i | *j );
                             Factor pietf = piet.marginal(v_i) * piet.marginal(*j);
-                            wg[UEdge(i,findVar(*j))] = real( KL_dist( piet, pietf ) );
+                            wg[UEdge(i,findVar(*j))] = KL_dist( piet, pietf );
                         } else
                             wg[UEdge(i,findVar(*j))] = 0;
                     }
@@ -455,14 +455,14 @@ double TreeEP::run() {
 }
 
 
-Complex TreeEP::logZ() const {
+Real TreeEP::logZ() const {
     double sum = 0.0;
 
     // entropy of the tree
     for( size_t beta = 0; beta < nrIRs(); beta++ )
-        sum -= real(_Qb[beta].entropy());
+        sum -= _Qb[beta].entropy();
     for( size_t alpha = 0; alpha < nrORs(); alpha++ )
-        sum += real(_Qa[alpha].entropy());
+        sum += _Qa[alpha].entropy();
 
     // energy of the on-tree factors
     for( size_t alpha = 0; alpha < nrORs(); alpha++ )
index 54294ed..fe42f43 100644 (file)
@@ -107,7 +107,7 @@ class TestAI {
                 obj->run();
                 time += toc() - tic;
                 try {
-                    logZ = real(obj->logZ());
+                    logZ = obj->logZ();
                     has_logZ = true;
                 } catch( Exception &e ) {
                     has_logZ = false;