help-octave
[Top][All Lists]

## Re: Computing distance matrix robustly

 From: Kevin Gross Subject: Re: Computing distance matrix robustly Date: Mon, 3 Oct 2005 09:17:16 -0400

```søn, 02 10 2005 kl. 18:50 -0400, skrev Tom Holroyd:

```
```On Sun, 2 Oct 2005, [utf-8] S?ren Hauberg wrote:

```
```Sorry about replying to my own post, but Tom Holroyd, contacted me
off-list with a solution. Simply replace
sX2 = sum(X.^2, 2);
with
sX2 = sumsq(X, 2);
in the vectorized version, and then things work.

```
```
I'm actually surprised that it worked.  But perhaps (in addition
to being faster) it is true that the FPU's registers are actually
wider than a double?  In that case, since X.^2 has to make
a temp array, it forces a loss of precision.  I wouldn't expect
1e-7 of roundoff anyway, though.  Unless the arrays in question
are floats?

```
```I must say I was quite surprised as well. The vectorized code was:

sX2 = sum(X.^2, 2); # use sumsq instead
XX = X*X';
D = sqrt( repmat(sX2, 1, n) + repmat(sX2', m, 1) - 2*XX );

That is, I also take the sqrt. So the round off error, is only about
1e-14, before sqrt.

/Søren

```
```Dr. Tom Holroyd
"A man of genius makes no mistakes. His errors are volitional and
are the portals of discovery." -- James Joyce

-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------

```
```

-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------

-------------------------------------------------------------
Octave is freely available under the terms of the GNU GPL.

Octave's home on the web:  http://www.octave.org
How to fund new projects:  http://www.octave.org/funding.html
Subscription information:  http://www.octave.org/archive.html
-------------------------------------------------------------

```