I need to adjust the following algorithm in Matlab, in particular I refer to the part for the construction of B
because it takes to much time when n
is large. The part for the construction of A
has been brilliantly suggested here Fast algorithm to get combinations in Matlab?
clear all
%% Construction of A: list of COMBINATIONS of the elements of the set {0,1} in n-1 places (organize
% 2 elements in n-1 places with repetitions and the order does NOT matter) repeated twice:
% once when the n-th element is 1, the other when the n-th element is 0.
n=3;
A=zeros(n,2*n);
A(:,1)=1;
for i=2:2:n*2-1
A(:,i-1)=circshift(A(:,i-1),[-1 0]);
A(:,i)=A(:,i-1);
A(end,i)=0;
A(:,i+1)=A(:,i);
end
A(:,end-1)=circshift(A(:,end-1),[-1 0]);
A=A';
groupindex=(linspace(1,size(A,1),size(A,1)))';
A=[A groupindex];
%% Construction of B: for each row of A -excluding the last column- (for example consider row 1 of A: 1 1 1) I want to
% (step 1) List all possible DISPOSITIONS of the elements of
% the set {0,1} in n-1 places
% (step 2) Combine row 1 with all
% possible dispositions listed in step 1 and
% get C=[1 1 1 1 1; 1 1 1 1 0; 1 1 1 0 1; 1 1 1 0 0]
% (step 3) Associate with the same number some rows of C according to this logic:
% [1 1 1 1 1]: the first two elements are ones; the third element is 1; the first two elements should be associated
% with the fourth and the fifth elements which are
% ones; so 1 1 1 1 can be represented as (1,1), (1,1) 1
% [1 1 1 1 0]: (1,1), (1,0), 1
% [1 1 1 0 1]: (1,0), (1,1), 1
% [1 1 1 0 0]: (1,0), (1,0), 1
% So row 2 and row 3 of C are equivalent
% if we exchange (1,0) and (0,1) and should
% be associated with the same number
%
% Consider now instead row 3 of A: 1 0 1
% (step 1) ...
% (step 2) C=[1 0 1 1 1; 1 0 1 1 0; 1 0 1 0 1; 1 0 1 0 0]
% (step 3) [1 0 1 1 1]: (1,1), (0,1) 1
% [1 0 1 1 0]: (1,1), (0,0) 1
% [1 0 1 0 1]: (1,0), (0,1) 1
% [1 0 1 0 0]: (1,0), (0,0) 1
% So there are no equivalent rows
%(step 1)
g=makeindex(n-1); %matrix 2^(n-1)x(n-1): list of dispositions of the elements of the set {0,1} in n-1 places
%(step 2)
repA=kron(A, ones(2^(n-1),1)); %repeat each row of A 2^(n-1) times
repg=repmat(g,size(A,1),1); %stack vertically the matrix g for size(A,1) times
B=[repA(:,1:n) repg repA(:,n+1)]; %matrix size(A,1)*2^(n-1) x (n+n-1+1)
%(step 3)
Xr=B(:,n);
m=sum(B(:,1:n-1),2);
lm=zeros(size(B,1),1);
for i=1:size(B,1)
count=0;
for j=1:n-1
if B(i,j)==1 && B(i,n+j)==1
lm(i)=count+1;
count=lm(i);
end
end
end
lf=zeros(size(B,1),1);
for i=1:size(B,1)
count=0;
for j=1:n-1
if B(i,j)==0 && B(i,n+j)==1
lf(i)=count+1;
count=lf(i);
end
end
end
classindex=zeros(size(B,1),1);
count=0;
for i=1:size(B,1)
old=classindex;
for j=1:size(B,1)
if j>=i && lm(i)==lm(j) && lf(i)==lf(j) && Xr(i)==Xr(j) && m(i)==m(j) && classindex(j)==0
classindex(j)=count+1;
end
end
diff=old-classindex;
if any(diff)==1;
count=count+1;
end
end
The function makeindex is
function [index]=makeindex(k) %
%
total=2^k; %
index=zeros(total,k); %
i=1; %
for i=1:k %
ii=1; %
cc=1; %
c=total/(2^i); %
while ii<=total %
if cc <=c %
index(ii,i)=1; %
cc=cc+1; %
ii=ii+1; %
else %
ii=ii+c; %
cc=1; %
end %
end %
end %
%
end