Prev: OuterPosition in absolute units and plotyy
Next: Error Embedded Matlab
From: Sukanta on 12 Jul 2010 22:19
How can I speed up the following code?

N = 256; x = randn(N,N,N);
tic;
y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N));
toc
From: James Tursa on 12 Jul 2010 23:13
"Sukanta " <sukanta.basu(a)ttu.edu> wrote in message <i1gieq$dog$1(a)fred.mathworks.com>...
> How can I speed up the following code?
>
> N = 256; x = randn(N,N,N);
> tic;
> y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N));
> toc

With a mex routine. That will avoid the unnecessary overhead of the array slice copies.

Here is an example. Might be able to optimize it even more for better memory access, but this shows the technique. Caution: there is no argument checking. The input must be a 3D double array (any size as long as the 3rd dimension is > 1).

To compile it, name the file averagen.c and do this:

>> mex averagen.c

Then you can test it like this:

>> mex averagen.c
>> N = 256; x = randn(N,N,N);
>> tic;y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N)); toc
Elapsed time is 1.230002 seconds.
>> tic;z=averagen(x);toc
Elapsed time is 0.335785 seconds.
>> isequal(y,z)
ans =
1


James Tursa

-----------------------------------------------------------------------------------------------------
// file averagen.c

#include "mex.h"
void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[])
{
mwSize *dims;
mwSize dimz[3];
mwSize m, n, p, mn, mnp, i;
double *a, *b, *c;

dims = mxGetDimensions(prhs[0]);
m = dims[0];
n = dims[1];
p = dims[2];
dimz[0] = m;
dimz[1] = n;
dimz[2] = p - 1;
mn = m*n;
mnp = mn * (p - 1);
plhs[0] = mxCreateNumericArray(3, dimz, mxDOUBLE_CLASS, mxREAL);
c = mxGetPr(plhs[0]);
a = mxGetPr(prhs[0]);
b = a + mn;
for( i=0; i<mnp; i++ ) {
c[i] = 0.5 * ( a[i] + b[i] );
}
}
From: Cris Luengo on 13 Jul 2010 04:40
"James Tursa" <aclassyguy_with_a_k_not_a_c(a)hotmail.com> wrote in message <i1glk0$rbs$1(a)fred.mathworks.com>...
> "Sukanta " <sukanta.basu(a)ttu.edu> wrote in message <i1gieq$dog$1(a)fred.mathworks.com>...
> > How can I speed up the following code?
> >
> > N = 256; x = randn(N,N,N);
> > tic;
> > y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N));
> > toc
>
> >> mex averagen.c
> >> N = 256; x = randn(N,N,N);
> >> tic;y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N)); toc
> Elapsed time is 1.230002 seconds.
> >> tic;z=averagen(x);toc
> Elapsed time is 0.335785 seconds.
> >> isequal(y,z)
> ans =
> 1

If you don't want to go the MEX-route, you can get a less impressive speed improvement using the convolution operation:

>> tic,y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N));toc
Elapsed time is 0.604676 seconds.
>> f = cat(3,0.5,0.5);
>> tic,z = convn(x,f,'valid');toc
Elapsed time is 0.285232 seconds.
>> isequal(y,z)
ans =
1

I'm frustrated looking for the right function, because MATLAB has a built-in function to compute the difference of subsequent elements, which is /almost/ the same as what you want to do:

>> tic,z = diff(x,[],3);toc
Elapsed time is 0.146744 seconds.

Almost the same, but not the same. It feels like it should be an easy problem! Oh, well...

Cheers,
Cris.
From: Oleg Komarov on 13 Jul 2010 05:07
"Sukanta " <sukanta.basu(a)ttu.edu> wrote in message <i1gieq$dog$1(a)fred.mathworks.com>...
> How can I speed up the following code?
>
> N = 256; x = randn(N,N,N);
> tic;
> y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N));
> toc

w/o mex:

N = 256; x = randn(N,N,N);

% Elapsed time is 0.546098 seconds.
tic; y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N)); toc

% Elapsed time is 0.301945 seconds.
tic; y1 = reshape((reshape(x(:,:,1:N-1),1,[]) + reshape(x(:,:,2:N),1,[]))*.5,256,256,[]); toc

isequal(y, y1) % 1

Oleg
From: Cris Luengo on 13 Jul 2010 05:43
"Oleg Komarov" <oleg.komarovRemove.this(a)hotmail.it> wrote in message <i1habo$dfn$1(a)fred.mathworks.com>...
> N = 256; x = randn(N,N,N);
>
> % Elapsed time is 0.546098 seconds.
> tic; y = 0.5*( x(:,:,1:N-1) + x(:,:,2:N)); toc
>
> % Elapsed time is 0.301945 seconds.
> tic; y1 = reshape((reshape(x(:,:,1:N-1),1,[]) + reshape(x(:,:,2:N),1,[]))*.5,256,256,[]); toc

Oleg,
I'm baffled as to why the reshaped arrays compute faster. I tried it, it is faster, I just don't see why. Please enlighten me!
Cheers,
Cris.
 |  Next  |  Last
Pages: 1 2
Prev: OuterPosition in absolute units and plotyy
Next: Error Embedded Matlab