wolffd@0: function Gs = mk_all_dags(N, order)
wolffd@0: % MK_ALL_DAGS generate all DAGs on N variables
wolffd@0: % G = mk_all_dags(N)
wolffd@0: %
wolffd@0: % G = mk_all_dags(N, order) only generates DAGs in which node i has parents from 
wolffd@0: % nodes in order(1:i-1). Default: order=[] (no constraints).
wolffd@0: %
wolffd@0: % G{i} is the i'th dag
wolffd@0: %
wolffd@0: % Note: the number of DAGs is super-exponential in N, so don't call this with N > 4.
wolffd@0: 
wolffd@0: if nargin < 2, order = []; end
wolffd@0: 
wolffd@0: use_file = 0;
wolffd@0: 
wolffd@0: global BNT_HOME
wolffd@0: fname = sprintf('%s/DAGS%d.mat', BNT_HOME, N);
wolffd@0: if use_file & exist(fname, 'file')
wolffd@0:   S = load(fname, '-mat');
wolffd@0:   fprintf('loading %s\n', fname);
wolffd@0:   Gs = S.Gs;
wolffd@0:   return;
wolffd@0: end
wolffd@0: 
wolffd@0: m = 2^(N*N);
wolffd@0: ind = ind2subv(2*ones(1,N^2), 1:m);
wolffd@0: Gs = {};
wolffd@0: j = 1;
wolffd@0: directed = 1;
wolffd@0: for i=1:m
wolffd@0:   dag = reshape(ind(i,:)-1, N, N);
wolffd@0:   if acyclic(dag, directed)
wolffd@0:     out_of_order = 0;
wolffd@0:     if ~isempty(order)
wolffd@0:       for k=1:N-1
wolffd@0: 	if any(dag(order(k+1:end), k))
wolffd@0: 	  out_of_order = 1;
wolffd@0: 	  break;
wolffd@0: 	end
wolffd@0:       end
wolffd@0:     end
wolffd@0:     if ~out_of_order
wolffd@0:       Gs{j} = dag;
wolffd@0:       j = j + 1;
wolffd@0:     end
wolffd@0:   end
wolffd@0: end
wolffd@0: 
wolffd@0: if use_file
wolffd@0:   disp(['mk_all_dags: saving to ' fname '!']);
wolffd@0:   save(fname, 'Gs');
wolffd@0: end