platform/upstream/lapack.git
14 years ago&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
langou [Fri, 2 Apr 2010 04:55:58 +0000 (04:55 +0000)]
&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&

"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

&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&

14 years agoFix CMakeLists.txt only blas is exported, lapack needs to be exported as well. (From...
julie [Tue, 23 Mar 2010 19:57:29 +0000 (19:57 +0000)]
Fix CMakeLists.txt only blas is exported, lapack needs to be exported as well. (From Bill Hoffman)

14 years agoFix comments to generate wrapper correctly
julie [Tue, 23 Mar 2010 19:55:21 +0000 (19:55 +0000)]
Fix comments to generate wrapper correctly

14 years agoCorrect bug0050 :: An extra argument after matrix A lines 348-249 in STFSM.f causes...
julie [Thu, 18 Mar 2010 19:46:42 +0000 (19:46 +0000)]
Correct bug0050 :: An extra argument after matrix A lines 348-249 in STFSM.f causes a crash of SGEMM call sent from Alexander Kobotov (Intel) on Wed, 3 Feb 2010 to lapack@cs.utk.edu

14 years agobug0049 :: output parameter GIVPTR stays uninitialized in case of quick exit (n.eq...
julie [Thu, 18 Mar 2010 19:41:50 +0000 (19:41 +0000)]
bug0049 :: output parameter GIVPTR stays uninitialized in case of quick exit (n.eq.0) or if the rank-1 modifier is small enough in *laed8 sent by Alexander Kobotov (Intel) on Wed, 3 Feb 2010 to lapack@cs.utk.edu

14 years agoFix bug0048-Hanging could occur in *gebal if a NaN is in input matrix sent by Alexand...
julie [Thu, 18 Mar 2010 19:31:57 +0000 (19:31 +0000)]
Fix bug0048-Hanging could occur in *gebal if a NaN is in input matrix sent by Alexander Kobotov(Intel) to lapack on Feb 3rd 2010

14 years agoFix bug0046-Incorrectly documented RWORK workspace in ZGESDD sent by user Zbigniew...
julie [Thu, 18 Mar 2010 18:33:47 +0000 (18:33 +0000)]
Fix bug0046-Incorrectly documented RWORK workspace in ZGESDD sent by user Zbigniew on Forum (topic 1779)

14 years agoTake off comment to put xerbla back in BLAS lib
julie [Wed, 10 Mar 2010 17:48:14 +0000 (17:48 +0000)]
Take off comment to put xerbla back in BLAS lib

14 years agoFix buggy comments reported by Intel folowing Jim's recommendation
julie [Tue, 9 Mar 2010 21:49:41 +0000 (21:49 +0000)]
Fix buggy comments reported by Intel folowing Jim's recommendation

14 years agoShorten long lines
julie [Thu, 4 Mar 2010 20:02:50 +0000 (20:02 +0000)]
Shorten long lines

14 years agoAdd testversion and testieee as executables
julie [Tue, 23 Feb 2010 22:58:21 +0000 (22:58 +0000)]
Add testversion and testieee as executables

14 years agoFollowing Jim and Sven report, sgeqpf and shseqr are computational routines. slanv2...
julie [Thu, 18 Feb 2010 15:40:14 +0000 (15:40 +0000)]
Following Jim and Sven report,  sgeqpf and shseqr are computational routines.  slanv2 is an auxiliary routine

14 years agoSome minor comments modifications
julie [Wed, 17 Feb 2010 17:33:47 +0000 (17:33 +0000)]
Some minor comments modifications

14 years agoAdd make.inc for ifort
julie [Wed, 10 Feb 2010 21:13:55 +0000 (21:13 +0000)]
Add make.inc for ifort

14 years agoFollowing a complain from John Tellefson (Salina, KS), added the make of the
langou [Wed, 10 Feb 2010 16:09:38 +0000 (16:09 +0000)]
Following a complain from John Tellefson (Salina, KS), added the make of the
variant library to the testing of the variants if the variant library is not
present.

14 years agoBug 0046 fixed.
langou [Fri, 29 Jan 2010 15:22:58 +0000 (15:22 +0000)]
Bug 0046 fixed.

author: Vasile Sima
committer: Julien Langou

Corrected LAPACK routines dlagv2 and slagv2 based on Vasile Sima (National
Institute for Research & Development in Informatics, Bucharest, Romania)'s
email to lapack@cs.utk.edu on Monday 25 January 2010.

> Specifically, the variable WI, which is used in the line 268 (close to the
> end), is not initialized (with ZERO) in the cases "A can be deflated" and "B
> is singular" (i.e., the code segments in the lines 138-144, 148-156, and
> 158-167).  The corrected version is included in the attached archive.

14 years agoPolish some comments, etc.. for the C wrapper
julie [Thu, 28 Jan 2010 18:44:38 +0000 (18:44 +0000)]
Polish some comments, etc.. for the C wrapper

14 years agoRemove a trailing blank line at the end of dlarfp.f following Andy May (Cardiff
langou [Thu, 31 Dec 2009 16:46:32 +0000 (16:46 +0000)]
Remove a trailing blank line at the end of dlarfp.f following Andy May (Cardiff
University) 12/31/09's comment.

15 years agoStill working on there subroutines ... I have added a variable SAVEALPHA, I am
langou [Mon, 21 Dec 2009 22:37:40 +0000 (22:37 +0000)]
Still working on there subroutines ...  I have added a variable SAVEALPHA, I am
not sure whether it is needed or not.  What I am pretty sure of, is that the
code is correct with, the previous code without may be correct as well ...
Anyway, I like it better like this.

15 years agoDarned!!! I messed up with my commit r709!!!!
langou [Mon, 21 Dec 2009 21:25:30 +0000 (21:25 +0000)]
Darned!!! I messed up with my commit r709!!!!

This is
            " ( 1 - CONJG( TAU ) ) * ( ALPHA ) = ABS( ALPHA ) "
that needs to be true so, indeed, the previous code was correct and r709 is a
mistake. (r710 is good.) So this commit rollbacks r709 and modifies a little
r710.

15 years agoThis is the relevant change. This is the bug fix for the MathWorks/PatQuillen's
langou [Mon, 21 Dec 2009 21:16:18 +0000 (21:16 +0000)]
This is the relevant change. This is the bug fix for the MathWorks/PatQuillen's
bug in xLARFP. The fix seems to work her in Denver ...

15 years agoAfter three irrelevant commit ( r706, r707, r708), this is the first relevant
langou [Mon, 21 Dec 2009 21:11:58 +0000 (21:11 +0000)]
After three irrelevant commit ( r706, r707, r708), this is the first relevant
modifications. I believe that in the complex case, when NORMX = ZERO but ALPHI
is not ZERO (i.e. ALPHA has a complex imaginary part), then TAU needs to be
such that:

            " ( 1 - TAU ) * ( ALPHA ) = ABS( ALPHA ) "

Since we have
            ALPHR = REAL( ALPHA )
            ALPHI = AIMAG( ALPHA )
            XNORM = DLAPY2( ALPHR, ALPHI )
The way to do this is to set TAU with

            TAU = CMPLX( ONE - ALPHR / XNORM, ALPHI / XNORM )

as opposed to

            TAU = CMPLX( ONE - ALPHR / XNORM, -ALPHI / XNORM )

(Note: XNORM is used as temporary variable here)

15 years agochange SAFMIN with SMLNUM
langou [Mon, 21 Dec 2009 20:57:28 +0000 (20:57 +0000)]
change SAFMIN with SMLNUM
change RSAFMN with BIGNUM

15 years agoremove the comments on the constraints of TAU
langou [Mon, 21 Dec 2009 20:54:29 +0000 (20:54 +0000)]
remove the comments on the constraints of TAU
this comment is true for xLARFG, it is not true anymore for xLARFP

15 years agostart some modification of the xLARFP routines
langou [Mon, 21 Dec 2009 20:52:57 +0000 (20:52 +0000)]
start some modification of the xLARFP routines
modification 1: remove the comments starting with ! and start them with *

15 years agoCorrect a bug found by Ashutosh Mahajan from Mathematics and Computer Science Divisio...
langou [Wed, 16 Dec 2009 01:58:23 +0000 (01:58 +0000)]
Correct a bug found by Ashutosh Mahajan from Mathematics and Computer Science Division of Argonne National Labs.
Email sent on Thu, 10 Dec 2009 16:53:57 to lapack@cs.utk.edu.

I think we have a bug in dstemr.f for the case when N is 2. The following
block:

392:                   IF (CS.NE.ZERO) THEN
393:                      ISUPPZ(2*M-1) = 1
394:                      ISUPPZ(2*M-1) = 2
395:                   ELSE
396:                      ISUPPZ(2*M-1) = 1
397:                      ISUPPZ(2*M-1) = 1
398:                   END IF

should really be (note lines 394, 397):

392:                   IF (CS.NE.ZERO) THEN
393:                      ISUPPZ(2*M-1) = 1
394:                      ISUPPZ(2*M) = 2
395:                   ELSE
396:                      ISUPPZ(2*M-1) = 1
397:                      ISUPPZ(2*M) = 1
398:                   END IF

similarly the block:
415:                IF (SN.NE.ZERO) THEN
416:                   IF (CS.NE.ZERO) THEN
417:                      ISUPPZ(2*M-1) = 1
418:                      ISUPPZ(2*M-1) = 2
419:                   ELSE
420:                      ISUPPZ(2*M-1) = 1
421:                      ISUPPZ(2*M-1) = 1
422:                   END IF
423:                ELSE

should really be (note lines 418, 421):

415:                IF (SN.NE.ZERO) THEN
416:                   IF (CS.NE.ZERO) THEN
417:                      ISUPPZ(2*M-1) = 1
418:                      ISUPPZ(2*M) = 2
419:                   ELSE
420:                      ISUPPZ(2*M-1) = 1
421:                      ISUPPZ(2*M) = 1
422:                   END IF
423:                ELSE

15 years agoCorrect a bug found by Ashutosh Mahajan from Mathematics and Computer Science Divisio...
langou [Wed, 16 Dec 2009 01:51:40 +0000 (01:51 +0000)]
Correct a bug found by Ashutosh Mahajan from Mathematics and Computer Science Division of Argonne National Labs.
Email sent on Fri, 11 Dec 2009 14:51:07 to lapack@cs.utk.edu.

Subject: [Lapack] bug in dsyevr.f when N = 1

when dsyevr is called with N=1, then it returns without setting up ISUPPZ values.

331:       IF( N.EQ.1 ) THEN
...
342:          IF( WANTZ )
343:      $      Z( 1, 1 ) = ONE
344:          RETURN
345:       END IF

It should rather do

IF( WANTZ )
    Z( 1, 1 ) = ONE
    ISUPPZ(1) = 1
    ISUPPZ(2) = 1
END IF

15 years ago some more typos
langou [Sat, 14 Nov 2009 19:22:15 +0000 (19:22 +0000)]
 some more typos

15 years agotypo in comments
langou [Sat, 14 Nov 2009 19:20:08 +0000 (19:20 +0000)]
typo in comments

15 years agoFix bug0038 - dgelsd waas not returning iwork size
julie [Tue, 20 Oct 2009 21:26:21 +0000 (21:26 +0000)]
Fix bug0038 - dgelsd waas not returning iwork size

15 years agoBug report from Michael Chuvelev from Intel on DSYEVR.
langou [Thu, 15 Oct 2009 14:21:05 +0000 (14:21 +0000)]
Bug report from Michael Chuvelev from Intel on DSYEVR.

======================================================================================
Date: Wed, 14 Oct 2009 02:42:06 -0600
From: "Chuvelev, Michael" <michael.chuvelev@intel.com>
To: "lapack@cs.utk.edu" <lapack@cs.utk.edu>
Subject: [Lapack] VL,   VU are referenced in lapack-3.2.1 dsyevr even if range.ne.'V'

Hello,

dsyevr contains this code (lines 360-361):

      VLL = VL
      VUU = VU

which means VL, VU is accessed regardless of RANGE.

Whereas this issue is fixed in ssyevr, for instance (lines 364-367):

      IF (VALEIG) THEN
         VLL = VL
         VUU = VU
      END IF

Best regards,

Michael.
======================================================================================

15 years agoCorrect problem reported by Kevin Wadleigh to keep compilers from optimizing NAN6...
julie [Wed, 30 Sep 2009 15:02:44 +0000 (15:02 +0000)]
Correct problem reported by Kevin Wadleigh to keep compilers from optimizing NAN6 = NAN5*0.0

15 years agoDate: Tue, 29 Sep 2009 06:17:06 -0600
langou [Tue, 29 Sep 2009 15:17:49 +0000 (15:17 +0000)]
Date: Tue, 29 Sep 2009 06:17:06 -0600
From: "Chuvelev, Michael"
Subject: [Lapack] Bug in lapack-3.2.1/TESTING/EIG/zdrvgbx.f

 Hi,

please pay attention to the printing format 9998 misuse in zdrvgbx.f:

                              WRITE( NOUT, FMT = 9998 )'ZGBSVXX', FACT,
     $                             TRANS, N, KL, KU, IMAT, 7,
     $                             RESULT( 7 )
...

9998 FORMAT( ' *** In ZDRVGB, LAFB=', I5, ' is too small for N=', I5,
     $      ', KU=', I5, ', KL=', I5, /
     $      ' ==> Increase LAFB to at least ', I5 )

This would possibly lead to a segmentation fault in case of errors in zgbsvxx.

The same issue concerns other precisions too.

Best regards,

Michael.

15 years ago Fix bug0038 see out of range element access in DLASET when using DGELSS (forum topic...
julie [Fri, 25 Sep 2009 18:48:09 +0000 (18:48 +0000)]
 Fix bug0038 see out of range element access in DLASET when using DGELSS (forum topic 1582)

15 years agoCommit Julien's proposed change for bug0019
julie [Fri, 11 Sep 2009 21:04:20 +0000 (21:04 +0000)]
Commit Julien's proposed change for bug0019

15 years agoremove extra definition
julie [Fri, 11 Sep 2009 21:01:17 +0000 (21:01 +0000)]
remove extra definition

15 years agoFix whitespace comments detected from parser
julie [Fri, 11 Sep 2009 20:28:33 +0000 (20:28 +0000)]
Fix whitespace comments detected from parser

15 years agoFix bug0023: SROTMG and DROTMG uses deprecated Fortran ASSIGN statement and assigned...
julie [Thu, 10 Sep 2009 22:59:31 +0000 (22:59 +0000)]
Fix bug0023: SROTMG and DROTMG uses deprecated Fortran ASSIGN statement and assigned GOTO statement, actually fixed ROTM also

15 years ago Fix ifort flags just for Unix
julie [Wed, 12 Aug 2009 19:47:39 +0000 (19:47 +0000)]
 Fix ifort flags just for Unix

15 years ago Fix variables used
julie [Wed, 12 Aug 2009 19:42:14 +0000 (19:42 +0000)]
 Fix variables used

15 years ago Use try_compile to look for the timing routine
julie [Wed, 12 Aug 2009 19:36:19 +0000 (19:36 +0000)]
 Use try_compile to look for the timing routine

15 years agoAdd conf file for submitting machine to LAPACK Dashboard
julie [Wed, 12 Aug 2009 17:26:34 +0000 (17:26 +0000)]
Add conf file for submitting machine to LAPACK Dashboard

15 years agoBetter solution as we set NONE as default...WINDOWS build run great...
julie [Tue, 11 Aug 2009 23:06:34 +0000 (23:06 +0000)]
Better solution as we set NONE as default...WINDOWS build run great...

15 years agoNeed the NONE for Windows...
julie [Tue, 11 Aug 2009 22:31:48 +0000 (22:31 +0000)]
Need the NONE for Windows...

15 years agoAdd CMAKE support to LAPACK
julie [Tue, 11 Aug 2009 21:58:21 +0000 (21:58 +0000)]
Add CMAKE support to LAPACK

15 years agoFixed comments on factorization in refinement routines. Found when coding in C++
deaglanhalligan [Fri, 17 Jul 2009 09:01:45 +0000 (09:01 +0000)]
Fixed comments on factorization in refinement routines. Found when coding in C++

15 years agoIn refinement routines, err_bnds initialization if-else statements should be if state...
deaglanhalligan [Fri, 17 Jul 2009 08:27:04 +0000 (08:27 +0000)]
In refinement routines, err_bnds initialization if-else statements should be if statements, found when coding in C++

15 years agothru -> through
langou [Wed, 15 Jul 2009 13:36:27 +0000 (13:36 +0000)]
thru -> through

15 years agofixed the variants_testing Makefile problem
du [Thu, 4 Jun 2009 20:24:24 +0000 (20:24 +0000)]
fixed the variants_testing Makefile problem

15 years agoBug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack...
langou [Sun, 10 May 2009 21:01:41 +0000 (21:01 +0000)]
Bug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack@cs.utk.edu".
"c(he/sy)equb: no slamch in external lists and lsame, intrinsics not described"

15 years ago[[ I did not really want to commit the previous commit, this is the follow-up commit...
langou [Sun, 10 May 2009 20:59:05 +0000 (20:59 +0000)]
[[ I did not really want to commit the previous commit, this is the follow-up commit ... ]]

Bug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack@cs.utk.edu".

"(d/s)tgsen.f: iwork(1) always referenced: line 455: iwork( 1 ) = LIWMIN,
documentation says that if IJOB=0 it shouldn't, so NULL pointer causes a
sigfault."

(There is indeed the same problem for the WORK array and the IWORK array in the
complex routines.)

I have changed the header of the routines ctgsen.f, dtgsen.f, stgsen.f, and
ztgsen.f. A Fortran array needs to be of size at least 1. So IWORK is of size
at least 1. It was indeed written in the header of the routine:
"IWORK (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))"
So since IWORK is of size at least 1, there is no reason not to reference it.
When there is a workspace query with IJOB.EQ.0, WORK(1) and IWORK(1) are both
set to 1.

15 years agoBug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack...
langou [Sun, 10 May 2009 20:51:42 +0000 (20:51 +0000)]
Bug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack@cs.utk.edu".

"(d/s)tgsen.f: iwork(1) always referenced: line 455: iwork( 1 ) = LIWMIN,
documentation says that if IJOB=0 it shouldn't, so NULL pointer causes a
sigfault."

There is indeed the same problem for array work.

I have changed the header of the routines ctgsen.f, dtgsen.f, stgsen.f, and
ztgsen.f. No matter what an array needs to be of size at least 1. So
IWORK is of size at least 1. It was written in the header of the routine:
"IWORK   (workspace/output) INTEGER array, dimension (MAX(1,LIWORK))"
So since IWORK is of size at least 1, there is no reason not to reference it.
When there is a workspace query with IJOB.EQ.0, IWORK(1) is set to 1.

15 years agoBug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack...
langou [Sun, 10 May 2009 20:45:07 +0000 (20:45 +0000)]
Bug report sent by Alexander V. Kobotov (from Intel) on Mon, 6 Apr 2009 to "lapack@cs.utk.edu".
"(c/z)geesx: no lwork=-1 branch at all, info = -15 is returned while doing lquery."
I noticed that the S and D version where not exiting after WORKSPACE query as well.
Thanks Alexander. Bug corrected.

15 years agoPatch from Lawrence Mulholland (NAG) on Wednesday May 6th 2009.
langou [Sun, 10 May 2009 18:10:23 +0000 (18:10 +0000)]
Patch from Lawrence Mulholland (NAG) on Wednesday May 6th 2009.
See forum: http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=2&t=1408.
Thanks Lawrence.

15 years agoBug reported by Julie Langou on Th. 04/30/09 thanks to the Intel Fortran
langou [Mon, 4 May 2009 03:08:48 +0000 (03:08 +0000)]
Bug reported by Julie Langou on Th. 04/30/09 thanks to the Intel Fortran
compiler on Windows.  I messed up the workspace types needed by CPOT01 and
ZPOT01, they need to be REAL (or DOUBLE PRECISION) and not COMPLEX (or DOUBLE
COMPLEX). Everything works fine on my machine with these changes. Julien.

15 years agoBig commit before 3.2.1 release.
julie [Thu, 16 Apr 2009 18:10:16 +0000 (18:10 +0000)]
Big commit before 3.2.1 release.
Those are just cosmetic changes to update version number and various other minor change.

15 years agoSome debug write statements were present in slafts.f and dlafts.f => removed.
langou [Wed, 8 Apr 2009 21:16:17 +0000 (21:16 +0000)]
Some debug write statements were present in slafts.f and dlafts.f => removed.

15 years agoUpdated documentation for EPIR routines. Changed ERRS_{N,C} variable names. Other...
deaglanhalligan [Wed, 8 Apr 2009 00:05:18 +0000 (00:05 +0000)]
Updated documentation for EPIR routines. Changed ERRS_{N,C} variable names. Other cosmetic changes.

15 years agofix bug :: in DGESDD, workspace query gives a value smaller than the minimal value...
julie [Tue, 17 Mar 2009 19:12:22 +0000 (19:12 +0000)]
fix bug :: in DGESDD, workspace query gives a value smaller than the minimal value given in the header to run the routine

    o reported by  Guy Bencteux on Sat Dec 06 2008
    o see forum topic 846 : https://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=2&t=846

remove a min(M,N) in the minimal workspace formula when JOBZ='O' and  JOBZ = 'S' or 'A'.

Note Julie on minimal value:
Case: M > N

IF( M > N*5.0D0  ) THEN
PATH 1 -> FOR JOBZ=N : MINWORK = 7*N + N
PATH 2 -> FOR JOBZ=O: MINWORK = 5*N*N + 7*N
PATH 3 -> FOR JOBZ=S: MINWORK = 4*N*N + 7*N
PATH 4 -> FOR JOBZ=A:  MINWORK = 4*N*N + 7*N

IF M>N (but not too large)
PATH 5 -> FOR JOBZ=N : MINWORK = MAX(7*N,M) + 3*N
PATH 6 -> FOR JOBZ=O: MINWORK = 3*N + MAX (M,4*N*N + 4*N)
PATH 7 -> FOR JOBZ=S : MINWORK = 3*N + MAX (M,3*N*N + 4*N)
PATH 8 -> FOR JOBZ=A :  MINWORK = 3*N + MAX (M,3*N*N + 4*N)

GENERAL FORMULA FROM THE TWO CASES WHEN M>N
FOR N LWORK >= 3*N + max (M,7*N)
FOR O LWORK >= 3*N + max (M,5*N*N + 4*N) [instead of 3*N*N + max (M,5*N*N + 4*N) ]
FOR S ET A LWORK >= 3*N + max(M,4*N*N + 4*N) [instead of 3*N*N + max(M,4*N*N + 4*N) ]

15 years agoAdd xSYEQUB bibliographic reference citing Livne & Golub.
jason [Sat, 14 Mar 2009 20:17:33 +0000 (20:17 +0000)]
Add xSYEQUB bibliographic reference citing Livne & Golub.

Reference: Livne, O.E. and Golub, G.H., "Scaling by Binormalization",
Numerical Algorithms, vol. 35, no. 1, pp. 97-120, January 2004.
DOI 10.1023/B:NUMA.0000016606.32820.69
Tech report version: http://ruready.utah.edu/archive/papers/bin.pdf

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
15 years agoFix ZLARFP and CLARFP optimizations when the vector is zero.
jason [Sat, 14 Mar 2009 20:06:19 +0000 (20:06 +0000)]
Fix ZLARFP and CLARFP optimizations when the vector is zero.

The x == zero branch ignored complex alphas.  The code still
functioned, but it scaled the entire zero vector.  Now that I think of
it, I should scan upwards for the scaling, too.  That will be a
separate enhancement patch; this is just the bug fix.

This patch also fixes an accidental precision shortening in ZLARFP,
which used CMPLX in the dead branch.

Reported by Igor Zhuravlov on the LAPACK web goo thingy, currently at
  http://icl.cs.utk.edu/lapack-forum/viewtopic.php?t=924
Also lead to finding a related (but different) error in the LAWN and
(accepted, still in editing) SISC paper.

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
15 years agoThe description of LDA in the header was number 8 in the header, while
langou [Wed, 11 Mar 2009 01:54:40 +0000 (01:54 +0000)]
The description of LDA in the header was number 8 in the header, while
LDA is the 4th argument in the interface.
=>
Rearrange the header accordingly.

15 years agoChange routines to xLARFP that were missed in the first round.
jason [Fri, 6 Mar 2009 18:22:09 +0000 (18:22 +0000)]
Change routines to xLARFP that were missed in the first round.

Thanks to Michael Chuvelev <michael.chuvelev@intel.com> for noting
that the following routines still called xLARFG rather than xLARFP:
  cgeqpf.f   slaqp2.f   claqp2.f  claqps.f   clatrz.f   ctzrqf.f

Bug reported in message
<A43BE3D8E247E440AF3477B3CDCB0A694A132639@irsmsx504.ger.corp.intel.com> .

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
15 years agoReturn SCALE = ONE in xTRSYL for zero-dimension matrices.
jason [Tue, 3 Mar 2009 17:16:37 +0000 (17:16 +0000)]
Return SCALE = ONE in xTRSYL for zero-dimension matrices.

I doubt if LAPACK is consistent about zero-dim matrices, but it's
worth a shot.  Plenty of discussion about their semantics on the
GNU Octave development lists.

Reported by Vasile Sima <vsima@ici.ro> in message
<1575.193.230.3.167.1235984292.squirrel@webmail.ici.ro>.

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
15 years agoOnly touch Q in {d,s}tgsen.f if WANTQ is true.
jason [Tue, 3 Mar 2009 17:12:04 +0000 (17:12 +0000)]
Only touch Q in {d,s}tgsen.f if WANTQ is true.

Reported by Daniel Waggoner <dwaggoner@frbatlanta.org> in message
<49AA2617.4060900@frbatlanta.org>.

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
15 years agoRestructured loops in {c,z}la_yyrcond_{x,c} and {s,c,d,z}la_yyamv to do less comparisons.
deaglanhalligan [Fri, 13 Feb 2009 23:58:46 +0000 (23:58 +0000)]
Restructured loops in {c,z}la_yyrcond_{x,c} and {s,c,d,z}la_yyamv to do less comparisons.

15 years agoShould call {c,z}la_heamv, not {c,z}la_syamv, from {c,z}la_porfsx_extended.
deaglanhalligan [Wed, 11 Feb 2009 09:20:58 +0000 (09:20 +0000)]
Should call {c,z}la_heamv, not {c,z}la_syamv, from {c,z}la_porfsx_extended.

15 years agoModified {c,z}syequb to use 2*N WORK rather than 3*N.
deaglanhalligan [Wed, 11 Feb 2009 07:13:27 +0000 (07:13 +0000)]
Modified {c,z}syequb to use 2*N WORK rather than 3*N.

15 years agoAdded Jasons patches to bring work and rwork down to 2*N in the complex refinement...
deaglanhalligan [Wed, 11 Feb 2009 06:06:54 +0000 (06:06 +0000)]
Added Jasons patches to bring work and rwork down to 2*N in the complex refinement routines.

15 years agocgesvxx and zgesvxx need at least 3*N WORK space after changes from last bug fix
deaglanhalligan [Tue, 10 Feb 2009 12:17:14 +0000 (12:17 +0000)]
cgesvxx and zgesvxx need at least 3*N WORK space after changes from last bug fix

15 years agoFixed WORK/RWORK bug in complex iterative refinement routines.
deaglanhalligan [Tue, 10 Feb 2009 03:14:40 +0000 (03:14 +0000)]
Fixed WORK/RWORK bug in complex iterative refinement routines.

15 years ago(no commit message)
langou [Mon, 9 Feb 2009 22:08:05 +0000 (22:08 +0000)]

15 years ago(no commit message)
langou [Mon, 9 Feb 2009 21:37:20 +0000 (21:37 +0000)]

15 years ago(no commit message)
langou [Sun, 1 Feb 2009 19:51:40 +0000 (19:51 +0000)]

15 years agoFix scaling in the transpose case of extra-precise refinement.
jason [Wed, 28 Jan 2009 20:34:18 +0000 (20:34 +0000)]
Fix scaling in the transpose case of extra-precise refinement.

The row-scaling factors need passed in the transpose case.  Apparently,
none of our release tests include ill-scaled matrices.

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
15 years ago==============================================================================
langou [Wed, 28 Jan 2009 16:11:32 +0000 (16:11 +0000)]
==============================================================================
Patch from Christof Voemel, ETH Zurich.
==============================================================================
I would like to make a small change  in dlarrd/slarrd/_larrd.f.base so that the
code deals better with certain matrices from Godunov.

Could you please replace in those files the lines, after label 40,
            GL = GL - FUDGE*SPDIAM*EPS*IN - FUDGE*PIVMIN
            GU = GU + FUDGE*SPDIAM*EPS*IN + FUDGE*PIVMIN

by
            GL = GL - FUDGE*TNORM*EPS*IN - FUDGE*PIVMIN
            GU = GU + FUDGE*TNORM*EPS*IN + FUDGE*PIVMIN

(Exchange SPDIAM to TNORM).
==============================================================================

15 years agoSet ithresh to 10, not 100. Removed la_linrx_max_n_errs.
deaglanhalligan [Mon, 26 Jan 2009 03:51:22 +0000 (03:51 +0000)]
Set ithresh to 10, not 100. Removed la_linrx_max_n_errs.

15 years agoFixed typos and variable types in testing files.
deaglanhalligan [Thu, 22 Jan 2009 23:28:13 +0000 (23:28 +0000)]
Fixed typos and variable types in testing files.

15 years agoAdded the rest of the extra precise iterative refinement testing.
deaglanhalligan [Thu, 22 Jan 2009 10:40:59 +0000 (10:40 +0000)]
Added the rest of the extra precise iterative refinement testing.

15 years agoFix out of bound access in LIN TESTING with XBLAS found with -fcheck-bounds with...
julie [Fri, 16 Jan 2009 22:27:20 +0000 (22:27 +0000)]
Fix out of bound access in LIN TESTING with XBLAS found with -fcheck-bounds with gfortran.
"Error in DLA_GBAMV.f line 226: forrtl: severe (408): fort: (2): Subscript #1 of the array AB has value 2 which is greater than the upper bound of 1"

Reported in http://icl.cs.utk.edu/trac/lapack-dev/ticket/45

15 years agoFollowing a out-of-bound complaint by gfortran
julie [Mon, 12 Jan 2009 22:16:18 +0000 (22:16 +0000)]
Following a out-of-bound complaint by gfortran

Modify size in declaration of DX and DY.
It was set to 1, and 5 or N in the comments !!!!!
Put * in the declaration and N in the comments

This routine may need to be double-checked.

15 years ago============================================================================
langou [Mon, 12 Jan 2009 05:17:21 +0000 (05:17 +0000)]
============================================================================
Bug in RFP routines xTFSM for N=1 independently reported by Jason and Julie.
============================================================================

See Jason's email DEC/27/2009: "Similar problem in xTFSM".
Jason used gfortran with -fbounds-check and got:
> At line 328 of file stfsm.f
> Fortran runtime error: Array reference out of bounds for array 'b', upper bound of dimension 1 exceeded (1 > 0)

Julie also observed the bug with MS Visual Studio.
See ticket http://icl.cs.utk.edu/trac/lapack-dev/ticket/46

The problem is the following. Sometimes in the RFP routines, we
call routines like for example
   DGEMM ( M, N, K, ... A( I, J), ...)
where M or N is 0 and A( I, J ) is out-of-bound.  The rationale was since M or
N is 0, DGEMM simply exits without doing anything, therefore the out-of-bound
A(I,J) is not an issue. Well, this is not Fortran correct.

In this commit, I have fix the problem for xTFSM. This is a fix. This is not
rocket-science and one can certainly do something more elegant ... This bug is
blocking Julie from porting on Windows so the matter was urgent.

Now RFP + ( gfortan -fbounds-check ) works fine.

TODO: find a cleaner way to fix this,

TODO: Unfortunately I believe there are other bugs like this one in other RFP
routines, always when N=1. See for example xPFTRF. I do not understand why the
TESTING with gfortran -fbounds-check do not trigger an error though ...

15 years agoFix some minors errors found while porting under Windows
julie [Thu, 8 Jan 2009 19:58:17 +0000 (19:58 +0000)]
Fix some minors errors found while porting under Windows

15 years agoLast round of modifications to the comments for the generation of the manpages
julie [Fri, 2 Jan 2009 21:57:02 +0000 (21:57 +0000)]
Last round of modifications to the comments for the generation of the manpages

15 years agoMerged revisions 609-614 via svnmerge from
jason [Tue, 30 Dec 2008 21:27:12 +0000 (21:27 +0000)]
Merged revisions 609-614 via svnmerge from
https://jason@icl.cs.utk.edu/svn/lapack-dev/lapack/branches/SC08-release

........
  r609 | julie | 2008-12-16 17:17:52 -0500 (Tue, 16 Dec 2008) | 1 line

  Polish routines to fit the LAPACK framework and allow manpages generation
........
  r610 | langou | 2008-12-19 12:12:38 -0500 (Fri, 19 Dec 2008) | 30 lines

  bug reported on the forum
  https://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=2&t=854

  the complete thread is available at
  http://groups.google.com/group/comp.lang.fortran/browse_thread/thread/635192e11beadb93#

  Tobias Burnus also sent us an email:

  > Hello,
  >
  > this was reported at
  > http://groups.google.com/group/comp.lang.fortran/browse_thread/thread/635192e11beadb93#
  >
  > The problem is the line 47:
  >
  > 47:       IF( M.EQ.0 .OR. A(M, 1).NE.ZERO .OR. A(M, N).NE.ZERO ) THEN
  >
  > If M == 0 the one accesses A(0,1) which is invalid as the lower bound is 1
  > and not 0.
  >
  > Note: Contrary to C there is no left-to-right evaluation of expressions in
  > Fortran; the order is left to the compiler. One might assume that a smart
  > compiler does not evaluate "A(M,1)" if "M==0", however, there is nothing in
  > the standard guarantees this.
  >
  > If bounds checks are turned on (see post at the URL above), gfortran aborts
  > with an out-of-bounds error.
........
  r611 | julie | 2008-12-19 15:00:58 -0500 (Fri, 19 Dec 2008) | 5 lines

  Modify the formatting of the comments.
  Replace Note and Notes section by Further Details
  This allow the manpages to be generated corectly.
........
  r612 | julie | 2008-12-19 16:29:21 -0500 (Fri, 19 Dec 2008) | 3 lines

  Reformat the xblas routines comments to be able to generate the manpages

........
  r613 | julie | 2008-12-19 16:30:31 -0500 (Fri, 19 Dec 2008) | 1 line

  Update version number
........
  r614 | jason | 2008-12-27 09:44:45 -0500 (Sat, 27 Dec 2008) | 13 lines

  Fix non-short-circuited tests in ILAxL{C,R}.

  Fortran doesn't short-circuit logical operators, so the check that the leading
  dimension /= 0 may not prevent indexing into a 0-length array.

  Reported by "hes selex" in
    http://groups.google.com/group/comp.lang.fortran/browse_thread/thread/635192e11beadb93
  and forwarded to the LAPACK maintainers by Tobias Burnus <burnus@net-b.de>.

  Chalk up more bugs found by gfortran's diagnostics!

Signed-off-by: Jason Riedy <ejr@cs.berkeley.edu>
Cc: Tobias Burnus <burnus@net-b.de>
........

16 years ago(no commit message)
julie [Tue, 16 Dec 2008 17:06:58 +0000 (17:06 +0000)]

16 years agoXBLAS: Added gbmv2 (including test codes).
yozo [Tue, 11 Nov 2008 19:57:27 +0000 (19:57 +0000)]
XBLAS: Added gbmv2 (including test codes).

This routines computes the matrix product:

    y  <-  alpha * op(A) * (x_head + x_tail) + beta * y

where A is a general banded matrix.

16 years agoXBLAS: Added hemv2 (including test codes).
yozo [Tue, 11 Nov 2008 19:56:54 +0000 (19:56 +0000)]
XBLAS: Added hemv2 (including test codes).

This routines computes the matrix product:

  y  <-  alpha * A * (x_head + x_tail) + beta * y

where A is a complex Hermitian matrix.

16 years agoXBLAS: Added symv2 (including tests).
yozo [Tue, 11 Nov 2008 19:56:25 +0000 (19:56 +0000)]
XBLAS: Added symv2 (including tests).

symv2 performs

  y  <-  alpha * A * (x_head + x_tail) + beta * y

where A is a symmetric matrix.

16 years agoXBLAS: Fix increment adjustment in do_test_axpby.
yozo [Tue, 11 Nov 2008 19:56:07 +0000 (19:56 +0000)]
XBLAS: Fix increment adjustment in do_test_axpby.

We need to double the increment for complex cases.

16 years agoXBLAS: Fix bug in testgen_BLAS_[cz]dot2.
yozo [Tue, 11 Nov 2008 19:55:57 +0000 (19:55 +0000)]
XBLAS: Fix bug in testgen_BLAS_[cz]dot2.

Wrong parameter was being passed to gen_y_to_cancel.
Should pass n_fix2+1 instead of k+1, since k and n_fix2
differ by factor of two for complex case.

16 years agoXBLAS: Fix bug in testing/test-dot2/testgen_BLAS_[sd]dot2.
yozo [Tue, 11 Nov 2008 19:55:50 +0000 (19:55 +0000)]
XBLAS: Fix bug in testing/test-dot2/testgen_BLAS_[sd]dot2.

Scaling factor to compute x_tail was being rounded to integer,
and thus was producing all zeros for x_tail (and hence didn't
test whether *mv2 routines were looking at x_tail at all).

16 years agoXBLAS: Removed unused variable x_vec in testing/test-symv/do_test_symv.c.
yozo [Tue, 11 Nov 2008 19:55:41 +0000 (19:55 +0000)]
XBLAS: Removed unused variable x_vec in testing/test-symv/do_test_symv.c.

16 years agoXBLAS: Make argument error return codes more specific in gbmv.
yozo [Tue, 11 Nov 2008 19:55:32 +0000 (19:55 +0000)]
XBLAS: Make argument error return codes more specific in gbmv.

Should return -k if the k-th argument is in error, instead of just 0.

16 years agoXBLAS: Fix lda check in gbmv.
yozo [Tue, 11 Nov 2008 19:55:11 +0000 (19:55 +0000)]
XBLAS: Fix lda check in gbmv.

For banded storage lda can be smaller than n or m.

16 years agoXBLAS: Fix lda check in gemv.
yozo [Tue, 11 Nov 2008 19:54:53 +0000 (19:54 +0000)]
XBLAS: Fix lda check in gemv.

We need to raise error if lda < m for column major format,
and lda < n for row major format.  Previously it was checking
for lda < leny.

Also test lda in do_test_gemv more carefully.

16 years agoMove LAPACK trunk into position.
jason [Tue, 28 Oct 2008 01:38:50 +0000 (01:38 +0000)]
Move LAPACK trunk into position.