diff options
| author | langou <langou@users.noreply.github.com> | 2009-01-12 05:17:21 +0000 |
|---|---|---|
| committer | langou <langou@users.noreply.github.com> | 2009-01-12 05:17:21 +0000 |
| commit | 3cff0c0c7d6f581bb14e645fbd5d035a53581446 (patch) | |
| tree | b16f6163f51de9ac12a691dcb407c55cb5086ed0 /SRC/stfsm.f | |
| parent | dc20e8382624a09011f4f57fac2f0dd543a0522d (diff) | |
| download | lapack-3cff0c0c7d6f581bb14e645fbd5d035a53581446.tar.gz lapack-3cff0c0c7d6f581bb14e645fbd5d035a53581446.tar.bz2 lapack-3cff0c0c7d6f581bb14e645fbd5d035a53581446.zip | |
============================================================================
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 ...
Diffstat (limited to 'SRC/stfsm.f')
| -rw-r--r-- | SRC/stfsm.f | 70 |
1 files changed, 46 insertions, 24 deletions
diff --git a/SRC/stfsm.f b/SRC/stfsm.f index 2ca36fb0..05ca537c 100644 --- a/SRC/stfsm.f +++ b/SRC/stfsm.f @@ -322,24 +322,34 @@ * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and * TRANS = 'N' * - CALL STRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, - + A( 0 ), M, B, LDB ) - CALL SGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ), M, - + B, LDB, ALPHA, B( M1, 0 ), LDB ) - CALL STRSM( 'L', 'U', 'T', DIAG, M2, N, ONE, - + A( M ), M, B( M1, 0 ), LDB ) + IF( M.EQ.1 ) THEN + CALL STRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, + + A, M, B, LDB ) + ELSE + CALL STRSM( 'L', 'L', 'N', DIAG, M1, N, ALPHA, + + A( 0 ), M, B, LDB ) + CALL SGEMM( 'N', 'N', M2, N, M1, -ONE, A( M1 ), + + M, B, LDB, ALPHA, B( M1, 0 ), LDB ) + CALL STRSM( 'L', 'U', 'T', DIAG, M2, N, ONE, + + A( M ), M, B( M1, 0 ), LDB ) + END IF * ELSE * * SIDE ='L', N is odd, TRANSR = 'N', UPLO = 'L', and * TRANS = 'T' * - CALL STRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA, - + A( M ), M, B( M1, 0 ), LDB ) - CALL SGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ), M, - + B( M1, 0 ), LDB, ALPHA, B, LDB ) - CALL STRSM( 'L', 'L', 'T', DIAG, M1, N, ONE, - + A( 0 ), M, B, LDB ) + IF( M.EQ.1 ) THEN + CALL STRSM( 'L', 'L', 'T', DIAG, M1, N, ALPHA, + + A( 0 ), M, B, LDB ) + ELSE + CALL STRSM( 'L', 'U', 'N', DIAG, M2, N, ALPHA, + + A( M ), M, B( M1, 0 ), LDB ) + CALL SGEMM( 'T', 'N', M1, N, M2, -ONE, A( M1 ), M, + + M, B( M1, 0 ), LDB, ALPHA, B, LDB ) + CALL STRSM( 'L', 'L', 'T', DIAG, M1, N, ONE, + + A( 0 ), M, B, LDB ) + END IF * END IF * @@ -388,24 +398,36 @@ * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and * TRANS = 'N' * - CALL STRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, - + A( 0 ), M1, B, LDB ) - CALL SGEMM( 'T', 'N', M2, N, M1, -ONE, A( M1*M1 ), - + M1, B, LDB, ALPHA, B( M1, 0 ), LDB ) - CALL STRSM( 'L', 'L', 'N', DIAG, M2, N, ONE, - + A( 1 ), M1, B( M1, 0 ), LDB ) + IF( M.EQ.1 ) THEN + CALL STRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, + + A( 0 ), M1, B, LDB ) + ELSE + CALL STRSM( 'L', 'U', 'T', DIAG, M1, N, ALPHA, + + A( 0 ), M1, B, LDB ) + CALL SGEMM( 'T', 'N', M2, N, M1, -ONE, + + A( M1*M1 ), M1, B, LDB, ALPHA, + + B( M1, 0 ), LDB ) + CALL STRSM( 'L', 'L', 'N', DIAG, M2, N, ONE, + + A( 1 ), M1, B( M1, 0 ), LDB ) + END IF * ELSE * * SIDE ='L', N is odd, TRANSR = 'T', UPLO = 'L', and * TRANS = 'T' * - CALL STRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA, - + A( 1 ), M1, B( M1, 0 ), LDB ) - CALL SGEMM( 'N', 'N', M1, N, M2, -ONE, A( M1*M1 ), - + M1, B( M1, 0 ), LDB, ALPHA, B, LDB ) - CALL STRSM( 'L', 'U', 'N', DIAG, M1, N, ONE, - + A( 0 ), M1, B, LDB ) + IF( M.EQ.1 ) THEN + CALL STRSM( 'L', 'U', 'N', DIAG, M1, N, ALPHA, + + A( 0 ), M1, B, LDB ) + ELSE + CALL STRSM( 'L', 'L', 'T', DIAG, M2, N, ALPHA, + + A( 1 ), M1, B( M1, 0 ), LDB ) + CALL SGEMM( 'N', 'N', M1, N, M2, -ONE, + + A( M1*M1 ), M1, B( M1, 0 ), LDB, + + ALPHA, B, LDB ) + CALL STRSM( 'L', 'U', 'N', DIAG, M1, N, ONE, + + A( 0 ), M1, B, LDB ) + END IF * END IF * |