Re: [Lapackpp-devel] LaLinearSolve not totally exact

[Christian S. <sti...@tu...>]

Am Mittwoch, 11. August 2004 02:09 schrieb Jacob \(Jack\) Gryn:
> I believe I solved the problem.
>
> The in the linslv.cc function, Xtmp should be a copy of B; however, at no
> point is an actual copy made.  Xtmp is left at 0's before making the dgels
> call.
>
> I replaced line 219 with the line Xtmp=B; and it seems to fix the problem.
>
> The output result appears to be correct, at least for an AX=A style
> problem, where X should return the identity matrix.

That's good :-))) . I think I really didn't work much with the QRSolve 
functions, so it might well be the case that you just discovered and fixed a 
bug. However, I thought in line 219 this assignment with LaIndex on the left 
and on the right should actually do some copying? Whatever. If you say it 
doesn't, then it probably doesn't and needs to be fixed. If you replace that 
line with a direct operator=, doesn't it run into some problems with the 
matrix dimension? The comment there seems to imply this. But of course that 
comment may well be wrong. If you see that it really is fixed by your 
proposal, then just go ahead and commit it.

> Although it actually returns I + some small values in the order of 10e-17;
> whereas CLAPACK returns 0's in these places for the same problem.  Is this
> only something that displays with higher precision, or is there actually a
> 10e-17 error intruduced somewhere?

I believe the 1e-17 is the actual epsilon, the machine precision. I see them 
more than once at places where analytically a zero (like, zero) should be 
found. But then, in real-world signals there is noise anyway, so these small 
numbers can really safely considered to be zero. Or in other words, I thought 
this is just the normal effect of finite-precision arithmetic.

Christian