Skip to content

Commit 43f9441

Browse files
martin-frbgmartin-frbg
authored
Merge pull request #5920 from martin-frbg/lapack1245
Avoid large intermediates in (C/Z)LARTG (Reference-LAPACK PR 1245)
2 parents 3e48730 + 53b1498 commit 43f9441

6 files changed

Lines changed: 20 additions & 24 deletions

File tree

‎lapack-netlib/SRC/cbbcsd.f‎

Lines changed: 3 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -351,9 +351,9 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
351351
* .. Parameters ..
352352
INTEGER MAXITR
353353
PARAMETER ( MAXITR = 6 )
354-
REAL HUNDRED, MEIGHTH, ONE, TEN, ZERO
354+
REAL HUNDRED, MEIGHTH, ZERO, ONE, TEN
355355
PARAMETER ( HUNDRED = 100.0E0, MEIGHTH = -0.125E0,
356-
$ ONE = 1.0E0, TEN = 10.0E0, ZERO = 0.0E0 )
356+
$ ZERO = 0.0E0, ONE = 1.0E0, TEN = 10.0E0 )
357357
COMPLEX NEGONECOMPLEX
358358
PARAMETER ( NEGONECOMPLEX = (-1.0E0,0.0E0) )
359359
REAL PIOVER2
@@ -576,7 +576,7 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
576576
END IF
577577
ELSE
578578
NU = SIGMA21
579-
MU = SQRT( 1.0 - NU**2 )
579+
MU = SQRT( ONE - NU**2 )
580580
IF( NU .LT. THRESH ) THEN
581581
MU = ONE
582582
NU = ZERO
@@ -1114,4 +1114,3 @@ SUBROUTINE CBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
11141114
* End of CBBCSD
11151115
*
11161116
END
1117-

‎lapack-netlib/SRC/clartg.f90‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -203,7 +203,7 @@ subroutine CLARTG( f, g, c, s, r )
203203
rtmax = rtmax * 2
204204
if( f2 > rtmin .and. h2 < rtmax ) then
205205
! safmin <= sqrt( f2*h2 ) <= safmax
206-
s = conjg( g ) * ( f / sqrt( f2*h2 ) )
206+
s = ( f / sqrt( f2 ) ) * ( conjg( g ) / sqrt( h2 ) )
207207
else
208208
s = conjg( g ) * ( r / h2 )
209209
end if
@@ -223,7 +223,7 @@ subroutine CLARTG( f, g, c, s, r )
223223
! sqrt(safmin) <= f2 * sqrt(safmax) <= h2 / sqrt(f2 * h2) <= h2 * (safmin / f2) <= h2 <= safmax
224224
r = f * ( h2 / d )
225225
end if
226-
s = conjg( g ) * ( f / d )
226+
s = ( f / sqrt( f2 ) ) * ( conjg( g ) / sqrt( h2 ) )
227227
end if
228228
else
229229
!
@@ -259,7 +259,7 @@ subroutine CLARTG( f, g, c, s, r )
259259
rtmax = rtmax * 2
260260
if( f2 > rtmin .and. h2 < rtmax ) then
261261
! safmin <= sqrt( f2*h2 ) <= safmax
262-
s = conjg( gs ) * ( fs / sqrt( f2*h2 ) )
262+
s = ( fs / sqrt( f2 ) ) * ( conjg( gs ) / sqrt( h2 ) )
263263
else
264264
s = conjg( gs ) * ( r / h2 )
265265
end if
@@ -279,7 +279,7 @@ subroutine CLARTG( f, g, c, s, r )
279279
! sqrt(safmin) <= f2 * sqrt(safmax) <= h2 / sqrt(f2 * h2) <= h2 * (safmin / f2) <= h2 <= safmax
280280
r = fs * ( h2 / d )
281281
end if
282-
s = conjg( gs ) * ( fs / d )
282+
s = ( fs / sqrt( f2 ) ) * ( conjg( gs ) / sqrt( h2 ) )
283283
end if
284284
! Rescale c and r
285285
c = c * w

‎lapack-netlib/SRC/dbbcsd.f‎

Lines changed: 3 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -351,9 +351,9 @@ SUBROUTINE DBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
351351
* .. Parameters ..
352352
INTEGER MAXITR
353353
PARAMETER ( MAXITR = 6 )
354-
DOUBLE PRECISION HUNDRED, MEIGHTH, ONE, TEN, ZERO
354+
DOUBLE PRECISION HUNDRED, MEIGHTH, ZERO, ONE, TEN
355355
PARAMETER ( HUNDRED = 100.0D0, MEIGHTH = -0.125D0,
356-
$ ONE = 1.0D0, TEN = 10.0D0, ZERO = 0.0D0 )
356+
$ ZERO = 0.0D0, ONE = 1.0D0, TEN = 10.0D0 )
357357
DOUBLE PRECISION NEGONE
358358
PARAMETER ( NEGONE = -1.0D0 )
359359
DOUBLE PRECISION PIOVER2
@@ -576,7 +576,7 @@ SUBROUTINE DBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
576576
END IF
577577
ELSE
578578
NU = SIGMA21
579-
MU = SQRT( 1.0 - NU**2 )
579+
MU = SQRT( ONE - NU**2 )
580580
IF( NU .LT. THRESH ) THEN
581581
MU = ONE
582582
NU = ZERO
@@ -1108,4 +1108,3 @@ SUBROUTINE DBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
11081108
* End of DBBCSD
11091109
*
11101110
END
1111-

‎lapack-netlib/SRC/sbbcsd.f‎

Lines changed: 3 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -351,9 +351,9 @@ SUBROUTINE SBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
351351
* .. Parameters ..
352352
INTEGER MAXITR
353353
PARAMETER ( MAXITR = 6 )
354-
REAL HUNDRED, MEIGHTH, ONE, TEN, ZERO
354+
REAL HUNDRED, MEIGHTH, ZERO, ONE, TEN
355355
PARAMETER ( HUNDRED = 100.0E0, MEIGHTH = -0.125E0,
356-
$ ONE = 1.0E0, TEN = 10.0E0, ZERO = 0.0E0 )
356+
$ ZERO = 0.0E0, ONE = 1.0E0, TEN = 10.0E0 )
357357
REAL NEGONE
358358
PARAMETER ( NEGONE = -1.0E0 )
359359
REAL PIOVER2
@@ -576,7 +576,7 @@ SUBROUTINE SBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
576576
END IF
577577
ELSE
578578
NU = SIGMA21
579-
MU = SQRT( 1.0 - NU**2 )
579+
MU = SQRT( ONE - NU**2 )
580580
IF( NU .LT. THRESH ) THEN
581581
MU = ONE
582582
NU = ZERO
@@ -1108,4 +1108,3 @@ SUBROUTINE SBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
11081108
* End of SBBCSD
11091109
*
11101110
END
1111-

‎lapack-netlib/SRC/zbbcsd.f‎

Lines changed: 3 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -351,9 +351,9 @@ SUBROUTINE ZBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
351351
* .. Parameters ..
352352
INTEGER MAXITR
353353
PARAMETER ( MAXITR = 6 )
354-
DOUBLE PRECISION HUNDRED, MEIGHTH, ONE, TEN, ZERO
354+
DOUBLE PRECISION HUNDRED, MEIGHTH, ZERO, ONE, TEN
355355
PARAMETER ( HUNDRED = 100.0D0, MEIGHTH = -0.125D0,
356-
$ ONE = 1.0D0, TEN = 10.0D0, ZERO = 0.0D0 )
356+
$ ZERO = 0.0D0, ONE = 1.0D0, TEN = 10.0D0 )
357357
COMPLEX*16 NEGONECOMPLEX
358358
PARAMETER ( NEGONECOMPLEX = (-1.0D0,0.0D0) )
359359
DOUBLE PRECISION PIOVER2
@@ -575,7 +575,7 @@ SUBROUTINE ZBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
575575
END IF
576576
ELSE
577577
NU = SIGMA21
578-
MU = SQRT( 1.0 - NU**2 )
578+
MU = SQRT( ONE - NU**2 )
579579
IF( NU .LT. THRESH ) THEN
580580
MU = ONE
581581
NU = ZERO
@@ -1113,4 +1113,3 @@ SUBROUTINE ZBBCSD( JOBU1, JOBU2, JOBV1T, JOBV2T, TRANS, M, P,
11131113
* End of ZBBCSD
11141114
*
11151115
END
1116-

‎lapack-netlib/SRC/zlartg.f90‎

Lines changed: 4 additions & 4 deletions
Original file line numberDiff line numberDiff line change
@@ -203,7 +203,7 @@ subroutine ZLARTG( f, g, c, s, r )
203203
rtmax = rtmax * 2
204204
if( f2 > rtmin .and. h2 < rtmax ) then
205205
! safmin <= sqrt( f2*h2 ) <= safmax
206-
s = conjg( g ) * ( f / sqrt( f2*h2 ) )
206+
s = ( f / sqrt( f2 ) ) * ( conjg( g ) / sqrt( h2 ) )
207207
else
208208
s = conjg( g ) * ( r / h2 )
209209
end if
@@ -223,7 +223,7 @@ subroutine ZLARTG( f, g, c, s, r )
223223
! sqrt(safmin) <= f2 * sqrt(safmax) <= h2 / sqrt(f2 * h2) <= h2 * (safmin / f2) <= h2 <= safmax
224224
r = f * ( h2 / d )
225225
end if
226-
s = conjg( g ) * ( f / d )
226+
s = ( f / sqrt( f2 ) ) * ( conjg( g ) / sqrt( h2 ) )
227227
end if
228228
else
229229
!
@@ -259,7 +259,7 @@ subroutine ZLARTG( f, g, c, s, r )
259259
rtmax = rtmax * 2
260260
if( f2 > rtmin .and. h2 < rtmax ) then
261261
! safmin <= sqrt( f2*h2 ) <= safmax
262-
s = conjg( gs ) * ( fs / sqrt( f2*h2 ) )
262+
s = ( fs / sqrt( f2 ) ) * ( conjg( gs ) / sqrt( h2 ) )
263263
else
264264
s = conjg( gs ) * ( r / h2 )
265265
end if
@@ -279,7 +279,7 @@ subroutine ZLARTG( f, g, c, s, r )
279279
! sqrt(safmin) <= f2 * sqrt(safmax) <= h2 / sqrt(f2 * h2) <= h2 * (safmin / f2) <= h2 <= safmax
280280
r = fs * ( h2 / d )
281281
end if
282-
s = conjg( gs ) * ( fs / d )
282+
s = ( fs / sqrt( f2 ) ) * ( conjg( gs ) / sqrt( h2 ) )
283283
end if
284284
! Rescale c and r
285285
c = c * w

0 commit comments

Comments
 (0)