This is a plain Octave (or Matlab) solution. The matrix Distances
is as in the question. The algorithm builds a 0-1 matrix a
in which each column encodes a set with sum of distances less than limit
(for example 10).
The matrix a
is initialized with identity, because all one-elements subsets are admissible (sum of distances is 0). Then each column is picked c = a(:,m);
and the possibility of adding another element is investigated (cand = c; cand(k) = 1;
means adding k-th element). Without loss of generality, it is enough to consider adding only the elements after the last current element of the set.
The product cand'*Distances*cand
is twice the sum of distances of the candidate set cand
. So, if it's less than twice the limit, the column is added: a = [a cand];
. At the end, the matrix a
is displayed in transposed form, for readability.
limit = 10;
n = length(Distances);
a = eye(n, n);
col1 = 1;
col2 = n;
while (col1 <= col2)
for m = col1:col2
c = a(:,m);
for k = max(find(c>0))+1:n
cand = c;
cand(k) = 1;
if cand'*Distances*cand < 2*limit
a = [a cand];
end
end
end
col1 = col2 + 1;
col2 = length(a);
end
disp(a')
Output:
1 0 0 0 0
0 1 0 0 0
0 0 1 0 0
0 0 0 1 0
0 0 0 0 1
1 1 0 0 0
1 0 1 0 0
1 0 0 1 0
0 1 1 0 0
0 1 0 1 0
0 0 0 1 1
With a 250 by 250 matrix, the performance will greatly depend on how large limit
is. If it is so large that all or most of 2^250 sets qualify, this is going to run out of memory. (2^250 is more than 10^75, so I don't think you'd ever want to see anywhere near that many sets).
And this is to have output in a more readable form:
for m = 1:length(a)
disp(find(a(:,m))')
end
Output:
1
2
3
4
5
1 2
1 3
1 4
2 3
2 4
4 5