How to vectorize this code to eliminate nested For loops

DaveF (view profile)

on 18 Jan 2013

Would like to know how this code be vectorized:

```X=zeros(nn);
for i=1:nn
for j=1:nn
X(i,j)=x(i)-x(j);
end
end
```

Products

No products are associated with this question.

Sean de Wolski (view profile)

on 18 Jan 2013

```X = bsxfun(@minus,x(1:nn).',x(1:nn));
```

Cedric Wannaz

Cedric Wannaz (view profile)

on 18 Jan 2013
` bsxfun(@minus,x(1:nn),x(1:nn).')`

Jan Simon (view profile)

on 18 Jan 2013
Edited by Jan Simon

Jan Simon (view profile)

on 19 Jan 2013

This cannot compete with the BSXFUN approach, but it is 32% faster than the original loop:

```nn = numel(x);
X  = zeros(nn);
for i = 1:nn
for j = i+1:nn
a      = x(i) - x(j);
X(i,j) = a;
X(j,i) = -a;
end
end
```

Using a temporary variable for x(i), which is updated in the outer loop only, is slightly slower.

Cedric Wannaz

Cedric Wannaz (view profile)

on 19 Jan 2013

Well then this one would be more efficient:

``` X  = zeros(nn);
for ii = 1:nn
X(:,ii) = x - x(ii);
end```

It almost competes with bsxfun for nn=1e4.

Jan Simon

Jan Simon (view profile)

on 19 Jan 2013

@Cedric: When you post this as an answer, you get my vote.

Your version does not exploit the symmetry as the other approaches. But this modification is much slower:

```X=zeros(nn);
for ii = 1:nn
a       = x(ii:nn) - x(ii);
X(ii:nn, ii) = a;
X(ii, ii:nn) = -a;
end
```

Obviously the overhead for explicit indexing is too high.

Cedric Wannaz

Cedric Wannaz (view profile)

on 19 Jan 2013

@Jan : thank you! Knowing that I would get your vote is enough, so I'll keep that as a comment ;-)

I just profiled your last version and it is very interesting indeed, because I would never have expected such an overhead. It might even lead me to review my code in a few projects, so thank you for the comment!

Jan Simon (view profile)

on 19 Jan 2013

And a C-Mex, which has the same speed as BSXFUN for a [1 x 1e3] vector under Matlab R2009a/64, Win7, MSVC 2008:

```#include "mex.h"
```
```void mexFunction(int nlhs, mxArray *plhs[], int nrhs, const mxArray *prhs[]) {
double *X, *Y;
mwSize i, j, n;
```
```     n = mxGetNumberOfElements(prhs[0]);
X = mxGetPr(prhs[0]);```
```     plhs[0] = mxCreateDoubleMatrix(n, n, mxREAL);
Y       = mxGetPr(plhs[0]);```
```     for (i = 0; i < n; i++) {
for (j = 0; j < n; j++) {
*Y++ = X[i] - X[j];
}
}
}```

Azzi Abdelmalek (view profile)

on 18 Jan 2013
Edited by Azzi Abdelmalek

Azzi Abdelmalek (view profile)

on 18 Jan 2013

```x=[1 3 4 10 20 30];
X=cell2mat(arrayfun(@(y) y-x,x','un',0))
```

%or

```X=bsxfun(@minus,repmat(x',1,numel(x)),x)
```

Cedric Wannaz

Cedric Wannaz (view profile)

on 18 Jan 2013

Or

` [X, Xt] = meshgrid(x, x) ;  Xt-X`

More efficient than cell2mat/arrayfun and less efficient than bsxfun. So the latter (bsxfun) should be the winner ;) But Sean's version without repmat beats them all.

Sean de Wolski

Sean de Wolski (view profile)

on 18 Jan 2013

Well the for-loops will smoke arrayfun/cell2mat like a Lambo would smoke bicycle...

bsxfun will be faster than the for-loops some of the time; meshgrid will only be faster if you need to reuse the gridded matrices.

Matt J (view profile)

on 19 Jan 2013

``` y=exp(x(1:nn));
X=log( y * (1./y.')  );```

Join the 15-year community celebration.

Play games and win prizes!

MATLAB and Simulink resources for Arduino, LEGO, and Raspberry Pi