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>
Reply-To: <HIDDEN>
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)
X-Newsreader: MATLAB Central Newsreader 688530
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