[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
Re: Computing distance matrix robustly
From: |
Søren Hauberg |
Subject: |
Re: Computing distance matrix robustly |
Date: |
Mon, 03 Oct 2005 09:17:20 +0200 |
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
-------------------------------------------------------------