julie [Wed, 6 May 2015 04:30:24 +0000 (04:30 +0000)]
Correct typo in variable name - reported by Elena (Oracle)
julie [Wed, 6 May 2015 04:21:45 +0000 (04:21 +0000)]
Propagate fix to bug111 to lapacke - minimum work space - Issue reported by Andrew James on LAPACK mailing list on April 30th 2015
langou [Tue, 5 May 2015 01:00:18 +0000 (01:00 +0000)]
Bug fix from Osni Marques, Beresford Parlett and Jim Demmel.
This fixes bugs 032 and 056.
Edits in cstein.f, dstein.f, sstein.f and zstein.f
From Osni on Monday, May 4:
As we discussed in our last conference call, I am attaching a new version of
_STEIN, which fixes bugs 032 and 056. (The bugs were reported for DSTEIN, but
I have propagated the fix to the other versions.) Also, in the process of
testing the fix for those two bugs with more difficult cases, we stumbled upon
a matrix that led SSTEIN to return NaNs in the eigenvectors ... We fixed that
too and, in summary, these are the changes in _STEIN:
1) The assignment after the GO TO 60 needs to be GPIND = J1 (instead of GPIND =
B1).
2) In 'Normalize and scale the righthand side vector Pb' (after GO TO 100 or GO
TO 120) we have replaced _ASUM with I_AMAX.
julie [Thu, 23 Apr 2015 03:25:41 +0000 (03:25 +0000)]
Correct info testing for SYSVXX routines
in ?SYSVXX
INFO = -10 needs to be changed to INFO = -11
The code in this section of the files is about to check and make sure all S
elements are positive. And S is argument 11 not 10.
in ?ERRVXX
Add info testing for input parameter 10 and 11
julie [Wed, 1 Apr 2015 14:56:02 +0000 (14:56 +0000)]
Problem reported by Elena (Oracle)
Some arguments return the wrong Info number
For example LDX is not the 14th argument, but 15th.
Correct corresponding testing routines
langou [Wed, 25 Mar 2015 22:45:19 +0000 (22:45 +0000)]
**** Correct bug 125 ****
Change the maximum number of QR iterations (ITMAX) for xLAHQR from ITMAX=30 to
ITMAX = 30 * MAX( 10, NH ).
Bug reported by Yoshihide Okimi on July 16, 2014
https://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=13&t=4566
and, independently, by Marco Caliari (Università di Verona) on February 11,
2015.
Involved in bug fix: Yoshihide Okimi, Marco Caliari (Università di Verona),
Meiyue Shao (LBL), Julien Langou, and Daniel Kressner (EPFL).
** Bug **
xLAHQR v3.5 does not converge on some matrices while xLAHQR v2 converges
** Note **
xLAHQR v2 is available through ARPACK so xLAHQR v2 is still in "wide" use
** Description **
xLAQR0 v3.5 and xLAQR4 v3.5 have ITMAX set to:
ITMAX = MAX( 30, 2*KEXSH )*MAX( 10, ( IHI-ILO+1 ) )
xLAHQR v2 has ITMAX set to:
ITMAX = 30 * NH
xLAHQR v3.5 has ITMAX set to:
ITMAX = 30
Yoshihide Okimi reports that, in his application, (1) xLAHQR v3.5 does not
converge, (2) xLAHQR v2 converges, (3) changing ITMAX to 30*NH in xLAHQR v3.5
enables convergence.
Marco Caliari reports that, in his application, using xLAHQR v2.0 converges,
while xLAHQR v3.5 does not converge.
** Bug fix **
Set ITMAX in xLAHQR to
ITMAX = 30 * MAX( 10, NH )
Email from Meiyue: I feel that setting ITMAX proportional to NH makes sense.
Perhaps we can change ITMAX to ITMAX*MAX( 10, NH ) so that it is also safe for
ill-conditioned tiny matrices.
langou [Wed, 25 Mar 2015 17:28:33 +0000 (17:28 +0000)]
**** Fix bug 124 ****
**** (I also fixed the documentation of the six subroutines ****
Bug reported by Victor_K on the forum.
See forum http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=13&t=4643
Bug fix by Julien Langou. Involved as well Nick Higham, Sven Hammarling and Craig Lucas.
** Bug **
xPSTF2 and xPSTRF returns NaNs when the matrix in input has N negative diagonal entries.
** Note **
A matrix with N negative diagonal entries is a highly non valid input, but, yet no reason to create some NaNs, better to have a cleaner solution.
** Description **
xPSTF2 and xPSTRF first scan the diagonal of the matrix to find the greatest element on the diagonal.
During this pass on the diagonal, xPSTF2/xPSTRF checks as well for ZERO or NAN on the diagonal and aborts if it finds a ZERO or a NAN on the diagonal.
*
* Compute stopping value
*
PVT = 1
AJJ = A( PVT, PVT )
DO I = 2, N
IF( A( I, I ).GT.AJJ ) THEN
PVT = I
AJJ = A( PVT, PVT )
END IF
END DO
IF( AJJ.EQ.ZERO.OR.DISNAN( AJJ ) ) THEN
RANK = 0
INFO = 1
GO TO 170
END IF
Since the input matrix is filled with a negative diagonal (and no ZERO and no NAN), xPSTF2/xPSTRF finds a maximum element (which is the largest negative elements on the diagonal), and then this (negative) element is used to compute DSTOP, the stopping value. But then we get a negative DSTOP,
IF( TOL.LT.ZERO ) THEN
DSTOP = N * DLAMCH( 'Epsilon' ) * AJJ
ELSE
DSTOP = TOL
END IF
and then I stop even thinking on what would happen next. It seems that the code continues and does produce NAN in the first column. Why not. But I think we should catch this case and try to handle it in a different way. Maybe.
** Bug fix **
Replace
IF( AJJ.EQ.ZERO.OR.DISNAN( AJJ ) ) THEN
With
IF( AJJ.LE.ZERO.OR.DISNAN( AJJ ) ) THEN
So that if there is any nonpositive element or any NAN on the diagonal of A, we return right away with RANK=0 and INFO=1.
** Email from Nick Higham **
[ Julien wrote ] So that if there is any nonpositive element or any NAN on the diagonal of A, we return right away with RANK=0 and INFO=1?
[ Nick wrote ]
I think so.
This is justified by the documentation, which says
INFO is INTEGER
If INFO = -K, the K-th argument had an illegal value,
= 0: algorithm completed successfully, and
> 0: the matrix A is either rank deficient with computed rank
as returned in RANK, or is indefinite.
We are in the indefinite case if any diagonal element is negative.
Actually, that is not quite true: we are in the "not positive definite"
case. We can only say "indefinite" if there are both positive *and* negative
diagonal entries! However, if that is going to be said then it would be
necessary to check that all other LAPACK routines are similarly precise in
their wording, so it's simplest, and not too misleading,
to stay with "indefinite".
julie [Wed, 18 Mar 2015 05:09:27 +0000 (05:09 +0000)]
Correcting Bug111 reported by Kevin Squire on LAPACK Mailing - August 15th 2013
From Kevin: We were following the recommendation for RWORK size calculation from
the ZGESDD header here ( http://netlib.org/lapack/complex16/zgesdd.f), and were
finding that when INFO=N, the size was insufficient, and generally led to a segmentation fault.
The current definition says that "If JOBZ = N, LRWORK >= 5*min(M,N)", but we found
that the minimum size needed to be at least 7*min(M,N). This was through experimentation
and comparison with code from Numpy and Octave, so you may wish to verify this through the code.
JL: bug confirmed and reproduced - correction in comments
If JOBZ = N, LRWORK >= 7*min(M,N)
julie [Tue, 17 Mar 2015 14:53:16 +0000 (14:53 +0000)]
Correct issue reported by Elena (Oracle)
Add SIGMA and EPS in a list of arguments and describe them in the "Arguments" section
julie [Fri, 6 Mar 2015 04:31:28 +0000 (04:31 +0000)]
Correction reported by Elena (Oracle)
LAPACK 3.5.0 [cz]hgeqz.f compute eigenvalues differently depending on whether
eigenvectors are requested or not. Since the computation is different, eigenvalues
are not bit-by-bit identical whether eigenvectors are requested or not.
This fix enables to have [cz]hgeqz.f behaving more of the same whether eigenvectors
are requested or not.
Note: The fact that eigenvalues are bit-by-bit identical whether eigenvectors are
requested or not is checked in the LAPACK test suite.
langou [Mon, 2 Mar 2015 16:04:35 +0000 (16:04 +0000)]
Correct two ". GT." to " .GT.". Bug reported by email ny Willem Vermin "Minor
issue in dorbdb.f and sorbdb.f" on Friday, February 27, 2015 at 4:03 AM.
philippe.theveny [Wed, 25 Feb 2015 23:40:47 +0000 (23:40 +0000)]
Fix wrong update when the right Schur vector is not computed.
julie [Wed, 25 Feb 2015 03:55:26 +0000 (03:55 +0000)]
Add back cggglm to Makefile...build goes fine
philippe.theveny [Tue, 24 Feb 2015 23:50:54 +0000 (23:50 +0000)]
Add xGGHD3: blocked Hessenberg reduction, code from Daniel Kressner.
Add xGGES3 and xGGEV3: computation of the Schur form, the Schur vectors, and
the generalized eigenvalues using the blocked Hessenberg reduction.
philippe.theveny [Tue, 24 Feb 2015 19:26:59 +0000 (19:26 +0000)]
Fix typo in comments.
philippe.theveny [Mon, 23 Feb 2015 18:06:06 +0000 (18:06 +0000)]
Fix typo in comments.
philippe.theveny [Fri, 6 Feb 2015 20:42:31 +0000 (20:42 +0000)]
Fix typo.
philippe.theveny [Fri, 6 Feb 2015 20:31:47 +0000 (20:31 +0000)]
Systematically return with INFO=1 when the root finder failed.
This corresponds to the documentation of [d/s]lasdx routines and fixes bug121
reported Justin Si.
julie [Fri, 30 Jan 2015 20:03:11 +0000 (20:03 +0000)]
Ed Anderson Bug Fix - Jan 26 2015
related to 1522
Correct call to [SD]chkhs to add new arguments.
julie [Fri, 30 Jan 2015 19:50:52 +0000 (19:50 +0000)]
Ed Anderson Improvement - Jan 26 2015
STGEX2 can call Level 3 BLAS and I did this in LAPACK3E.
I replaced STGEX2 with the LAPACK3E version (when I say STGEX2, I mean DTGEX2 too)
and also fixed an incorrect leading dimension in DTGEX2.
julie [Fri, 30 Jan 2015 19:24:26 +0000 (19:24 +0000)]
Ed Anderson Bug Fix - Jan 26 2015
Also in TESTING/EIG, I fixed two bugs in s/d chkhs; the diff of schkhs is:
790c790
< CALL SLACPY( ' ', N, N, U, LDU, UZ, LDU )
---
> CALL SLACPY( ' ', N, N, U, LDU, UZ, LDA )
831c831
< $ ABS( WI1( J )-WI3( J ) ) )
---
> $ ABS( WR1( J )-WR3( J ) ) )
In further testing, I was able to produce a test case where Test 8 fails comparing
the eigenvalues from xHSEQR(‘E’,’N’) with the eigenvalues from xHSEQR(‘S’,’V’).
The eigenvalues were the same, but in a different order. It appears that xHSEQR
does some rearrangement of the Schur matrix and so we should not expect the
eigenvalues to be in the same order, although I am trying to confirm this
with the authors. In the meantime, I have modified xCHKEE and xCHKHS to
compute the eigenvalues from xHSEQR(‘S’,’N’) in a separate array W2 (or
WR2 and WI2) so that test 8 can compare the eigenvalues from xHSEQR(‘S’,’N’)
with the eigenvalues from xHSEQR(‘S’,’V’).
julie [Fri, 30 Jan 2015 07:29:56 +0000 (07:29 +0000)]
Ed Anderson Bug Fix - Jan 26 2015
In TESTING/LIN/schkrq.f and schklq.f (only those, not the other precisions), the loop to look for test failures runs from 1 to NTESTS, it should be 1 to NT.
julie [Fri, 30 Jan 2015 07:26:22 +0000 (07:26 +0000)]
Ed Anderson Bug Fix - Jan 26 2015
In TESTING/LIN, the test routines xchkqr kept track of the number of tests performed in the variable NT, which could be 2, 6, or 7. Recently another test was added at the end and the loop that checks for failing tests was changed to run from 1 to NT to 1 to NTESTS, including tests that were never performed. As a result, with THRESH = 0., the number of tests failing exceeds the number of tests performed (for example, with n = 1 as the only N value and THRESH = 0.0 you get “604 out of 514 tests failed to pass the threshold”.) The right way to fix this would be to renumber the tests so that the new test is test #3 and the loop can run from 1 to NT again, but I just changed the count of tests performed from NT to NTESTS. For consistency between schkqr and schkrq and other test routines in the QR family, it would have been better for the new tests to go in a separate test module.
julie [Fri, 30 Jan 2015 07:22:01 +0000 (07:22 +0000)]
Ed Anderson Bug Fix - Jan 26 2015 - I fixed an old bug in the count of the number of tests in xdrvge and xdrvgb
julie [Sat, 10 Jan 2015 05:57:14 +0000 (05:57 +0000)]
Applied Jakub's CBLAS patch sent on Jan 9th
Netlib CBLAS can be put on top of BLAS with 64-bit integers.
It can be accomplished using the flag "WeirdNEC" :)
It basically redefines int to long, but at the same time, it had multiple bugs.
langou [Mon, 1 Dec 2014 05:07:52 +0000 (05:07 +0000)]
Apply patch from Pavel Holoborodko.
See http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=13&t=4624
Please thanks to
Pavel Holoborodko. Multiprecision Computing Toolbox for MATLAB. Advanpix LLC.
Below is a copy-paste from Pavel's post on the forum.
ZLARFT:
ZGEMM is not listed in EXTERNAL section (but extensively used in the code).
DLADIV:
This file has three different routines - DLADIV, DLADIV1, DLADIV2.
DLADIV2 has "if" in lower case:
if( BR.NE.ZERO ) THEN
DGGBAL/BAK, ZGGBAL/BAK:
There are few places in these routines where DOUBLE PRECISION variables are implicitly converted to INTEGER.
For example:
IR = LSCALE( I ) + SIGN( HALF, LSCALE( I ) )
JC=RSCALE(I)+SIGN(HALF,RSCALE(I))
K=RSCALE(I)
K=LSCALE(I)
R/LSCALE(I) is double precision, IR, JC and K are integers. For the sake of
correct semantic and F77 rules we have to use INT to show conversion
explicitly. Also this would make conversion to other languages much easier.
Interestingly, only those routines are missing the INT.
julie [Sat, 15 Nov 2014 05:06:16 +0000 (05:06 +0000)]
Fix typo reported by Elena (Oracle)
julie [Sat, 15 Nov 2014 05:01:47 +0000 (05:01 +0000)]
Bug Fix form Elena (Oracle) - Correct test of error exits
julie [Sat, 15 Nov 2014 04:52:36 +0000 (04:52 +0000)]
Bux fix from elena Ivanova (Oracle) Rearrange If statements for not overwriting info value
julie [Fri, 31 Oct 2014 15:51:06 +0000 (15:51 +0000)]
Change INFO value when testing N - N is first parameter - INFO should be -1 - bug reported by Elena Ivanova (oracle)
julie [Fri, 31 Oct 2014 15:48:04 +0000 (15:48 +0000)]
Set default value for INFO - bug reported by Elena Ivanova (oracle) and David Hough
julie [Mon, 20 Oct 2014 16:31:31 +0000 (16:31 +0000)]
Comment Fix from Jim: Update of Purpose section of Jacobi SVD about accuracy
julie [Mon, 20 Oct 2014 16:23:51 +0000 (16:23 +0000)]
Comment fix from Jim: argument A in xsyconv.f should be in,out
julie [Fri, 17 Oct 2014 16:46:43 +0000 (16:46 +0000)]
From Jim: Routine description improvement - Adding guarantees that each R(i,i) >= 0 in comment and reference to lawn 203 in Furhter details section
julie [Thu, 9 Oct 2014 01:12:17 +0000 (01:12 +0000)]
correct typo while editing Makefile with VI
julie [Wed, 8 Oct 2014 18:40:57 +0000 (18:40 +0000)]
Fixing folder uppercase / lower case issue - Thank you Don
julie [Tue, 7 Oct 2014 04:21:49 +0000 (04:21 +0000)]
Commit fix from Osni Marquez for Bug 115 reported from Duncan Po (Mathworks)
Description of bug http://icl.utk.edu/lapack-forum/viewtopic.php?f=13&t=4391,
It turns out to be related to bug http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=5&t=4498 as well. Further, I suspect that the problem reported in an e-mail sent by Justin Weiguang Si to lapack@cs.utk.edu on Sep 22 is also related to those, although Justin has not yet replied to my request for the offending matrix...
The bug has been traced to [d/s]laed6, which computes the root closest to the origin of a secular equation and is used in the D&C tridiagonal eigensolver (DSTEDC) and D&C least squares solver (DGELSD). We have interacted with Ren-Cang Li (the original developer of [d/s]laed6) about possible fixes, and I am attaching a new version of [d/s]laed6. This version has been tested with our 'torture cases' (more on this below) in:
- single and double precision
- Intel Xeon Westmere with gnu (gfortran) and intel (ifort) compilers
- AMD MagnyCours with pgi, cray, intel and gnu compilers
The bug was related to a too stringent tolerance convergence criterion, line 390 in laed6. The fix is
390 IF( ( ABS( F ).LE.FOUR*EPS*ERRETM ) .OR.
391 $ ( (UBD-LBD).LE.FOUR*EPS*ABS(TAU) ) )
to replace
390 IF( ABS( F ).LE.EPS*ERRETM )
corrected bugs: 115 and 121
julie [Tue, 7 Oct 2014 04:13:57 +0000 (04:13 +0000)]
Fix problem in Info return value when checking input parameter
Definition of xLAED7( ICOMPQ, N, QSIZ, TLVLS, CURLVL, CURPBM, D, Q, [ etc..]
QSIZ is the third argument, not the fourth
ELSE IF( ICOMPQ.EQ.1 .AND. QSIZ.LT.N ) THEN
- INFO = -4
+ INFO = -3
Bug reported by Elena Ivanova (Oracle) directly to Julie on Oct 2nd 2014
Bug labelled as Bug122
julie [Mon, 22 Sep 2014 15:38:39 +0000 (15:38 +0000)]
Correction to bug120
Fix workspace definition in comment following bug report from Rasmus Munk Larsen, Google Knowledge & Research
See forum https://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=13&t=4587
See revision 1349, and 649 for history
langou [Tue, 16 Sep 2014 03:51:58 +0000 (03:51 +0000)]
Typos in the comments of [s/d]tpmqrt.f reported by Andreas Noack Jensen on Sept. 15 2014.
See LAPACK forum: http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=13&t=4583
julie [Sat, 23 Aug 2014 01:53:37 +0000 (01:53 +0000)]
Clean CMAKE build for includes
julie [Sat, 23 Aug 2014 01:31:15 +0000 (01:31 +0000)]
Update make.inc's following the integration of CBLAS
julie [Sat, 23 Aug 2014 01:21:57 +0000 (01:21 +0000)]
Commit to integrate CMAKE and MAKEFILE build system for CBLAS
Also modify the way the CBLAS mangling was generated...now CBLAS follows LAPACKE model
Details on MAKEFILE build
To compile CBLAS library from the LAPACK folder:
- make cblaslib --> create the cblas lib
- make cblas_testing --> compile BLAS Ref if necessary, build exe, then run exe
Details on CMAKE build
Now you have an option to select CBLAS
If CBLAS is check the Library is built and if BUILD_TESTING is enable, test and example are built and run
following
--> modification of all make.inc
--> more clean up
--> debugging on other platforms
** PLEASE TRY and SEND FEEDBACK - Thanks - Julie **
julie [Fri, 15 Aug 2014 20:14:26 +0000 (20:14 +0000)]
langou [Sat, 12 Jul 2014 20:30:57 +0000 (20:30 +0000)]
See post on the forum from Pauli Virtanen.
http://icl.cs.utk.edu/lapack-forum/viewtopic.php?t=4560
Documentation of xSTEVR says that VL and VU are not referenced if RANGE = 'A'
or 'I', however the code does reference VL and VU for any value of RANGE:
VLL = VL
VUU = VU
So simply protect these two lines of code:
IF( VALEIG ) THEN
VLL = VL
VUU = VU
END IF
Julien.
langou [Sat, 12 Jul 2014 19:58:39 +0000 (19:58 +0000)]
Edited the code of xGBMV, xGEMV and xGEMM by removing the test for an entry of
B (or X) being zero to skip a loop so as to enable better NaN (or Inf)
propagation.
Taking DGEMM as an example, if you use notation:
C(I,J) = C(I,J) + alpha * A(I,L) * B(L,J),
the reference BLAS DGEMM is using the JLI version of matrix matrix multiply and
is checking, for each J (from 1 to N) and for each L (from 1 to K), whether
B(L,J) is zero (or not) to save (or not) the 2M following operations.
(See the "IF (B(L,J).NE.ZERO) THEN" in the code below.)
The snippets of code is as follows
DO 90 J = 1,N
DO 80 L = 1,K
IF (B(L,J).NE.ZERO) THEN
TEMP = ALPHA*B(L,J)
DO 70 I = 1,M
C(I,J) = C(I,J) + TEMP*A(I,L)
70 CONTINUE
END IF
80 CONTINUE
90 CONTINUE
This induces some non NaN-propagation in a pretty ad-hoc way. For better NaN
propagation, this patch removes the above IF statement.
The snippet of code now becomes
DO 90 J = 1,N
DO 80 L = 1,K
TEMP = ALPHA*B(L,J)
DO 70 I = 1,M
C(I,J) = C(I,J) + TEMP*A(I,L)
70 CONTINUE
80 CONTINUE
90 CONTINUE
This enables correct NaN propagation for this piece of code.
Rationale: BLAS does not correctly propagate all NaNs (and Infs). We still have
no NaN propagation where for example ALPHA=0, etc. The goal of this commit is
to have correct NaN propagation no matter what the entries of the input
matrices/vectors (A, B, C, X, etc.) are. BLAS do not correctly propagate NaNs
and Infs based on some values of the scalars (ALPHA, BETA, etc.).
See below the email from Tom Callaway from RedHat, sent on July 9th to
lapack@cs.utk.edu.
Hello LAPACK people,
Martyn & Lejeczek (on CC) reported an issue to Fedora relating to R using our
system copy of BLAS (from LAPACK).
As noted in the R administration and Installation Manual, "R relies on ISO/IEC
60559 compliance of an external BLAS. This can be broken if for example the
code assumes that terms with a zero factor are always zero and do not need to
be computed - whereas x*0 can be NaN. This is checked in the test suite."
In the stock BLAS, DGBMV, DGEMM, and DGEMV fail this. R has been patching their
bundled BLAS to resolve this issue since 2010, but Fedora now uses the system
BLAS.
Attached is a patch (from upstream R) to fix this issue in the LAPACK BLAS.
Please consider applying it.
Thanks,
~tom
julie [Wed, 9 Jul 2014 16:43:58 +0000 (16:43 +0000)]
fix lapack_testing.py for windows platforms - Fix provided by amran.illoul@ensam.eu on July 9th
langou [Fri, 27 Jun 2014 11:46:34 +0000 (11:46 +0000)]
Bug in comments on array size reported by a Sergey Kuznetsov from Intel.
The size of array RWORK in {c,z}dbsqr was given incorrectly.
Previous version for cbdsqr.f read
RWORK is REAL array, dimension (2*N)
if NCVT = NRU = NCC = 0, (max(1, 4*N-4)) otherwise
Now reads
RWORK is REAL array, dimension (4*N)
langou [Mon, 16 Jun 2014 12:16:58 +0000 (12:16 +0000)]
bug found by Themos Tsikas fron NAG Ltd
LAPACK FORUM topic #4558
Hello
In zlanhf.f, first executable statements
Code: Select all
IF( N.EQ.0 ) THEN
ZLANHF = ZERO
RETURN
ELSE IF( N.EQ.1 ) THEN
ZLANHF = ABS(A(0))
RETURN
END IF
When N=1, we should be ignoring any imaginary part (as we do for general diagonal elements later on) and so it should read
Code: Select all
IF( N.EQ.0 ) THEN
ZLANHF = ZERO
RETURN
ELSE IF( N.EQ.1 ) THEN
ZLANHF = ABS(DBLE(A(0)))
RETURN
END IF
Themos Tsikas
NAG Ltd
julie [Sat, 14 Jun 2014 17:04:05 +0000 (17:04 +0000)]
Applied patch provided from vitaut on LAPACK forum on June 13th.
The patch fixes the warnings (Policy CMP0026) by replacing the deprecated LOCATION target
property with the generator expression $<TARGET_FILE>
See http://icl.cs.utk.edu/lapack-forum/viewtopic.php?f=13&t=4556&p=10939#p10939
julie [Tue, 15 Apr 2014 20:19:45 +0000 (20:19 +0000)]
Update version number in CMAKE - Thanks Allan Wollaber for the bug report (April 15th 2014)
julie [Tue, 15 Apr 2014 20:18:55 +0000 (20:18 +0000)]
Update version number in CMAKE - Thanks Allan Wollaber for the bug report (April 15th 2014)
julie [Mon, 7 Apr 2014 15:40:04 +0000 (15:40 +0000)]
Adding missing routine in repository
julie [Wed, 19 Mar 2014 01:15:50 +0000 (01:15 +0000)]
Corrected minor bug submitted by Elena Ivanova, Oracle on March 18th 2014.
From Elena:
there is a minor bug in /LIN/dchkql.f:
Array IWORK (a parameter of this routine) is not used at all.
It would be better to remove IWORK from the list of parameters,
remove its declaration and description in the code of /LIN/dchkql.f:
Note: only dchkql had that IWORK parameter, [c,s,z]chkql don't have it.
julie [Thu, 20 Feb 2014 15:58:23 +0000 (15:58 +0000)]
Fix typos found by Jeremy Du Croz
julie [Tue, 18 Feb 2014 16:47:10 +0000 (16:47 +0000)]
INTEGRATED CHANGES PROVIDED BY INTEL TEAM FOR LAPACKE -- SEE Eugene Chereshnev email on Feb 18th
--> THANK YO INTEL <--
We performed auto-replacement matrix_order -> matrix_layout in the LAPACKE C Interface.
Also we applied several bug fixes which were presented on our side and added two LAPACKE interfaces for
function LAPACK_stpqrt which weren't represented on your side (although interfaces for d/c/z/pqrt were).
All changes are reflected in the archive in attachment (only modified files in the archive). All changes were performed for revision 1471.
Summary of changes:
1. Files:
lapacke_stprfb.c
lapacke_dtprfb.c
lapacke_ctprfb.c
lapacke_ztprfb.c
lapacke_stprfb_work.c
lapacke_dtprfb_work.c
lapacke_ctprfb_work.c
lapacke_ztprfb_work.c
lapacke.h
Expression MAX(1,ldwork)*MAX(n,k) was replaced by 'work_size'. It's not a bug but apparently 'work_size' is right expression in this place.
Also in functions lapacke_(c/z)tprfb type of 'work' array was changed from float/double to
lapack_complex_float/lapack_complex_double. Array 'work' is passed to lapacke_(c/z)tprfb_work and it calls
LAPACK_(c/z)tprfb (and passes array 'work') but LAPACK_(c/z)tprfb has argument 'work' of complex type:
SRC/ctprfb.f:
*> WORK is COMPLEX array, dimension
*> (LDWORK,N) if SIDE = 'L',
*> (LDWORK,K) if SIDE = 'R'.
and SRC/ztprfb.f:
*> WORK is COMPLEX*16 array, dimension
*> (LDWORK,N) if SIDE = 'L',
*> (LDWORK,K) if SIDE = 'R'.
Types of 'work' array in lapacke.h for LAPACK_(c/z)tprfb are float/double so they aren't correct.
The archive includes changes in lapacke.h and source files:
types of 'work' array in LAPACK_(c/z)tprfb: float/double -> lapack_complex_float/lapack_complex_double
types of 'work' array in LAPACKE_(c/z)tprfb_work: float/double -> lapack_complex_float/lapack_complex_double
In a similar way types of allocated 'work' arrays in LAPACKE_(c/z)tprfb were replaced.
Also you can see here why MAX(1,ldwork)*MAX(n,k) isn't fully correct expression for work_size.
2. Files:
lapacke_cuncsd.c
lapacke_dorcsd.c
lapacke_sorcsd.c
lapacke_zuncsd.c
From documentation of LAPACK_cuncsd (for other three LAPACK functions IWORK must have same size):
*> IWORK is INTEGER array, dimension (M-MIN(P,M-P,Q,M-Q))
So it's not correct here to allocate working array with MAX(1,m-q) size.
Changes:
lapack_int r;
r=MIN(p,m-p);
r=MIN(r,q);
r=MIN(r,m-q);
MAX(1,m-q) -> MAX(1,m-r)
3. Files:
lapacke_cgbsvxx_work.c
lapacke_cgesvxx_work.c
lapacke_chesvxx_work.c
lapacke_cposvxx_work.c
lapacke_csysvxx_work.c
lapacke_dgbsvxx_work.c
lapacke_dgesvxx_work.c
lapacke_dposvxx_work.c
lapacke_dsysvxx_work.c
lapacke_sgbsvxx_work.c
lapacke_sgesvxx_work.c
lapacke_sposvxx_work.c
lapacke_ssysvxx_work.c
lapacke_zgbsvxx_work.c
lapacke_zgesvxx_work.c
lapacke_zhesvxx_work.c
lapacke_zposvxx_work.c
lapacke_zsysvxx_work.c
Functions doesn't perform right transposition for 'err_bnds_norm' and 'err_bnds_comp' matrices:
From documentation of these functions:
*> ERR_BNDS_NORM is (SINGLE/DOUBLE) PRECISION array, dimension (NRHS, N_ERR_BNDS)
*> ERR_BNDS_COMP is (SINGLE/DOUBLE) PRECISION array, dimension (NRHS, N_ERR_BNDS)
So for interface functions right 'lda' for 'err_bnds_norm' and 'err_bnds_comp' matrices is 'n_err_bnds' rather than 'nhrs' in the case of LAPACK_ROW_MAJOR.
4. Files:
lapacke_stpqrt.c
lapacke_stpqrt_work.c
lapacke/src/Makefile
lapacke/src/CMakeLists.txt
Interfaces for LAPACK_spqrt were added: lapacke_stpqrt.c and lapacke_stpqrt_work.c and appropriate changes were performed in Makefile and CMakeLists.txt.
Declaration of LAPACK_spqrt was added to lapacke.h.
julie [Fri, 14 Feb 2014 21:40:56 +0000 (21:40 +0000)]
Correct a bug in the description of the WORK parameter in (c,z)la_syrpvgrw.f
WORK is real
Bug report sent from Elena Ivanova (Oracle) on Feb 14th 2014
julie [Fri, 14 Feb 2014 21:13:20 +0000 (21:13 +0000)]
Fix potential initiation issue in (c,d,s,z)lansy.f
Bug report sent from Elena Ivanova on Feb 2014
Array WORK is initialized with ZERO when 'L' and is not initialized when 'U'.
There can be some garbage in WORK. Move up the initialization of WORK before the line
IF( LSAME( UPLO, 'U' ) ) THEN.
Rearrange loop ordering numbers.
langou [Mon, 10 Feb 2014 19:55:39 +0000 (19:55 +0000)]
( Let me know if I broke anything. I double checked. We should be fine but )
( please let me know if concerns/problems. )
1) The comments in ICMAX1 (resp. IZMAX1) were not correct. The comments read:
"ICMAX1 finds the index of the element whose real part has maximum absolute
value." It should have read: "ICMAX1 finds the index of the first vector
element of maximum absolute value." This is corrected. The problem was reported
Eloy Romero Alcalde from the SLEPc team, Universitat Politècnica de València
and confirmed by Nick Higham, Sven Hammarling and Julien Langou.
2) The routine ICAMAX (resp. IZAMAX) from the BLAS evolved without ICMAX1
(resp. IZMAX1). Reconcialated both version. So essentially took the ICAMAX
currently in the BLAS and changed it appropriately for an ICMAX1.
3) Remove the use of statement function in ICMAX1 (resp. IZMAX1). Now obsolete
in FORTRAN.
4) Change comments in BLAS routines: SCABS1, DCABS1, ICAMAX, IZAMAX, SCASUM,
and DZASUM. Remove the use of "absolute value of a complex number" for the
quantity "| Re(.) | + | Im(.) |". For example:
before
DCABS1 computes absolute value of a double complex number
after
DCABS1 computes |Re(.)| + |Im(.)| of a double complex number
julie [Sun, 9 Feb 2014 00:37:45 +0000 (00:37 +0000)]
== Patch provided by Brad King from Kitware - brad.king@kitware.com ==
Provide CMake packages for both LAPACK and LAPACKE
Teach "lapack-config.cmake" to provide variables
LAPACK_blas_LIBRARIES
LAPACK_lapack_LIBRARIES
that contain either the target names when using the reference
implementation or the system libraries found for them.
Configure a "lapacke-config.cmake" file for the build and install
trees to package LAPACKE. Teach it to load the LAPACK package
installed with it. Provide variables
LAPACKE_INCLUDE_DIRS
LAPACKE_LIBRARIES
containing the header file search path for lapacke headers
and the list of lapacke library targets.
This requires CMake 2.8.10 to separate the installation export
for the lapacke library from the other targets.
julie [Thu, 9 Jan 2014 20:34:39 +0000 (20:34 +0000)]
Problem reported by zhaowei@sccas.cn on Jan 9th 2014
Corrected comments to list Z as input/output
> Email
i just found a mistake in the comments of 'dlasq3.c' while reading the code:
/* Z (input) DOUBLE PRECISION array, dimension ( 4*N ) */
/* Z holds the qd array. */
but Z can be changed in that subroutine...i just think '(output)' should be added to avoid misleading user...
julie [Mon, 9 Dec 2013 17:55:10 +0000 (17:55 +0000)]
Do not hide the CMAKE_GNUtoMS variable to users (Thanks Brad)
julie [Sat, 7 Dec 2013 00:32:08 +0000 (00:32 +0000)]
Fix typo in comment reported by elena.x.ivanova@oracle.com - zhetri2x fix was committed in r1465 (oops)
julie [Sat, 7 Dec 2013 00:25:51 +0000 (00:25 +0000)]
Fixed typos reported by elena.x.ivanova@oracle.com on Dec 5th
julie [Fri, 6 Dec 2013 23:13:50 +0000 (23:13 +0000)]
Remove unused leftover GNUtoMS support file
julie [Fri, 6 Dec 2013 03:59:31 +0000 (03:59 +0000)]
Remove GNUtoMS code in favor of CMake builtin version - Thanks Brad (Kitware)
julie [Fri, 6 Dec 2013 03:58:29 +0000 (03:58 +0000)]
Require CMake >= 2.8.7 so we can use CMAKE_GNUtoMS
lawrence.mulholland [Fri, 29 Nov 2013 16:08:47 +0000 (16:08 +0000)]
fix to bug report 4331
julie [Wed, 27 Nov 2013 16:07:34 +0000 (16:07 +0000)]
Split lines to make INtel Windows Compiler happy...
langou [Tue, 26 Nov 2013 19:36:29 +0000 (19:36 +0000)]
Correct bug reported by Alex Zotkevich from Intel
See: http://icl.utk.edu/lapack-forum/viewtopic.php?f=13&t=4392
Bug: During workspace computation, LAPACK code CGESVD was calling other
subroutines (e.g. CGEQRF) with REAL DUM variable as COMPLEX WORK variable. DUM
(in CGESVD) is REAL while WORK (in called subroutines) is COMPLEX. This
corrupts the stack when a value is set in WORK.
Fix: In CGESVD, use the COMPLEX CDUM variable (already present in the code)
instead of the REAL DUM variable. Since I was at it, the COMPLEX "TAU"
variables (not referenced anyway) were passed the REAL DUM variable, I changed
the code so that the COMPLEX CDUM variable is passed. This is cleaner like
this.
Same problem with ZGESVD. Same fix.
Alex's post:
Hi, We recently found a stack corruption issue in (C,Z)GESVD that potentially
could even lead to incorrect xerbla error message. In ZGESVD array DUM which
is used in LWORK query is a double precision array of size 1 allocated on
stack:
DOUBLE PRECISION DUM( 1 )
DUM comes to ( ZGEQRF, ZUNGQR, ... ) as a WORK array to return an optimal LWORK
value. But in ( ZGEQRF, ZUNGQR, ... ) array WORK is declared as a COMPLEX*16
array. So WORK(1) = 1 corrupts the stack as it deals with complex value while
pointer on input of the function is a pointer to double: (oooooooo|xxxxxxxx),
oooooooo fills with LWORK value, xxxxxxxx corrupts. Let compiler use xxxxxxxx
to hold some value. After LWORK query the value will turn to be a zero. "Hacky
fix" would be to allocate DUM array of size 2.
W.B.R.
Alex Zotkevich
langou [Mon, 25 Nov 2013 23:47:55 +0000 (23:47 +0000)]
Apply patch from Nadezhda Mozartova from Intel. (Sent to Julien on November
20th.) Thanks Nadezhda! Thanks Intel!
julie [Wed, 20 Nov 2013 18:57:47 +0000 (18:57 +0000)]
Fix info value problem for (c,d,s,z)tgsy2.f reported by Elena Ivanova (ORACLE) - Sent Nov 20th
julie [Tue, 19 Nov 2013 10:18:52 +0000 (10:18 +0000)]
Applying 2nd part of Mathieu's patch - modification in the include file
julie [Tue, 19 Nov 2013 09:52:10 +0000 (09:52 +0000)]
Upadte version and date
julie [Tue, 19 Nov 2013 09:43:31 +0000 (09:43 +0000)]
Apply patch submitted by Matthew Faverge, INRIA on Nov 19th 2013
There are problems in LAPACKE complex lacn2 interfaces. Those functions don't have a "isgn" parmeter in LAPACK. This exists only for real interface. The problem is present in the four files of C and Z functions.
In Lapack:
SUBROUTINE ZLACN2( N, V, X, EST, KASE, ISAVE )
In Lapacke:
LAPACK_zlacn2( &n, v, x, isgn, est, kase, isave );
igor175 [Mon, 18 Nov 2013 23:51:04 +0000 (23:51 +0000)]
LAPACK bug fix for XBLAS: bug revealed in testing of linear equation routines, when XBLAS is linked
julie [Sun, 17 Nov 2013 01:17:17 +0000 (01:17 +0000)]
Update version and date on make.inc templates
julie [Sun, 17 Nov 2013 01:15:04 +0000 (01:15 +0000)]
fix Illegal grep pb in lapack_testing
julie [Sun, 17 Nov 2013 01:08:52 +0000 (01:08 +0000)]
Slight modification to Lawrence patch for xTPQRT - Update conditions for testing parameter L - MIN(M,N) can be 0
julie [Sun, 17 Nov 2013 01:06:41 +0000 (01:06 +0000)]
update ilaver
julie [Sun, 17 Nov 2013 00:25:52 +0000 (00:25 +0000)]
Update release number and date
julie [Sun, 17 Nov 2013 00:01:25 +0000 (00:01 +0000)]
Patch provided by Lawrence Mulholland from Nag on Nov 1st 2013
Email below:
============
I have been incorporating some routines into the NAG Library,
which means some automatic code translation and writing some
example and test programs.
The routines I have been adding are:
?geqrt, ?gemqrt, ?tpqrt, ?tpmqrt, ?orcsd, ?uncsd
At the end of this message I will give you my current svn status
and svn diff for consideration and approval before I commit.
In each case, when testing immediate exits, my tests failed because
constraints were mutually exclusive for the immediate return case.
I have already committed changes to the constraints for some of
the above to allow immediate exit.
I have completed this for the remainder of this set.
Less importantly, there are things in the code that trip up a checking
compiler:
a) an
IF ( clause1(i) .AND. clause2(array(i)) ) THEN
where array(i) is either not initialized or is out of bounds if
clause1(i) is .FALSE.
This is wrong since a Fortran compiler is at liberty to test clause2 first.
In my changes this has been split into two as best suits the case.
b) an
CALL SUB (i, array(N-i+2))
with i = 1 and array(N+1) either not initialized or out of bounds, but
internally array(N+1) is not referenced.
In this case I don't think the Fortran standard is clear, but it trips up the
nagfor compiler with checking on. So in the NAG incorporated versions
of Lapack routines such calls are protected and/or
a special i=1 call is made.
The changes I want to commit also do this.
c) workspace queries passing zero instead of array references
e.g.
lwork = -1
call barf(n,m,0,0,0,0,0,-1,info)
a checking compiler won't like this.
I have changed cases like this to pass available arrays of sufficient size
and the right shape in place of the zeros.
langou [Tue, 15 Oct 2013 02:57:28 +0000 (02:57 +0000)]
Typo in the comment on complex xHSEQR.
The formula:
(QZ)*H*(QZ)**H.
should read:
(QZ)*T*(QZ)**H.
james [Sat, 5 Oct 2013 21:19:17 +0000 (21:19 +0000)]
Added check for NaN in input parameter G to avoid infinite loop -- fixes bug110.
langou [Fri, 4 Oct 2013 19:55:20 +0000 (19:55 +0000)]
FIX BUG 112.
(Bug reported on Aug 18 2013 by Daniel Strobusch on LAPACK mailing list)
(Bug fixed by Rodney on Aug 19 2013)
This is a problem with gfortran. One need to force gfortran to allocate all
local arrays on the stack. One way to do this is to have the option
-frecursive
julie [Fri, 4 Oct 2013 03:20:34 +0000 (03:20 +0000)]
Remove unused code
julie [Mon, 16 Sep 2013 16:44:21 +0000 (16:44 +0000)]
Correct comment
Problem reported by Elena Ivanova <Elena.x.Ivanova@Oracle.com> on Sept 16th
"Both RES, and AYB have to be described as REAL."
julie [Fri, 13 Sep 2013 21:32:40 +0000 (21:32 +0000)]
Apply fix sent by Elena Ivanova (Oracle)
There is FPE in the 457th line of CDRVRFP:
cdrvrfp.f
(dbx) where
=>[1] cdrvrfp(nout = 6, nn = 9, nval = ARRAY, nns = 3, nsval = ARRAY, nnt = 9, ntval = ARRAY, thresh = 30.0, a = ARRAY, asav = ARRAY, afac = ARRAY, ainv = ARRAY, b = ARRAY, bsav = ARRAY, xact = ARRAY, x = ARRAY, arf = ARRAY, arfinv = ARRAY, c_work_clatms = ARRAY, c_work_cpot02 = ARRAY, c_work_cpot03 = ARRAY, s_work_clatms = ARRAY, s_work_clanhe = ARRAY, s_work_cpot01 = ARRAY, s_work_cpot02 = ARRAY, s_work_cpot03 = ARRAY), line 457 in "cdrvrfp.f"
MAIN(), line 244 in "cchkrfp.f"
When N=0, attempt to calculate a condition number → FPE, because we try to delete by zero (ANORM = 0.0, AINVNM = 0.0) in the 457th line:
Line 457:
RCONDC = ( ONE / ANORM ) / AINVNM
julie [Fri, 13 Sep 2013 21:24:08 +0000 (21:24 +0000)]
Add small fixes sent by Ake Sandgren <ake.sandgren@hpc2n.umu.se>
julie [Thu, 12 Sep 2013 15:42:59 +0000 (15:42 +0000)]
Fix "if test" semantic in xdrvpox
Bug report sent by Elena Ivanova <elena.x.ivanova@oracle.com> on Sept 10th
lawrence.mulholland [Wed, 21 Aug 2013 10:39:32 +0000 (10:39 +0000)]
typo in previous fix: MIN(M,N)>0 NOT MIN(M,N)<0 when checking NB>MIN(M,N)
langou [Fri, 16 Aug 2013 15:03:03 +0000 (15:03 +0000)]
Fix from Lawrence Mulholland (NAG).
From Lawrence:
The *GEQRT routines are supposed to return immediately when N=0 or M=0.
However the code is such that an INFO = -3 is triggered since we are
trying to satisfy
ELSE IF( NB.LT.1 .OR. NB.GT.MIN(M,N) )THEN
INFO = -3
should be
ELSE IF( NB.LT.1 .OR. ( NB.GT.MIN(M,N) .AND. MIN(M,N).LT.0 ) )THEN
INFO = -3
Also on a side note, I think we could allow NB to be greater than MIN(M,N) in
the interface and then reset NB as needed if the user's NB is too large.
Another day. Good enough for now.
Julien.
langou [Mon, 22 Jul 2013 23:01:15 +0000 (23:01 +0000)]
apply patch from Hong Xu
patch sent to lapack mailing list on Sunday July 21st
In the input parameter checking, LDZ was not checked correctly
LDZ was checked against M instead of N essentially
julie [Sat, 20 Jul 2013 19:29:54 +0000 (19:29 +0000)]
Add missing file in CMakeLists.txt
igor175 [Tue, 9 Jul 2013 04:45:46 +0000 (04:45 +0000)]
Added LAPACK tests for 'rook' routines and drivers to LIN/Makefile LIN/CMakeLists.txt
igor175 [Tue, 9 Jul 2013 04:42:19 +0000 (04:42 +0000)]
Modified test files for 'rook' pivoting LAPACK routines: LIN/cchkaa.f LIN/cchkhe_rook.f LIN/cdrvhe_rook.f LIN/cerrhe.f LIN/cerrvx.f LIN/zchkaa.f LIN/zchkhe_rook.f LIN/zdrvhe_rook.f LIN/zerrhe.f LIN/zerrvx.f ctest.in ztest.in
igor175 [Tue, 9 Jul 2013 04:06:39 +0000 (04:06 +0000)]
Added 'rook' LAPACK routines and drivers to SRC/Makefile SRC/CMakeLists.txt
igor175 [Tue, 9 Jul 2013 04:02:51 +0000 (04:02 +0000)]
Modified new 'rook' pivoting routines and drivers for Hermitian indefinite matrices, Complex precision: zhetf2_rook.f zlahef_rook.f zhetri_rook.f
igor175 [Tue, 9 Jul 2013 01:55:07 +0000 (01:55 +0000)]
Modified new 'rook' pivoting routines and drivers for Hermitian indefinite matrices, Complex precision: chetf2_rook.f clahef_rook.f chetri_rook.f
igor175 [Thu, 4 Jul 2013 00:52:02 +0000 (00:52 +0000)]
modified error hadling routines for LIN testing of new Hermitian 'rook' pivoting code: aladhd.f