Array indexing question for vectorization
1 visualización (últimos 30 días)
Mostrar comentarios más antiguos
I have three arrays A, B, and IND, with A containing some values that I want to update, B containing the values that I want to use to update A, and IND containing the indeces to update in A. IND may reference the same index in A several times. For example, I may have:
A = ones(1, 10);
B = 1:10;
IND = [1 1 2 5 5 5 5 6 8 9];
and the update I want to do is multiplication. I can achieve this by doing
for ii = 1:10
A(IND(ii)) = A(IND(ii)) * B(ii);
end
which gives me my desired answer of
A =
2 3 1 1 840 8 1 9 10 1
Is there a way I can do this in a vectorized operation, avoiding the for loop (I'm potentially doing this kind of operation thousands of times on large arrays)? Doing
A(IND) = A(IND) .* B
results in
A =
2 3 1 1 7 8 1 9 10 1
Any tips greatly appreciated!
2 comentarios
Matt Fig
el 13 de Nov. de 2012
Are A and B always like that? Could they be more like:
A = randi(N,1,M);
B = randi(K,1,J);
or not?
Respuesta aceptada
Honglei Chen
el 13 de Nov. de 2012
Editada: Honglei Chen
el 13 de Nov. de 2012
You can try the following trick
c = unique(IND)
d = accumarray(IND',B',[],@prod)
A(c) = A(c).*d(c)'
6 comentarios
José-Luis
el 13 de Nov. de 2012
True, but it still is slower than a for loop:
Elapsed time is 0.000003 seconds. %for
Elapsed time is 0.000828 seconds. %unique+accumarray
Elapsed time is 0.000515 seconds. %accumarray
Más respuestas (1)
Jan
el 13 de Nov. de 2012
Editada: Jan
el 13 de Nov. de 2012
Is your A pre-allocated? If IND is sorted, running the loop backwards wil be faster:
for ii = 10:-1:1
A(IND(ii)) = A(IND(ii)) * B(ii);
end
If IND is not sorted, pre-allocate explicitly:
A = zeros(1, max(IND));
For more speed create a C-Mex function. If you are interested and have a compiler installed, I can post the few required lines of code.
2 comentarios
Jan
el 13 de Nov. de 2012
@Richard: A = zeros(...); ... A(k) = A(k) * ... leads to zeros. As in your original question your need ones().
Ver también
Categorías
Más información sobre Logical en Help Center y File Exchange.
Productos
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!