```Path: news.mathworks.com!not-for-mail
From: <HIDDEN>
Newsgroups: comp.soft-sys.matlab
Subject: Re: find vector r<=p (p: a given vector)
Date: Mon, 24 Nov 2008 21:55:04 +0000 (UTC)
Organization: Battelle Energy Alliance (INL)
Lines: 35
Message-ID: <ggf7ro\$fal\$1@fred.mathworks.com>
References: <gf2ppf\$q6l\$1@fred.mathworks.com> <ggcs55\$aiv\$1@fred.mathworks.com>
NNTP-Posting-Host: webapp-03-blr.mathworks.com
Content-Type: text/plain; charset="ISO-8859-1"
Content-Transfer-Encoding: 8bit
X-Trace: fred.mathworks.com 1227563704 15701 172.30.248.38 (24 Nov 2008 21:55:04 GMT)
X-Complaints-To: news@mathworks.com
NNTP-Posting-Date: Mon, 24 Nov 2008 21:55:04 +0000 (UTC)
Xref: news.mathworks.com comp.soft-sys.matlab:502953

"Roger Stafford" <ellieandrogerxyzzy@mindspring.com.invalid> wrote in message
>   Just for variety's sake here is a for-loop method, Oriole.  It should be much faster
> than your double loop.  Bruno's vectorized solution is certainly more elegant-looking > and could well be even faster.
>
>  n = length(p);
>  d = max(floor(p),0);
>  R = (0:d(n))';
>  s = d(n)+1;
>  for k = n-1:-1:1
>   R = repmat(R,d(k)+1,1);
>   t = repmat(0:d(k),s,1);
>   s = s*(d(k)+1);
>   R = [t(:),R];
>  end
>
> Roger Stafford

I liked Roger's approach, except that the R is not pre-allocated and is built-up through repmatting.  I fixed this.   Using timeit this version shows a 2.5 to 3.2 times gain in speed depending on p.  Here is my for-loop approach:

sz = prod(p+1); % The final array size.
[idx,idx] = sort(p,'descend'); % Largest to smallest.
R = zeros(sz,length(p)); % Preallocation.
y = zeros(sz,1); % Preallocation.

for kk = idx % Fill R in this order.
div = sz;
sz = sz/(p(kk)+1);
y(1:sz:end) = 1; % Poke in ones.
y(div+1:div:end) = -p(kk); %  Poke in drops.
R(:,kk) = cumsum(y) - 1;
y(:) = 0;  % Reset for the next pass.
end

```