Problem with generating all binary vectors of constant Hamming weight

Hi, all.

I need to precompute a matrix containing all binary vectors with the length n and Hamming weight k, listed in a lexicographical order. It is known, that total number of such vectors is defined by appropriate binomial coefficient C(n,k).

Anyway, I wrote the following function, that do the job:

function vectors = constant_weight(n, k)

total = factorial(n)/(factorial(k)*factorial(n-k));
vectors = zeros(total, n);

next = [ones(1,k) zeros(1,n-k)];

num=0;
while 1

      num = num + 1;
      vectors(num,:) = next;

      if num == total
          break
      end

      for i = n-1:-1:1
          if (next(i) == 1) && (next(i+1) == 0)
              next(i) = 0;
              next(i+1) = 1;
              break;
          end
      end

      if i+2 < n
          next(i+2:n) = next(n:-1:i+2);
      end

end

end

Inside the while loop I just used Narayana's algoritm to find lexicographically next vector.

My fuction works pretty fast, for example, with n=36 and k=4. But already with n=37 and k=4 it dramatically slows down: from few seconds to hours.

Can anyone suggest an idea, why this slowdown happens? Is it my fault, or MATLAB feature?

Or maybe, someone could propose another solution of the original task, which will be faster?

Answer 1

[Alexander Shishkin]

It seems that I've found the solution. I just need to replace

total = factorial(n)/(factorial(k)*factorial(n-k));

with

total = round(factorial(n)/(factorial(k)*factorial(n-k)));

But any ideas about further performance improvement are welcome.

Answer 2

[Michal Kuba]

This solution is very similar to yours but it seems that the recursion is more efficient (but this is the question):

function [] = constant_weight_demo()

% ---------------------------------------------------------------------------
% Input:

% Vector length:
n = 9;

% Vector weight:
k = 5;

% ---------------------------------------------------------------------------
% Inicialization:

% First vector memory allocation:
x = zeros([1, n]);

% First vector memory inicialization:
x(1 : k) = 1;

% Last vector memory allocation:
y = zeros([1, n]);

% Last vector memory initialization:
y((n - k + 1) : n) = 1;

% ---------------------------------------------------------------------------
% Recursion:

while 1

     % ----------------------------------------------------------------------
     % Do it something with current vector x ...

     % ... for example printing the current vector to the command window:
     fprintf(1, '  x = %s\n', vect2str(x, 1));

     % Note: vect2str() - see below.
     % ----------------------------------------------------------------------

     % Termination condition:
     if isequal(x, y)

        break

     end

     % Next vector generation from current vector (see below):
     x = constant_weight(x);

     % ----------------------------------------------------------------------

end

% ----------------------------------------------------------------------------

end

% ----------------------------------------------------------------------------
% ----------------------------------------------------------------------------

function next = constant_weight(prev)

% ----------------------------------------------------------------------------

next = [];

if isempty(prev)

     return

end

n = numel(prev);

next = double(prev ~= 0);

k = sum(next(:));

if k == 0

     return

end

if k == n

     return

end

% ----------------------------------------------------------------------------

i = find((next(1 : n - 1) - next(2 : n)) == 1, 1, 'first');

if isempty(i)

     return

end

x = next(1);

next(i) = 0;

next(i + 1) = 1;

if x == 1

     return

end

if i < 3

     return

end

next(1 : (i - 1)) = next((i - 1) : -1 : 1);

% ----------------------------------------------------------------------------

end

% ----------------------------------------------------------------------------
% ----------------------------------------------------------------------------

function s = vect2str(x, n)

I = numel(x);

n = length(sprintf('%d', n));

if I == 0

     s = '[]';

     return

end

s = sprintf('[%*d', n, x(1));

for i = 2 : I

     s = sprintf('%s, %*d', s, n, x(i));

end

s = strcat(s, ']');

end

% ----------------------------------------------------------------------------
% ----------------------------------------------------------------------------