How to avoid for loops when generating index arrays?

6 次查看(过去 30 天)
I often find myself coding nested for loops to generate vectors of integer indices. For example:
n = 4;
i = 1;
for L = 0:n
for M = -L:L
l(i) = L;
m(i) = M;
i = i+1;
end
end
All I need are the vectors "l" and "m". I can preallocate to save some speed, but my real problem is having to use the for loops as sometimes the index vectors I need to create have many more nested for loops whose (note that the inner loop index depends on the outer loop index).
Is there a simple way to avoid using loops to generate index vectors like these?

采纳的回答

Matt J
Matt J 2013-4-30
Here's another method, less memory consuming than NDGRID
mm=sparse(-n:n);
ll=sparse(0:n);
map=bsxfun(@le,abs(mm.'),ll);
idx=nonzeros(bsxfun(@times,map,1:length(ll) ));
l=full(ll(idx));
idx=nonzeros(bsxfun(@times,map,(1:length(mm)).')) ;
m=full(mm(idx));

更多回答(5 个)

Roger Stafford
Roger Stafford 2013-4-30
编辑:Matt J 2013-4-30
For your particular problem you can do this:
M = (0:n*(n+2))';
L = floor(sqrt(M));
M = M-L.*(L+1);
(I've used uppercase letters, 'L' and 'M', in place of your lowercase 'l' and 'm'.)
As with Matt Kindig, I am not sure this will be any faster than your for-loops. Time it with a large value for n and see.

Matt J
Matt J 2013-4-30
编辑:Matt J 2013-4-30
Here's a way to do it using NDGRID. It's not apriori obvious whether for loops would or would not be faster. It depends what you plan to reuse.
[mg,lg]=ndgrid(-n:n,0:n);
idx=abs(mg)<=lg;
l=lg(idx).',
m=mg(idx).',
  6 个评论
Oliver
Oliver 2013-5-1
The reason that I care is that I have only given you a simplified example. In practice I have things like:
for N = 0:2:30
for L = 0:N
for M = -L:L
for P = 0:M
for Q = -P:P
...
end
end
end
end
end
So making the intermediate arrays in this case would require a 5-dimensional array that takes up a lot of space.
Matt J
Matt J 2013-5-1
编辑:Matt J 2013-5-2
I'm starting to think Sean's advice about sticking with for-loops is the best one. There can definitely be ways to cut down on the loop nesting (see my newest Answer based on cell arrays), but the required form would depend on the body of the original set of for-loops.

请先登录,再进行评论。


Sean de Wolski
Sean de Wolski 2013-4-30
编辑:Sean de Wolski 2013-4-30
doc meshgrid
doc ndgrid %?
:)
And of course, depending on your application, two nested for-loops or bsxfun() might be better.
  2 个评论
Oliver
Oliver 2013-4-30
I can't use ndgrid or meshgrid because the inner loop depends on the outer loop. If i had something like
i = 1;
for L = 3:6
for M = 1:5
l(i) = L;
m(i) = M;
i = i+1;
end
end
then I could use:
[l,m] = meshgrid(3:6,1:5);
But, this won't work with the dependency.
Sean de Wolski
Sean de Wolski 2013-4-30
Just use the for-loops, they'll be the fastest by far. If you want to disguise it, write a function that takes L and M and returns l and m.
cellfun and arrayfun are slow and converting between cells and numeric types is slow. The above with preallocation will be pretty quick.

请先登录,再进行评论。


Matt Kindig
Matt Kindig 2013-4-30
It's kind of hack-y, but it gives the same output as your original posting:
n=4;
l = cell2mat(arrayfun(@(x) x*ones(1,2*x+1), 0:n, 'uni', false));
m= cell2mat( arrayfun(@(x) (-x:1:x), 0:n, 'uni', false));
Keep in mind that this may very well be slower than for-loops--I haven't done any timing comparisons.

Matt J
Matt J 2013-5-1
编辑:Matt J 2013-5-1
Here's a way to eliminate one nested loop
l=cell(1,n+1);
m=l;
for L=0:n
i=L+1;
m{i}=-L:L;
l{i}=m{i};
l{i}(:)=L;
end
l=[l{:}],
m=[m{:}],
  2 个评论
Sean de Wolski
Sean de Wolski 2013-5-1
I'd be surprised if this is faster due to the cell array conversions. I guess one of us will have to run a timing test.
Matt J
Matt J 2013-5-1
编辑:Matt J 2013-5-2
For n=1000 I get this,
Original Approach:
Elapsed time is 0.093130 seconds.
Cell-Based Approach
Elapsed time is 0.027393 seconds.
I think the vectorization inherent in
m{i}=-L:L;
l{i}=m{i};
l{i}(:)=L;
trumps the overhead from the cell conversion.

请先登录,再进行评论。

类别

Help CenterFile Exchange 中查找有关 Matrix Indexing 的更多信息

产品

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by