1 /* This file is part of libDAI - http://www.libdai.org/
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.
7 * Copyright (C) 2006-2010 Joris Mooij [joris dot mooij at libdai dot org]
8 * Copyright (C) 2006-2007 Radboud University Nijmegen, The Netherlands
17 #include <boost/program_options.hpp>
18 #include <dai/factorgraph.h>
19 #include <dai/weightedgraph.h>
21 #include <dai/graph.h>
23 #include <dai/properties.h>
28 namespace po
= boost::program_options
;
31 /// Possible factor types
32 DAI_ENUM(FactorType
,ISING
,EXPGAUSS
,POTTS
);
35 /// Creates a factor graph from a pairwise interactions graph
36 /** \param G Graph describing interactions between variables
37 * \param ft Type of factors to use for interactions
38 * \param states Number of states of the variables
39 * \param props Additional properties for generating the interactions
41 FactorGraph
createFG( const GraphAL
&G
, FactorType ft
, size_t states
, const PropertySet
&props
) {
42 size_t N
= G
.nrNodes();
45 DAI_ASSERT( ft
!= FactorType::ISING
);
47 // Get inverse temperature
49 if( ft
!= FactorType::ISING
)
50 beta
= props
.GetAs
<Real
>("beta");
52 // Get properties for Ising factors
57 if( ft
== FactorType::ISING
) {
58 mean_h
= props
.GetAs
<Real
>("mean_th");
59 sigma_h
= props
.GetAs
<Real
>("sigma_th");
60 mean_J
= props
.GetAs
<Real
>("mean_w");
61 sigma_J
= props
.GetAs
<Real
>("sigma_w");
67 for( size_t i
= 0; i
< N
; i
++ )
68 vars
.push_back( Var( i
, states
) );
71 vector
<Factor
> factors
;
72 factors
.reserve( G
.nrEdges() + N
);
74 for( size_t i
= 0; i
< N
; i
++ )
75 foreach( const GraphAL::Neighbor
&j
, G
.nb(i
) )
77 if( ft
== FactorType::POTTS
)
78 factors
.push_back( createFactorPotts( vars
[i
], vars
[j
], beta
) );
79 else if( ft
== FactorType::EXPGAUSS
)
80 factors
.push_back( createFactorExpGauss( VarSet( vars
[i
], vars
[j
] ), beta
) );
81 else if( ft
== FactorType::ISING
) {
82 Real J
= rnd_stdnormal() * sigma_J
+ mean_J
;
83 factors
.push_back( createFactorIsing( vars
[i
], vars
[j
], J
) );
87 if( ft
== FactorType::ISING
)
88 for( size_t i
= 0; i
< N
; i
++ ) {
89 Real h
= rnd_stdnormal() * sigma_h
+ mean_h
;
90 factors
.push_back( createFactorIsing( vars
[i
], h
) );
93 return FactorGraph( factors
.begin(), factors
.end(), vars
.begin(), vars
.end(), factors
.size(), vars
.size() );
97 /// Return a random factor graph with higher-order interactions
98 /** \param N number of variables
99 * \param M number of factors
100 * \param k number of variables that each factor depends on
101 * \param beta standard-deviation of Gaussian log-factor entries
103 FactorGraph
createHOIFG( size_t N
, size_t M
, size_t k
, Real beta
) {
105 vector
<Factor
> factors
;
108 for( size_t i
= 0; i
< N
; i
++ )
109 vars
.push_back(Var(i
,2));
111 for( size_t I
= 0; I
< M
; I
++ ) {
113 while( vars
.size() < k
) {
115 size_t newind
= (size_t)(N
* rnd_uniform());
116 Var newvar
= Var(newind
, 2);
117 if( !vars
.contains( newvar
) ) {
123 factors
.push_back( createFactorExpGauss( vars
, beta
) );
126 return FactorGraph( factors
.begin(), factors
.end(), vars
.begin(), vars
.end(), factors
.size(), vars
.size() );
130 /// Creates a regular random bipartite graph
131 /** \param N1 = number of nodes of type 1
132 * \param d1 = size of neighborhoods of nodes of type 1
133 * \param N2 = number of nodes of type 2
134 * \param d2 = size of neighborhoods of nodes of type 2
135 * \note asserts that N1 * d1 == N2 * d2
137 BipartiteGraph
createRandomBipartiteGraph( size_t N1
, size_t N2
, size_t d1
, size_t d2
) {
140 DAI_ASSERT( N1
* d1
== N2
* d2
);
142 // build lists of degree-repeated vertex numbers
143 std::vector
<size_t> stubs1( N1
*d1
, 0 );
144 for( size_t n1
= 0; n1
< N1
; n1
++ )
145 for( size_t t
= 0; t
< d1
; t
++ )
146 stubs1
[n1
*d1
+ t
] = n1
;
148 // build lists of degree-repeated vertex numbers
149 std::vector
<size_t> stubs2( N2
*d2
, 0 );
150 for( size_t n2
= 0; n2
< N2
; n2
++ )
151 for( size_t t
= 0; t
< d2
; t
++ )
152 stubs2
[n2
*d2
+ t
] = n2
;
155 random_shuffle( stubs1
.begin(), stubs1
.end() );
156 random_shuffle( stubs2
.begin(), stubs2
.end() );
159 vector
<BipartiteGraph::Edge
> edges
;
160 edges
.reserve( N1
*d1
);
161 for( size_t e
= 0; e
< N1
*d1
; e
++ )
162 edges
.push_back( BipartiteGraph::Edge(stubs1
[e
], stubs2
[e
]) );
164 // finish construction
165 G
.construct( N1
, N2
, edges
.begin(), edges
.end() );
171 /// Returns x**n % p, assuming p is prime
172 int powmod (int x
, int n
, int p
) {
174 for( int m
= 0; m
< n
; m
++ )
180 /// Returns order of x in GF(p) with p prime
181 size_t order (int x
, int p
) {
183 DAI_ASSERT( x
!= 0 );
187 prod
= (prod
* x
) % p
;
189 } while( prod
!= 1 );
194 /// Returns whether n is a prime number
195 bool isPrime (size_t n
) {
197 for( size_t k
= 2; (k
< n
) && result
; k
++ )
204 /// Constructs a regular LDPC graph with N=6, j=2, K=4, k=3
205 BipartiteGraph
createSmallLDPCGraph() {
207 size_t N
=4, j
=3, K
=4; // k=3;
209 typedef BipartiteGraph::Edge Edge
;
211 edges
.reserve( N
*j
);
212 edges
.push_back( Edge(0,0) ); edges
.push_back( Edge(1,0) ); edges
.push_back( Edge(2,0) );
213 edges
.push_back( Edge(0,1) ); edges
.push_back( Edge(1,1) ); edges
.push_back( Edge(3,1) );
214 edges
.push_back( Edge(0,2) ); edges
.push_back( Edge(2,2) ); edges
.push_back( Edge(3,2) );
215 edges
.push_back( Edge(1,3) ); edges
.push_back( Edge(2,3) ); edges
.push_back( Edge(3,3) );
217 // finish construction
218 G
.construct( N
, K
, edges
.begin(), edges
.end() );
224 /// Creates group-structured LDPC code
225 /** Use construction described in "A Class of Group-Structured LDPC Codes"
226 * by R. M. Tanner, D. Sridhara and T. Fuja
227 * Proceedings of ICSTA, 2001
229 * Example parameters: (p,j,k) = (31,3,5)
234 * j and k must be divisors of p-1
236 BipartiteGraph
createGroupStructuredLDPCGraph( size_t p
, size_t j
, size_t k
) {
244 for( a
= 2; a
< p
; a
++ )
245 if( order(a
,p
) == k
)
247 DAI_ASSERT( a
!= p
);
248 for( b
= 2; b
< p
; b
++ )
249 if( order(b
,p
) == j
)
251 DAI_ASSERT( b
!= p
);
252 // cout << "# order(a=" << a << ") = " << order(a,p) << endl;
253 // cout << "# order(b=" << b << ") = " << order(b,p) << endl;
255 DAI_ASSERT( N
* n
== K
* k
);
257 typedef BipartiteGraph::Edge Edge
;
259 edges
.reserve( N
* n
);
261 for( size_t s
= 0; s
< j
; s
++ )
262 for( size_t t
= 0; t
< k
; t
++ ) {
263 size_t P
= (powmod(b
,s
,p
) * powmod(a
,t
,p
)) % p
;
264 for( size_t m
= 0; m
< p
; m
++ )
265 edges
.push_back( Edge(t
*p
+ m
, s
*p
+ ((m
+ P
) % p
)) );
268 // finish construction
269 G
.construct( N
, K
, edges
.begin(), edges
.end() );
275 // Constructs a parity check table
276 void createParityCheck( Real
*result
, size_t n
, Real eps
) {
278 for( size_t i
= 0; i
< N
; i
++ ) {
280 for( size_t t
= 0; t
< n
; t
++ )
286 result
[i
] = 1.0 - eps
;
292 /// Predefined names of various factor graph types
293 const char *FULL_TYPE
= "FULL";
294 const char *DREG_TYPE
= "DREG";
295 const char *LOOP_TYPE
= "LOOP";
296 const char *TREE_TYPE
= "TREE";
297 const char *GRID_TYPE
= "GRID";
298 const char *GRID3D_TYPE
= "GRID3D";
299 const char *HOI_TYPE
= "HOI";
300 const char *LDPC_TYPE
= "LDPC";
303 /// Possible LDPC structures
304 DAI_ENUM(LDPCType
,SMALL
,GROUP
,RANDOM
);
308 int main( int argc
, char *argv
[] ) {
310 // Variables for storing command line arguments
314 size_t d
, N
, K
, k
, j
, n1
, n2
, n3
, prime
;
315 bool periodic
= false;
318 Real beta
, sigma_w
, sigma_th
, mean_w
, mean_th
, noise
;
320 // Declare the supported options.
321 po::options_description
opts("General command line options");
323 ("help", "produce help message")
324 ("seed", po::value
<size_t>(&seed
), "random number seed (tries to read from /dev/urandom if not specified)")
325 ("states", po::value
<size_t>(&states
), "number of states of each variable (default=2 for binary variables)")
328 // Graph structure options
329 po::options_description
opts_graph("Options for specifying graph structure");
330 opts_graph
.add_options()
331 ("type", po::value
<string
>(&type
), "factor graph type (one of 'FULL', 'DREG', 'LOOP', 'TREE', 'GRID', 'GRID3D', 'HOI', 'LDPC')")
332 ("d", po::value
<size_t>(&d
), "variable connectivity (only for type=='DREG');\n\t<d><N> should be even")
333 ("N", po::value
<size_t>(&N
), "number of variables (not for type=='GRID','GRID3D')")
334 ("n1", po::value
<size_t>(&n1
), "width of grid (only for type=='GRID','GRID3D')")
335 ("n2", po::value
<size_t>(&n2
), "height of grid (only for type=='GRID','GRID3D')")
336 ("n3", po::value
<size_t>(&n3
), "length of grid (only for type=='GRID3D')")
337 ("periodic", po::value
<bool>(&periodic
), "periodic grid? (only for type=='GRID','GRID3D'; default=0)")
338 ("K", po::value
<size_t>(&K
), "number of factors (only for type=='HOI','LDPC')")
339 ("k", po::value
<size_t>(&k
), "number of variables per factor (only for type=='HOI','LDPC')")
343 po::options_description
opts_factors("Options for specifying factors");
344 opts_factors
.add_options()
345 ("factors", po::value
<FactorType
>(&ft
), "factor type (one of 'EXPGAUSS','POTTS','ISING')")
346 ("beta", po::value
<Real
>(&beta
), "inverse temperature (ignored for factors=='ISING')")
347 ("mean_w", po::value
<Real
>(&mean_w
), "mean of pairwise interactions w_{ij} (only for factors=='ISING')")
348 ("mean_th", po::value
<Real
>(&mean_th
), "mean of unary interactions th_i (only for factors=='ISING')")
349 ("sigma_w", po::value
<Real
>(&sigma_w
), "stddev of pairwise interactions w_{ij} (only for factors=='ISING')")
350 ("sigma_th", po::value
<Real
>(&sigma_th
), "stddev of unary interactions th_i (only for factors=='ISING'")
354 po::options_description
opts_ldpc("Options for specifying LDPC code factor graphs");
355 opts_ldpc
.add_options()
356 ("ldpc", po::value
<LDPCType
>(&ldpc
), "type of LDPC code (one of 'SMALL','GROUP','RANDOM')")
357 ("j", po::value
<size_t>(&j
), "number of parity checks per bit (only for type=='LDPC')")
358 ("noise", po::value
<Real
>(&noise
), "bitflip probability for binary symmetric channel (only for type=='LDPC')")
359 ("prime", po::value
<size_t>(&prime
), "prime number for construction of LDPC code (only for type=='LDPC' with ldpc='GROUP'))")
363 opts
.add(opts_graph
).add(opts_factors
).add(opts_ldpc
);
365 // Parse command line arguments
366 po::variables_map vm
;
367 po::store(po::parse_command_line(argc
, argv
, opts
), vm
);
370 // Display help message if necessary
371 if( vm
.count("help") || !vm
.count("type") ) {
372 cout
<< "This program is part of libDAI - http://www.libdai.org/" << endl
<< endl
;
373 cout
<< "Usage: ./createfg [options]" << endl
<< endl
;
374 cout
<< "Creates a factor graph according to the specified options." << endl
<< endl
;
376 cout
<< endl
<< opts
<< endl
;
378 cout
<< "The following factor graph types with pairwise interactions can be created:" << endl
;
379 cout
<< "\t'FULL': fully connected graph of <N> variables" << endl
;
380 cout
<< "\t'DREG': random regular graph of <N> variables where each variable is connected with <d> others" << endl
;
381 cout
<< "\t'LOOP': a single loop of <N> variables" << endl
;
382 cout
<< "\t'TREE': random tree-structured (acyclic, connected) graph of <N> variables" << endl
;
383 cout
<< "\t'GRID': 2D grid of <n1>x<n2> variables" << endl
;
384 cout
<< "\t'GRID3D': 3D grid of <n1>x<n2>x<n3> variables" << endl
;
385 cout
<< "The following higher-order interactions factor graphs can be created:" << endl
;
386 cout
<< "\t'HOI': random factor graph consisting of <N> variables and <K> factors," << endl
;
387 cout
<< "\t each factor being an interaction of <k> variables." << endl
;
388 cout
<< "The following LDPC code factor graphs can be created:" << endl
;
389 cout
<< "\t'LDPC': simulates LDPC decoding problem, using an LDPC code of <N> bits and <K>" << endl
;
390 cout
<< "\t parity checks, with <k> bits per check and <j> checks per bit, transmitted" << endl
;
391 cout
<< "\t on a binary symmetric channel with probability <noise> of flipping a bit." << endl
;
392 cout
<< "\t The transmitted codeword has all bits set to zero. The argument 'ldpc'" << endl
;
393 cout
<< "\t determines how the LDPC code is constructed: either using a group structure," << endl
;
394 cout
<< "\t or randomly, or a fixed small code with (N,K,k,j) = (4,4,3,3)." << endl
<< endl
;
396 cout
<< "For all types except type=='LDPC', the factors have to be specified as well." << endl
<< endl
;
398 cout
<< "EXPGAUSS factors (the default) are created by drawing all log-factor entries" << endl
;
399 cout
<< "independently from a Gaussian with mean 0 and standard deviation <beta>." << endl
<< endl
;
401 cout
<< "In case of pairwise interactions, one can also choose POTTS factors, for which" << endl
;
402 cout
<< "the log-factors are simply delta functions multiplied by the strength <beta>." << endl
<< endl
;
404 cout
<< "For pairwise interactions and binary variables, one can also use ISING factors." << endl
;
405 cout
<< "Here variables x1...xN are assumed to be +1/-1--valued, and unary interactions" << endl
;
406 cout
<< "are of the form exp(th*xi) with th drawn from a Gaussian distribution with mean" << endl
;
407 cout
<< "<mean_th> and standard deviation <sigma_th>, and pairwise interactions are of the" << endl
;
408 cout
<< "form exp(w*xi*xj) with w drawn from a Gaussian distribution with mean <mean_w>" << endl
;
409 cout
<< "and standard deviation <sigma_w>." << endl
;
413 // Set default number of states
414 if( !vm
.count("states") )
417 // Set default factor type
418 if( !vm
.count("factors") )
419 ft
= FactorType::EXPGAUSS
;
420 // Check validness of factor type
421 if( ft
== FactorType::POTTS
)
422 if( type
== HOI_TYPE
)
423 throw "For factors=='POTTS', interactions should be pairwise (type!='HOI')";
424 if( ft
== FactorType::ISING
)
425 if( ((states
!= 2) || (type
== HOI_TYPE
)) )
426 throw "For factors=='ISING', variables should be binary (states==2) and interactions should be pairwise (type!='HOI')";
429 if( !vm
.count("seed") ) {
432 infile
.open( "/dev/urandom" );
433 success
= infile
.is_open();
435 infile
.read( (char *)&seed
, sizeof(size_t) / sizeof(char) );
436 success
= infile
.good();
440 throw "Please specify random number seed.";
444 // Set default periodicity
445 if( !vm
.count("periodic") )
448 // Store some options in a PropertySet object
450 if( vm
.count("mean_th") )
451 options
.Set("mean_th", mean_th
);
452 if( vm
.count("sigma_th") )
453 options
.Set("sigma_th", sigma_th
);
454 if( vm
.count("mean_w") )
455 options
.Set("mean_w", mean_w
);
456 if( vm
.count("sigma_w") )
457 options
.Set("sigma_w", sigma_w
);
458 if( vm
.count("beta") )
459 options
.Set("beta", beta
);
461 // Output some comments
462 cout
<< "# Factor graph made by " << argv
[0] << endl
;
463 cout
<< "# type = " << type
<< endl
;
464 cout
<< "# states = " << states
<< endl
;
466 // The factor graph to be constructed
469 #define NEED_ARG(name, desc) do { if(!vm.count(name)) throw "Please specify " desc " with --" name; } while(0);
471 if( type
== FULL_TYPE
|| type
== DREG_TYPE
|| type
== LOOP_TYPE
|| type
== TREE_TYPE
|| type
== GRID_TYPE
|| type
== GRID3D_TYPE
) {
472 // Pairwise interactions
474 // Check command line options
475 if( type
== GRID_TYPE
) {
476 NEED_ARG("n1", "width of grid");
477 NEED_ARG("n2", "height of grid");
479 } else if( type
== GRID3D_TYPE
) {
480 NEED_ARG("n1", "width of grid");
481 NEED_ARG("n2", "height of grid");
482 NEED_ARG("n3", "depth of grid");
485 NEED_ARG("N", "number of variables");
487 if( states
> 2 || ft
== FactorType::POTTS
) {
488 NEED_ARG("beta", "stddev of log-factor entries");
490 NEED_ARG("mean_w", "mean of pairwise interactions");
491 NEED_ARG("mean_th", "mean of unary interactions");
492 NEED_ARG("sigma_w", "stddev of pairwise interactions");
493 NEED_ARG("sigma_th", "stddev of unary interactions");
496 if( type
== DREG_TYPE
)
497 NEED_ARG("d", "connectivity (number of neighboring variables of each variable)");
499 // Build pairwise interaction graph
501 if( type
== FULL_TYPE
)
502 G
= createGraphFull( N
);
503 else if( type
== DREG_TYPE
)
504 G
= createGraphRegular( N
, d
);
505 else if( type
== LOOP_TYPE
)
506 G
= createGraphLoop( N
);
507 else if( type
== TREE_TYPE
)
508 G
= createGraphTree( N
);
509 else if( type
== GRID_TYPE
)
510 G
= createGraphGrid( n1
, n2
, periodic
);
511 else if( type
== GRID3D_TYPE
)
512 G
= createGraphGrid3D( n1
, n2
, n3
, periodic
);
514 // Construct factor graph from pairwise interaction graph
515 fg
= createFG( G
, ft
, states
, options
);
517 // Output some additional comments
518 if( type
== GRID_TYPE
|| type
== GRID3D_TYPE
) {
519 cout
<< "# n1 = " << n1
<< endl
;
520 cout
<< "# n2 = " << n2
<< endl
;
521 if( type
== GRID3D_TYPE
)
522 cout
<< "# n3 = " << n3
<< endl
;
524 if( type
== DREG_TYPE
)
525 cout
<< "# d = " << d
<< endl
;
526 cout
<< "# options = " << options
<< endl
;
527 } else if( type
== HOI_TYPE
) {
528 // Higher order interactions
530 // Check command line arguments
531 NEED_ARG("N", "number of variables");
532 NEED_ARG("K", "number of factors");
533 NEED_ARG("k", "number of variables per factor");
534 NEED_ARG("beta", "stddev of log-factor entries");
536 // Create higher-order interactions factor graph
538 fg
= createHOIFG( N
, K
, k
, beta
);
539 } while( !fg
.isConnected() );
541 // Output some additional comments
542 cout
<< "# K = " << K
<< endl
;
543 cout
<< "# k = " << k
<< endl
;
544 cout
<< "# beta = " << beta
<< endl
;
545 } else if( type
== LDPC_TYPE
) {
548 // Check command line arguments
549 NEED_ARG("ldpc", "type of LDPC code");
550 NEED_ARG("noise", "bitflip probability for binary symmetric channel");
552 // Check more command line arguments (seperately for each LDPC type)
553 if( ldpc
== LDPCType::RANDOM
) {
554 NEED_ARG("N", "number of variables");
555 NEED_ARG("K", "number of factors");
556 NEED_ARG("k", "number of variables per factor");
557 NEED_ARG("j", "number of parity checks per bit");
559 throw "Parameters should satisfy N * j == K * k";
560 } else if( ldpc
== LDPCType::GROUP
) {
561 NEED_ARG("prime", "prime number");
562 NEED_ARG("k", "number of variables per factor");
563 NEED_ARG("j", "number of parity checks per bit");
565 if( !isPrime(prime
) )
566 throw "Parameter <prime> should be prime";
567 if( !((prime
-1) % j
== 0 ) )
568 throw "Parameters should satisfy (prime-1) % j == 0";
569 if( !((prime
-1) % k
== 0 ) )
570 throw "Parameters should satisfy (prime-1) % k == 0";
574 } else if( ldpc
== LDPCType::SMALL
) {
581 // Output some additional comments
582 cout
<< "# N = " << N
<< endl
;
583 cout
<< "# K = " << K
<< endl
;
584 cout
<< "# j = " << j
<< endl
;
585 cout
<< "# k = " << k
<< endl
;
586 if( ldpc
== LDPCType::GROUP
)
587 cout
<< "# prime = " << prime
<< endl
;
588 cout
<< "# noise = " << noise
<< endl
;
590 // Construct likelihood and paritycheck factors
591 Real likelihood
[4] = {1.0 - noise
, noise
, noise
, 1.0 - noise
};
592 Real
*paritycheck
= new Real
[1 << k
];
593 createParityCheck(paritycheck
, k
, 0.0);
595 // Create LDPC structure
596 BipartiteGraph ldpcG
;
599 if( ldpc
== LDPCType::GROUP
)
600 ldpcG
= createGroupStructuredLDPCGraph( prime
, j
, k
);
601 else if( ldpc
== LDPCType::RANDOM
)
602 ldpcG
= createRandomBipartiteGraph( N
, K
, j
, k
);
603 else if( ldpc
== LDPCType::SMALL
)
604 ldpcG
= createSmallLDPCGraph();
607 for( size_t i
= 0; i
< N
; i
++ )
608 if( ldpcG
.nb1(i
).size() != j
)
610 for( size_t I
= 0; I
< K
; I
++ )
611 if( ldpcG
.nb2(I
).size() != k
)
613 } while( !regular
&& !ldpcG
.isConnected() );
615 // Convert to FactorGraph
616 vector
<Factor
> factors
;
617 for( size_t I
= 0; I
< K
; I
++ ) {
619 for( size_t _i
= 0; _i
< k
; _i
++ ) {
620 size_t i
= ldpcG
.nb2(I
)[_i
];
623 factors
.push_back( Factor( vs
, paritycheck
) );
627 // Generate noise vector
628 vector
<char> noisebits(N
,0);
630 for( size_t i
= 0; i
< N
; i
++ ) {
631 if( rnd_uniform() < noise
) {
636 cout
<< "# bitflips = " << bitflips
<< endl
;
638 // Simulate transmission of all-zero codeword
639 vector
<char> input(N
,0);
640 vector
<char> output(N
,0);
641 for( size_t i
= 0; i
< N
; i
++ )
642 output
[i
] = (input
[i
] + noisebits
[i
]) & 1;
645 for( size_t i
= 0; i
< N
; i
++ )
646 factors
.push_back( Factor(Var(i
,2), likelihood
+ output
[i
]*2) );
648 // Construct Factor Graph
649 fg
= FactorGraph( factors
);
651 throw "Invalid type";
653 // Output additional comments
654 cout
<< "# N = " << fg
.nrVars() << endl
;
655 cout
<< "# seed = " << seed
<< endl
;
657 // Output factor graph
659 } catch( const char *e
) {
660 /// Display error message
661 cerr
<< "Error: " << e
<< endl
;