|
&&&&&
"Bug" reported by Ondra Kamenik, confirmed by Julie Langou, fixed by Jim Demmel
and Guillaume Revy. See forum post 1783.
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
Fix problem in xTGEX2. The threshold value is too stringent and some matrices
are failing. Relax the threshold by a factor 2. More below.
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
See: http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=2&t=1783
dgges issue >> by OndraKamenik >> Tue Mar 16, 2010 6:01 am
All,
I have the following problem with dgges. For version 3.1.1 and sooner, I get a
reasonable result, for version 3.2 and 3.2.1 I get info=n+2.
I am separating eigenvalues in the unit circle from one outside the unique
circle.
The two D and E matrices have relatively well separated null spaces, the
minimum angle is acos(0.97). However, if i calculate condition numbers of
E-lambda*D of lambda=[-1:0.01:1], they are quite bad.
The matrices are attached with a small c++ program which calls dgges and sorts
eigenvalues for a better comparison. There is also a Makefile, which links with
different version of lapack.
I use the reference blas.
My question is if it is a bug in Lapack introduced between 3.2. and 3.1.1 or
the matrix is just very bad and in version 3.1.1 I was just lucky to get a
reasonable solution.
Many thanks for any help.
Ondra Kamenik
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
Much more conversation on mailing list and forum [ skipped ]
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
Date: Wed, 31 Mar 2010 19:53:11 -0600
From: James Demmel
This bug was introduced by changing the routine that computes Householder
transformations in the last release of LAPACK from dlarfg to the new dlarfp,
which makes the diagonal of the R factor in QR nonnegative. The bug was a
failure (INFO = N+3) in the dgges routine for computing selected deflating
subspaces of a matrix pencil A - lambda*B. I will describe the bug, a quick
fix, and implications for the floating point debugging project that Guillaume
and others of us are working on.
The bug occurred in dtgex2, in the swapping of eigenvalues (adjacent 1x1 and/or
2x2 blocks on the diagonal of the generalized Schur form) in order to compute
a selected subspace. The swapping involves QR decompositions of small (n <= 4)
matrices.
The code in dtgex2 performs two of its own internal correctness tests (a "weak"
one and a "strong" one) to see if the swapping has been performed stably. The
two tests compute a residual in slightly different ways. The weak test passed,
but the strong test (which can be commented out by setting the internal
parameter WANDS to be .false.) failed, leading to returning INFO = N+3.
However, it failed by exceeding the threshold THRESH only by a factor like 1.2
or less. THRESH is set to 10*macheps*dnorm, so if we changed the (probably
somewhat arbitrary) factor from 10 to 20, it would work. Or we could set WANDS
= false.
Since Bo's name is on this routine, his comments are particularly welcome.
Guillaume and I spent a while tracking this down, and will continue to find out
why dlarfp led to a (slightly!) larger residual than the old dlarfg.
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
Date: Thu, 1 Apr 2010 20:01:21 -0600
From: James Demmel <demmel@cs.berkeley.edu>
Actually, on reexamining the data, this example would be fixed by
changing 10 to 11,
but let's go with 20 to be on the safe side :) . I think it would be
most efficient
if Julie or you appropriately change the one line of code in DTGEX2:
THRESH = MAX( TEN*EPS*DNORM, SMLNUM )
All versions (S/D/C/Z) have an analogous line of code with the same
constant TEN that I would change.
Sorry, I've lost track of the Mathworks bug report on dlarfp. Can you
remind me?
Thanks,
Jim
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
|
