wolffd@0: function gdl = mk_gdl_graph(G, domains, node_sizes, kernels, varargin)
wolffd@0: % MK_GDL_GRAPH Make a GDL (generalized distributed law) graph
wolffd@0: % gdl = mk_gdl_graph(G, domains, node_sizes, kernels, ...)
wolffd@0: %
wolffd@0: % A GDL graph is like a moralized, but untriangulated, Bayes net:
wolffd@0: % each "node" represents a domain with a corresponding kernel function.
wolffd@0: % For details, see "The Generalized Distributive Law", Aji and McEliece,
wolffd@0: % IEEE Trans. Info. Theory, 46(2): 325--343, 2000
wolffd@0: % 
wolffd@0: % G(i,j) = 1 if there is an (undirected) edge between domains i,j
wolffd@0: %
wolffd@0: % domains{i} is the domain of node i
wolffd@0: %
wolffd@0: % node_sizes(i) is the number of values node i can take on,
wolffd@0: %   or the length of node i if i is a continuous-valued vector.
wolffd@0: % node_sizes(i) = 1 if i is a utility node.
wolffd@0: %
wolffd@0: % kernels is the list of kernel functions
wolffd@0: %
wolffd@0: % The list below gives optional arguments [default value in brackets].
wolffd@0: % 
wolffd@0: % equiv_class - equiv_class(i)=j  means factor node i gets its params from factors{j} [1:F]
wolffd@0: % discrete - the list of nodes which are discrete random variables [1:N]
wolffd@0: % chance   - the list of nodes which are random variables [1:N]
wolffd@0: % decision - the list of nodes which are decision nodes [ [] ]
wolffd@0: % utility  - the list of nodes which are utility nodes [ [] ]
wolffd@0: 
wolffd@0: 
wolffd@0: ns = node_sizes;
wolffd@0: N = length(domains);
wolffd@0: vars = [];
wolffd@0: for i=1:N
wolffd@0:   vars = myunion(vars, domains{i});
wolffd@0: end
wolffd@0: Nvars  = length(vars);
wolffd@0: 
wolffd@0: gdl.equiv_class = 1:length(kernels);
wolffd@0: gdl.chance_nodes = 1:Nvars;
wolffd@0: gdl.utility_nodes = [];
wolffd@0: gdl.decision_nodes = [];
wolffd@0: gdl.dnodes = 1:Nvars;
wolffd@0: 
wolffd@0: if nargin >= 5
wolffd@0:   args = varargin;
wolffd@0:   nargs = length(args);
wolffd@0:   for i=1:2:nargs
wolffd@0:     switch args{i},
wolffd@0:      case 'equiv_class', bnet.equiv_class = args{i+1}; 
wolffd@0:      case 'chance',      bnet.chance_nodes = args{i+1}; 
wolffd@0:      case 'utility',     bnet.utility_nodes = args{i+1}; 
wolffd@0:      case 'decision',    bnet.decision_nodes = args{i+1}; 
wolffd@0:      case 'discrete',    bnet.dnodes = args{i+1}; 
wolffd@0:      otherwise,  
wolffd@0:       error(['invalid argument name ' args{i}]);       
wolffd@0:     end
wolffd@0:   end
wolffd@0: end
wolffd@0: 
wolffd@0: 
wolffd@0: gdl.G = G;
wolffd@0: gdl.vars = vars;
wolffd@0: gdl.doms = domains;
wolffd@0: gdl.node_sizes = node_sizes;
wolffd@0: gdl.cnodes = mysetdiff(vars, gdl.dnodes);
wolffd@0: gdl.kernels = kernels;
wolffd@0: gdl.type = 'gdl';
wolffd@0: 
wolffd@0: % Compute a bit vector representation of the set of domains
wolffd@0: % dom_bitv(i,j) = 1 iff variable j occurs in domain i
wolffd@0: gdl.dom_bitv = zeros(N, length(vars));
wolffd@0: for i=1:N
wolffd@0:   gdl.dom_bitv(i, domains{i}) = 1;
wolffd@0: end
wolffd@0: 
wolffd@0: % compute the interesection of the domains on either side of each edge (separating set)
wolffd@0: gdl.sepset = cell(N, N);
wolffd@0: gdl.nbrs = cell(1,N);
wolffd@0: for i=1:N
wolffd@0:   nbrs = neighbors(G, i);
wolffd@0:   gdl.nbrs{i} = nbrs;
wolffd@0:   for j = nbrs(:)'
wolffd@0:     gdl.sepset{i,j} = myintersect(domains{i}, domains{j});
wolffd@0:   end
wolffd@0: end
wolffd@0: 
wolffd@0: