Replaced doubles by Reals, fixed two bugs
[libdai.git] / src / mf.cpp
1 /* This file is part of libDAI - http://www.libdai.org/
2 *
3 * libDAI is licensed under the terms of the GNU General Public License version
4 * 2, or (at your option) any later version. libDAI is distributed without any
5 * warranty. See the file COPYING for more details.
6 *
7 * Copyright (C) 2006-2009 Joris Mooij [joris dot mooij at libdai dot org]
8 * Copyright (C) 2006-2007 Radboud University Nijmegen, The Netherlands
9 */
10
11
12 #include <iostream>
13 #include <sstream>
14 #include <map>
15 #include <set>
16 #include <dai/mf.h>
17 #include <dai/util.h>
18
19
20 namespace dai {
21
22
23 using namespace std;
24
25
26 const char *MF::Name = "MF";
27
28
29 void MF::setProperties( const PropertySet &opts ) {
30 DAI_ASSERT( opts.hasKey("tol") );
31 DAI_ASSERT( opts.hasKey("maxiter") );
32
33 props.tol = opts.getStringAs<Real>("tol");
34 props.maxiter = opts.getStringAs<size_t>("maxiter");
35 if( opts.hasKey("verbose") )
36 props.verbose = opts.getStringAs<size_t>("verbose");
37 else
38 props.verbose = 0U;
39 if( opts.hasKey("damping") )
40 props.damping = opts.getStringAs<Real>("damping");
41 else
42 props.damping = 0.0;
43 }
44
45
46 PropertySet MF::getProperties() const {
47 PropertySet opts;
48 opts.Set( "tol", props.tol );
49 opts.Set( "maxiter", props.maxiter );
50 opts.Set( "verbose", props.verbose );
51 opts.Set( "damping", props.damping );
52 return opts;
53 }
54
55
56 string MF::printProperties() const {
57 stringstream s( stringstream::out );
58 s << "[";
59 s << "tol=" << props.tol << ",";
60 s << "maxiter=" << props.maxiter << ",";
61 s << "verbose=" << props.verbose << ",";
62 s << "damping=" << props.damping << "]";
63 return s.str();
64 }
65
66
67 void MF::construct() {
68 // create beliefs
69 _beliefs.clear();
70 _beliefs.reserve( nrVars() );
71 for( size_t i = 0; i < nrVars(); ++i )
72 _beliefs.push_back( Factor( var(i) ) );
73 }
74
75
76 string MF::identify() const {
77 return string(Name) + printProperties();
78 }
79
80
81 void MF::init() {
82 for( vector<Factor>::iterator qi = _beliefs.begin(); qi != _beliefs.end(); qi++ )
83 qi->fill(1.0);
84 }
85
86
87 Real MF::run() {
88 double tic = toc();
89
90 if( props.verbose >= 1 )
91 cerr << "Starting " << identify() << "...";
92
93 size_t pass_size = _beliefs.size();
94 Diffs diffs(pass_size * 3, 1.0);
95
96 size_t t=0;
97 for( t=0; t < (props.maxiter*pass_size) && diffs.maxDiff() > props.tol; t++ ) {
98 // choose random Var i
99 size_t i = (size_t) (nrVars() * rnd_uniform());
100
101 Factor jan;
102 Factor piet;
103 foreach( const Neighbor &I, nbV(i) ) {
104 Factor henk;
105 foreach( const Neighbor &j, nbF(I) ) // for all j in I \ i
106 if( j != i )
107 henk *= _beliefs[j];
108 piet = factor(I).log(true);
109 piet *= henk;
110 piet = piet.marginal(var(i), false);
111 piet = piet.exp();
112 jan *= piet;
113 }
114
115 jan.normalize();
116
117 if( jan.hasNaNs() ) {
118 cerr << Name << "::run(): ERROR: jan has NaNs!" << endl;
119 return 1.0;
120 }
121
122 if( props.damping != 0.0 )
123 jan = (jan^(1.0 - props.damping)) * (_beliefs[i]^props.damping);
124 diffs.push( dist( jan, _beliefs[i], Prob::DISTLINF ) );
125
126 _beliefs[i] = jan;
127 }
128
129 _iters = t / pass_size;
130 if( diffs.maxDiff() > _maxdiff )
131 _maxdiff = diffs.maxDiff();
132
133 if( props.verbose >= 1 ) {
134 if( diffs.maxDiff() > props.tol ) {
135 if( props.verbose == 1 )
136 cerr << endl;
137 cerr << Name << "::run: WARNING: not converged within " << props.maxiter << " passes (" << toc() - tic << " seconds)...final maxdiff:" << diffs.maxDiff() << endl;
138 } else {
139 if( props.verbose >= 2 )
140 cerr << Name << "::run: ";
141 cerr << "converged in " << t / pass_size << " passes (" << toc() - tic << " seconds)." << endl;
142 }
143 }
144
145 return diffs.maxDiff();
146 }
147
148
149 Factor MF::beliefV( size_t i ) const {
150 Factor piet;
151 piet = _beliefs[i];
152 piet.normalize();
153 return(piet);
154 }
155
156
157 Factor MF::belief (const VarSet &ns) const {
158 if( ns.size() == 1 )
159 return belief( *(ns.begin()) );
160 else {
161 DAI_ASSERT( ns.size() == 1 );
162 return Factor();
163 }
164 }
165
166
167 Factor MF::belief (const Var &n) const {
168 return( beliefV( findVar( n ) ) );
169 }
170
171
172 vector<Factor> MF::beliefs() const {
173 vector<Factor> result;
174 for( size_t i = 0; i < nrVars(); i++ )
175 result.push_back( beliefV(i) );
176 return result;
177 }
178
179
180 Real MF::logZ() const {
181 Real s = 0.0;
182
183 for(size_t i=0; i < nrVars(); i++ )
184 s -= beliefV(i).entropy();
185 for(size_t I=0; I < nrFactors(); I++ ) {
186 Factor henk;
187 foreach( const Neighbor &j, nbF(I) ) // for all j in I
188 henk *= _beliefs[j];
189 henk.normalize();
190 Factor piet;
191 piet = factor(I).log(true);
192 piet *= henk;
193 s -= piet.sum();
194 }
195
196 return -s;
197 }
198
199
200 void MF::init( const VarSet &ns ) {
201 for( size_t i = 0; i < nrVars(); i++ ) {
202 if( ns.contains(var(i) ) )
203 _beliefs[i].fill( 1.0 );
204 }
205 }
206
207
208 } // end of namespace dai