Algorithm optimization: do not use for loop

Hi, I have the following code which I'd like to optimize:

% Creation of all scenarios
AoA__ = [0 1 2 3];
Flight_Phase__ = [0 1];
AoA___ = allcomb(AoA__,AoA__,AoA__,Flight_Phase__)
s = size(AoA___);
n = s(1);
m = s(2);

% Calculation of sum and norm of each row

somme_ligne = sum(AoA___')';
norme_ligne = sqrt(sum(AoA___'.*AoA___'))';

somme_norme = [somme_ligne,norme_ligne];

unicite(1) = 1;

% Finding the redudancies

 for i = 2:n
    if ismember(somme_norme(i,:),somme_norme(1:i-1,:),'rows')
        unicite(i,1) = 0;
    else
        unicite(i,1) = 1;    
    end
end

nbr_scenarios = sum(unicite);

% Saving the corresponding indices 
for i = 1:n
    if unicite(i,1) == 1
        indices(i,1) = i;
    end
end

index = indices(indices~=0);

AoA_without_order = AoA___(index,:)

Let me explain the interest. I am creating cell vectors AoA and Flight Phase. Then, I encode them in underlying values (Simulink objective here). Then we have two vectors: AoA = [0 1 2 3] and Flight Phase = [0 1];

I create all combinations with order and repetitions with 'allcomb' function. But I only need combinations with répétitions and without order. For that, I compute the sum and norm of my allcomb matrix and gather them in one matrix 'somme_norme' with 'somme_norme(:,1) = sum of each row' and 'somme_norme(:,2) = norm of each row'. If two lines have the same norm and sum, normally, the combination is redundant, then unicite = 0. Else, unicite = 1. Afterwards, I extract the rows number of the first apparition of one combination and run it in my allcomb matrix.

With this example, my all combination matrix is 124*4. Without considering the order, we must have only 40 combinations (formula used: C(n+p-1,p)).

I would like to find an algorithm which does not use a for loop. With that value of n, the total runtim is more than 256h, which is too much. Any ideas for that? I am not used to vectorize my algorithm, so bad knowledge of Matlab functions.

Thanks

[Adam]

ismember is a function to avoid in a big for loop. It is not at all a fast function and includes a lot of checks that can usually be avoided.

[Antoine Guilleux]

Good to know, thanks. It was the only solution I found yet. Do you know what could replace this function?

[Jan]

How many columns does your input have? Which type is it?

[Antoine Guilleux]

The input 'somme_norme' has two columns of type double.

Answer 1

[Walter Roberson]

[~, ia] = unique(somme_norme, 'rows');
unicite = zeros(size(somme_norme,1),1);
unicite(ia) = 1;

[Antoine Guilleux]

Thanks for this answer. However, I've tested on my short example but I don't have the good result at the end. May be I was not clear enough? But I will try to adapt it anyway!

[Jan]

@Antoine: Please post some example data which shows, that this approach is not the solution. It looks perfect.

I dare to vote for it in the hope, that it works and there way a problem in your observation of "not a good result" only.

[Walter Roberson]

I created some random data and my code seemed to work fine.

[Antoine Guilleux]

Walter, yes finally it works. My use of it seems to be wrong. I'm gonna explain the meaning of this small part of my whole code.

The objective is to create combinations of the following vectors without order but with repetitions: [0 1 2 3], [0 1 2 3],[0 1 2 3],[0 1]. First, I compute all the combinations, the 'allcomb' function gives me a 4^3*2 = 124 by 4 matrix, that's correct. But for me, the combinations [1 0 0 1] and [0 1 0 1] are the same, and I only need one of them. That's why, for each row of my 'somme_norme' column, I compute the sum and the norm. Then, the values obtained are put in the 'somme_norme' matrix: first column: the sum of each row, second column, the norm of each row. Normally, two combinations with order must have the same sum and norm.

Then, if I sum the vector 'unicite', I must obtain the total number of combinations with répétitions and without order. With those inputs, I must obtain sum(unicite) = 40. However, your code returns 29 instead. But when my inputs vector are only [0 1 2 3], [0 1 2 3], [0 1 2 3] then it works and I have 20 combinations, and it is right.

Conclusion: your code works well, but I need to adapt it. For other cases, we can have 12 vectors as inputs with different number of elements. And here it does not work directly, I have to think about it.

But thanks anyway and sorry for my mistake!

[Walter Roberson]

Different approach than what you describe:

Generate your data. sort it along the second dimension: since you do not care what order it is in, you can leave it sorted. Now you can use unique() like I showed, without bothering with the sum or norm.

[Antoine Guilleux]

Yes that's exactly what I found out yesterday and it works well. However, the problem was coming from the inputs: you cannot do this with my short example:

allcomb([0 1 2 3],[0 1 2 3],[0 1 2 3],[0 1])

Here physically the first three vectors represent one parameter three times, and the last vector another one. I cannot mix them if I want to have all combinations without order. I think there's no other way than use it on the first parameter and then build the matrix "with your hands":

A = sort(allcomb([0 1 2 3],[0 1 2 3],[0 1 2 3]),2)
[~, ia] = unique(A, 'rows');
A_without_order = A(ia,:)
m = size(A_without_order,1)

% Building final matrix
Combinations(1:m,:) = A_without_order(:,:);
Combinations(m+1:2*m,:) = A_without_order(:,:);
Combinations(1:m,4) = 0;
Combinations(m+1:2*m,4) = 1;

This is may be not pretty, but quite efficient and without for loops, and I have my 40 combinations!

Answer 2

[Jan]

Creating the temporary arrays somme_norme(1:i-1,:) is cruel! Remember how much work this is for n=16384000. If I guess, that somme_norme has 10 columns and it is a double array, Matlab requests this number of bytes from the OS to create these arrays:

sum(1:16384000) * 10 * 8 Byte = 
10737418895360000 Byte = 
10'737 TByte

Of course this will take a long time, even without considering the work in ismember.

The same problem occurs for letting the output unicite grow iteratively: this requests sum(1:n)*8 bytes from the OS, because in each iteration a new vector is created and the old data are copied. Pre-allocation is the simple solution: Create the output before the loop with its final size.

It would be much more efficient to operate on columns instead of rows, because the element neighboring in columns are store side by side in the memory also.

S = somme_norme.';
unicite = ones(1, n);  % Pre-allocate!!!
for i = 2:n
  Si = S(:, i);
  for k = 1:i-1
    if all(Si == S(:, k))
      unicite(i) = 0;  % Stop the inner loop when the first match is found
      break;
    end
  end
end

This is a brute-force attack. For a unique data set, the inner loop runs again sum(1:16384000) times.

There are smarter solutions: sortrows can sort the array and then you have to check the neighboring rows only. But this is exactly what Walter's suggestion does in unique('rows'). So please post a small example to demonstrate, what you want and what Walter's code does instead.

[Antoine Guilleux]

Thanks for that. I obtain the same number of combinations than with Walter's, so definitely the problem is coming from somewhere else, my bad. I will put my whole code in the first post explaining what I want it to do. Thanks!

[Antoine Guilleux]

EDIT: I realized that computing the sum and norm of each row is not a good idea. Another example: [1 1 4] and [0 3 3] have the same sum = 6 and norm = sqrt(18). But there's no link between those combinations -> reason why I cannot find the right answer even with your codes... ><"