phynexis-deps/lapack
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Fix DROTMG/SROTMG rescaling overwriting DH12

Browse Source 
When DFLAG=0 and both DD1 and DD2 need rescaling, the first
scaling loop transitions DFLAG from 0 to -1 and correctly
scales DH11 and DH12. But the second scaling loop then hits
the ELSE branch (matching DFLAG=-1) which unconditionally
resets DH21=-1 and DH12=1, overwriting the scaled DH12.

Change ELSE to ELSE IF (DFLAG.EQ.ONE) so the implied-element
initialization (DH21=-1, DH12=1) only fires when DFLAG=1,
not when DFLAG=-1 (where all elements are already explicit).

Fixes #244
This commit is contained in:
Julien Schueller
2026-06-11 15:05:19 +02:00
parent 526a4fbbd3
commit 6dc1a78eae
4 changed files with 102 additions and 100 deletions
Download Patch File Download Diff File Expand all files Collapse all files
BLAS/SRC/drotmg.f
+40 -46
View File
@@ -193,54 +193,48 @@
END IF
END IF
* Normalize flag to -ONE so all H elements are explicit
*
IF( DFLAG.EQ.ZERO ) THEN
DH11 = ONE
DH22 = ONE
DFLAG = -ONE
ELSE IF( DFLAG.EQ.ONE ) THEN
DH21 = -ONE
DH12 = ONE
DFLAG = -ONE
END IF
*
* PROCEDURE..SCALE-CHECK
IF (DD1.NE.ZERO) THEN
DO WHILE ((DD1.LE.RGAMSQ) .OR. (DD1.GE.GAMSQ))
IF (DFLAG.EQ.ZERO) THEN
DH11 = ONE
DH22 = ONE
DFLAG = -ONE
ELSE
DH21 = -ONE
DH12 = ONE
DFLAG = -ONE
END IF
IF (DD1.LE.RGAMSQ) THEN
DD1 = DD1*GAM**2
DX1 = DX1/GAM
DH11 = DH11/GAM
DH12 = DH12/GAM
ELSE
DD1 = DD1/GAM**2
DX1 = DX1*GAM
DH11 = DH11*GAM
DH12 = DH12*GAM
END IF
ENDDO
END IF
IF (DD1.NE.ZERO) THEN
DO WHILE ((DD1.LE.RGAMSQ) .OR. (DD1.GE.GAMSQ))
IF (DD1.LE.RGAMSQ) THEN
DD1 = DD1*GAM**2
DX1 = DX1/GAM
DH11 = DH11/GAM
DH12 = DH12/GAM
ELSE
DD1 = DD1/GAM**2
DX1 = DX1*GAM
DH11 = DH11*GAM
DH12 = DH12*GAM
END IF
ENDDO
END IF
IF (DD2.NE.ZERO) THEN
DO WHILE ( (DABS(DD2).LE.RGAMSQ) .OR. (DABS(DD2).GE.GAMSQ) )
IF (DFLAG.EQ.ZERO) THEN
DH11 = ONE
DH22 = ONE
DFLAG = -ONE
ELSE
DH21 = -ONE
DH12 = ONE
DFLAG = -ONE
END IF
IF (DABS(DD2).LE.RGAMSQ) THEN
DD2 = DD2*GAM**2
DH21 = DH21/GAM
DH22 = DH22/GAM
ELSE
DD2 = DD2/GAM**2
DH21 = DH21*GAM
DH22 = DH22*GAM
END IF
END DO
END IF
IF (DD2.NE.ZERO) THEN
DO WHILE ((DABS(DD2).LE.RGAMSQ) .OR. (DABS(DD2).GE.GAMSQ))
IF (DABS(DD2).LE.RGAMSQ) THEN
DD2 = DD2*GAM**2
DH21 = DH21/GAM
DH22 = DH22/GAM
ELSE
DD2 = DD2/GAM**2
DH21 = DH21*GAM
DH22 = DH22*GAM
END IF
END DO
END IF
END IF
BLAS/SRC/srotmg.f
+40 -46
View File
@@ -193,54 +193,48 @@
END IF
END IF
* Normalize flag to -ONE so all H elements are explicit
*
IF( SFLAG.EQ.ZERO ) THEN
SH11 = ONE
SH22 = ONE
SFLAG = -ONE
ELSE IF( SFLAG.EQ.ONE ) THEN
SH21 = -ONE
SH12 = ONE
SFLAG = -ONE
END IF
*
* PROCEDURE..SCALE-CHECK
IF (SD1.NE.ZERO) THEN
DO WHILE ((SD1.LE.RGAMSQ) .OR. (SD1.GE.GAMSQ))
IF (SFLAG.EQ.ZERO) THEN
SH11 = ONE
SH22 = ONE
SFLAG = -ONE
ELSE
SH21 = -ONE
SH12 = ONE
SFLAG = -ONE
END IF
IF (SD1.LE.RGAMSQ) THEN
SD1 = SD1*GAM**2
SX1 = SX1/GAM
SH11 = SH11/GAM
SH12 = SH12/GAM
ELSE
SD1 = SD1/GAM**2
SX1 = SX1*GAM
SH11 = SH11*GAM
SH12 = SH12*GAM
END IF
ENDDO
END IF
IF (SD1.NE.ZERO) THEN
DO WHILE ((SD1.LE.RGAMSQ) .OR. (SD1.GE.GAMSQ))
IF (SD1.LE.RGAMSQ) THEN
SD1 = SD1*GAM**2
SX1 = SX1/GAM
SH11 = SH11/GAM
SH12 = SH12/GAM
ELSE
SD1 = SD1/GAM**2
SX1 = SX1*GAM
SH11 = SH11*GAM
SH12 = SH12*GAM
END IF
ENDDO
END IF
IF (SD2.NE.ZERO) THEN
DO WHILE ( (ABS(SD2).LE.RGAMSQ) .OR. (ABS(SD2).GE.GAMSQ) )
IF (SFLAG.EQ.ZERO) THEN
SH11 = ONE
SH22 = ONE
SFLAG = -ONE
ELSE
SH21 = -ONE
SH12 = ONE
SFLAG = -ONE
END IF
IF (ABS(SD2).LE.RGAMSQ) THEN
SD2 = SD2*GAM**2
SH21 = SH21/GAM
SH22 = SH22/GAM
ELSE
SD2 = SD2/GAM**2
SH21 = SH21*GAM
SH22 = SH22*GAM
END IF
END DO
END IF
IF (SD2.NE.ZERO) THEN
DO WHILE ( (ABS(SD2).LE.RGAMSQ) .OR. (ABS(SD2).GE.GAMSQ) )
IF (ABS(SD2).LE.RGAMSQ) THEN
SD2 = SD2*GAM**2
SH21 = SH21/GAM
SH22 = SH22/GAM
ELSE
SD2 = SD2/GAM**2
SH21 = SH21*GAM
SH22 = SH22*GAM
END IF
END DO
END IF
END IF
BLAS/TESTING/dblat1.f
+11 -4
View File
@@ -172,7 +172,7 @@
E 4.D10, 2.D-2, 1.D-5, 10.D0,
F 2.D-10, 4.D-2, 1.D5, 10.D0,
G 2.D10, 4.D-2, 1.D-5, 10.D0,
H 4.D0, -2.D0, 8.D0, 4.D0 /
H 4.D-9, 2.D-9, 2.D0, 1.D0/
* TRUE RESULTS FOR MODIFIED GIVENS
DATA DTRUE/0.D0,0.D0, 1.3D0, .2D0, 0.D0,0.D0,0.D0, .5D0, 0.D0,
A 0.D0,0.D0, 4.5D0, 4.2D0, 1.D0, .5D0, 0.D0,0.D0,0.D0,
@@ -207,8 +207,15 @@
DTRUE(9,7) = 1.D4 / D12
DTRUE(1,8) = DTRUE(1,7)
DTRUE(2,8) = 2.D10 / (1.5D0 * D12 * D12)
DTRUE(1,9) = 32.D0 / 7.D0
DTRUE(2,9) = -16.D0 / 7.D0
DTRUE(1,9) = 5.9652323555555560D-02
DTRUE(2,9) = 2.9826161777777780D-02
DTRUE(3,9) = 5.4931640625000000D-04
DTRUE(4,9) = 1.D0
DTRUE(5,9) = -1.D0
DTRUE(6,9) = 2.4414062500000000D-04
DTRUE(7,9) = -1.2207031250000000D-04
DTRUE(8,9) = 6.1035156250000000D-05
DTRUE(9,9) = 2.4414062500000000D-04
* .. Executable Statements ..
*
* Compute true values which cannot be prestored
@@ -218,7 +225,7 @@
DBTRUE(3) = -1.0D0/0.6D0
DBTRUE(5) = 1.0D0/0.6D0
*
DO 20 K = 1, 8
DO 20 K = 1, 9
* .. Set N=K for identification in output if any ..
N = K
IF (ICASE.EQ.3) THEN
BLAS/TESTING/sblat1.f
+11 -4
View File
@@ -171,7 +171,7 @@
E 4.E10, 2.E-2, 1.E-5, 10.E0,
F 2.E-10, 4.E-2, 1.E5, 10.E0,
G 2.E10, 4.E-2, 1.E-5, 10.E0,
H 4.E0, -2.E0, 8.E0, 4.E0 /
H 4.E-9, 2.E-9, 2.E0, 1.E0/
* TRUE RESULTS FOR MODIFIED GIVENS
DATA DTRUE/0.E0,0.E0, 1.3E0, .2E0, 0.E0,0.E0,0.E0, .5E0, 0.E0,
A 0.E0,0.E0, 4.5E0, 4.2E0, 1.E0, .5E0, 0.E0,0.E0,0.E0,
@@ -206,8 +206,15 @@
DTRUE(9,7) = 1.E4 / D12
DTRUE(1,8) = DTRUE(1,7)
DTRUE(2,8) = 2.E10 / (1.5E0 * D12 * D12)
DTRUE(1,9) = 32.E0 / 7.E0
DTRUE(2,9) = -16.E0 / 7.E0
DTRUE(1,9) = 5.9652323555555560E-02
DTRUE(2,9) = 2.9826161777777780E-02
DTRUE(3,9) = 5.4931640625000000E-04
DTRUE(4,9) = 1.E0
DTRUE(5,9) = -1.E0
DTRUE(6,9) = 2.4414062500000000E-04
DTRUE(7,9) = -1.2207031250000000E-04
DTRUE(8,9) = 6.1035156250000000E-05
DTRUE(9,9) = 2.4414062500000000E-04
* .. Executable Statements ..
*
* Compute true values which cannot be prestored
@@ -217,7 +224,7 @@
DBTRUE(3) = -1.0E0/0.6E0
DBTRUE(5) = 1.0E0/0.6E0
*
DO 20 K = 1, 8
DO 20 K = 1, 9
* .. Set N=K for identification in output if any ..
N = K
IF (ICASE.EQ.3) THEN