Skip to content

Commit 8d23489

Browse files
julielangoujulielangou
authored
Merge pull request #290 from mgates3/norms
Update Frobenius norms for better accuracy.
2 parents c3b03d8 + 0d9289b commit 8d23489

Some content is hidden

Large Commits have some content hidden by default. Use the searchbox below for content that may be hidden.

46 files changed

+1458
-496
lines changed

‎SRC/CMakeLists.txt

Lines changed: 4 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -150,7 +150,8 @@ set(SLASRC
150150
stplqt.f stplqt2.f stpmlqt.f
151151
ssytrd_2stage.f ssytrd_sy2sb.f ssytrd_sb2st.F ssb2st_kernels.f
152152
ssyevd_2stage.f ssyev_2stage.f ssyevx_2stage.f ssyevr_2stage.f
153-
ssbev_2stage.f ssbevx_2stage.f ssbevd_2stage.f ssygv_2stage.f)
153+
ssbev_2stage.f ssbevx_2stage.f ssbevd_2stage.f ssygv_2stage.f
154+
scombssq.f)
154155

155156
set(DSLASRC spotrs.f sgetrs.f spotrf.f sgetrf.f)
156157

@@ -341,7 +342,8 @@ set(DLASRC
341342
dtplqt.f dtplqt2.f dtpmlqt.f
342343
dsytrd_2stage.f dsytrd_sy2sb.f dsytrd_sb2st.F dsb2st_kernels.f
343344
dsyevd_2stage.f dsyev_2stage.f dsyevx_2stage.f dsyevr_2stage.f
344-
dsbev_2stage.f dsbevx_2stage.f dsbevd_2stage.f dsygv_2stage.f)
345+
dsbev_2stage.f dsbevx_2stage.f dsbevd_2stage.f dsygv_2stage.f
346+
dcombssq.f)
345347

346348
set(DXLASRC dgesvxx.f dgerfsx.f dla_gerfsx_extended.f dla_geamv.f
347349
dla_gercond.f dla_gerpvgrw.f dsysvxx.f dsyrfsx.f

‎SRC/Makefile

Lines changed: 2 additions & 2 deletions
Original file line numberDiff line numberDiff line change
@@ -174,7 +174,7 @@ SLASRC = \
174174
ssytrd_2stage.o ssytrd_sy2sb.o ssytrd_sb2st.o ssb2st_kernels.o \
175175
ssyevd_2stage.o ssyev_2stage.o ssyevx_2stage.o ssyevr_2stage.o \
176176
ssbev_2stage.o ssbevx_2stage.o ssbevd_2stage.o ssygv_2stage.o \
177-
sgesvdq.o
177+
sgesvdq.o scombssq.o
178178

179179
DSLASRC = spotrs.o sgetrs.o spotrf.o sgetrf.o
180180

@@ -373,7 +373,7 @@ DLASRC = \
373373
dsytrd_2stage.o dsytrd_sy2sb.o dsytrd_sb2st.o dsb2st_kernels.o \
374374
dsyevd_2stage.o dsyev_2stage.o dsyevx_2stage.o dsyevr_2stage.o \
375375
dsbev_2stage.o dsbevx_2stage.o dsbevd_2stage.o dsygv_2stage.o \
376-
dgesvdq.o
376+
dgesvdq.o dcombssq.o
377377

378378
ifdef USEXBLAS
379379
DXLASRC = dgesvxx.o dgerfsx.o dla_gerfsx_extended.o dla_geamv.o \

‎SRC/clangb.f

Lines changed: 17 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -130,6 +130,7 @@ REAL FUNCTION CLANGB( NORM, N, KL, KU, AB, LDAB,
130130
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
131131
* December 2016
132132
*
133+
IMPLICIT NONE
133134
* .. Scalar Arguments ..
134135
CHARACTER NORM
135136
INTEGER KL, KU, LDAB, N
@@ -147,14 +148,17 @@ REAL FUNCTION CLANGB( NORM, N, KL, KU, AB, LDAB,
147148
* ..
148149
* .. Local Scalars ..
149150
INTEGER I, J, K, L
150-
REAL SCALE, SUM, VALUE, TEMP
151+
REAL SUM, VALUE, TEMP
152+
* ..
153+
* .. Local Arrays ..
154+
REAL SSQ( 2 ), COLSSQ( 2 )
151155
* ..
152156
* .. External Functions ..
153157
LOGICAL LSAME, SISNAN
154158
EXTERNAL LSAME, SISNAN
155159
* ..
156160
* .. External Subroutines ..
157-
EXTERNAL CLASSQ
161+
EXTERNAL CLASSQ, SCOMBSSQ
158162
* ..
159163
* .. Intrinsic Functions ..
160164
INTRINSIC ABS, MAX, MIN, SQRT
@@ -207,15 +211,22 @@ REAL FUNCTION CLANGB( NORM, N, KL, KU, AB, LDAB,
207211
ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
208212
*
209213
* Find normF(A).
214+
* SSQ(1) is scale
215+
* SSQ(2) is sum-of-squares
216+
* For better accuracy, sum each column separately.
210217
*
211-
SCALE = ZERO
212-
SUM = ONE
218+
SSQ( 1 ) = ZERO
219+
SSQ( 2 ) = ONE
213220
DO 90 J = 1, N
214221
L = MAX( 1, J-KU )
215222
K = KU + 1 - J + L
216-
CALL CLASSQ( MIN( N, J+KL )-L+1, AB( K, J ), 1, SCALE, SUM )
223+
COLSSQ( 1 ) = ZERO
224+
COLSSQ( 2 ) = ONE
225+
CALL CLASSQ( MIN( N, J+KL )-L+1, AB( K, J ), 1,
226+
$ COLSSQ( 1 ), COLSSQ( 2 ) )
227+
CALL SCOMBSSQ( SSQ, COLSSQ )
217228
90 CONTINUE
218-
VALUE = SCALE*SQRT( SUM )
229+
VALUE = SSQ( 1 )*SQRT( SSQ( 2 ) )
219230
END IF
220231
*
221232
CLANGB = VALUE

‎SRC/clange.f

Lines changed: 16 additions & 6 deletions
Original file line numberDiff line numberDiff line change
@@ -120,6 +120,7 @@ REAL FUNCTION CLANGE( NORM, M, N, A, LDA, WORK )
120120
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
121121
* December 2016
122122
*
123+
IMPLICIT NONE
123124
* .. Scalar Arguments ..
124125
CHARACTER NORM
125126
INTEGER LDA, M, N
@@ -137,14 +138,17 @@ REAL FUNCTION CLANGE( NORM, M, N, A, LDA, WORK )
137138
* ..
138139
* .. Local Scalars ..
139140
INTEGER I, J
140-
REAL SCALE, SUM, VALUE, TEMP
141+
REAL SUM, VALUE, TEMP
142+
* ..
143+
* .. Local Arrays ..
144+
REAL SSQ( 2 ), COLSSQ( 2 )
141145
* ..
142146
* .. External Functions ..
143147
LOGICAL LSAME, SISNAN
144148
EXTERNAL LSAME, SISNAN
145149
* ..
146150
* .. External Subroutines ..
147-
EXTERNAL CLASSQ
151+
EXTERNAL CLASSQ, SCOMBSSQ
148152
* ..
149153
* .. Intrinsic Functions ..
150154
INTRINSIC ABS, MIN, SQRT
@@ -196,13 +200,19 @@ REAL FUNCTION CLANGE( NORM, M, N, A, LDA, WORK )
196200
ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
197201
*
198202
* Find normF(A).
203+
* SSQ(1) is scale
204+
* SSQ(2) is sum-of-squares
205+
* For better accuracy, sum each column separately.
199206
*
200-
SCALE = ZERO
201-
SUM = ONE
207+
SSQ( 1 ) = ZERO
208+
SSQ( 2 ) = ONE
202209
DO 90 J = 1, N
203-
CALL CLASSQ( M, A( 1, J ), 1, SCALE, SUM )
210+
COLSSQ( 1 ) = ZERO
211+
COLSSQ( 2 ) = ONE
212+
CALL CLASSQ( M, A( 1, J ), 1, COLSSQ( 1 ), COLSSQ( 2 ) )
213+
CALL SCOMBSSQ( SSQ, COLSSQ )
204214
90 CONTINUE
205-
VALUE = SCALE*SQRT( SUM )
215+
VALUE = SSQ( 1 )*SQRT( SSQ( 2 ) )
206216
END IF
207217
*
208218
CLANGE = VALUE

‎SRC/clanhb.f

Lines changed: 35 additions & 13 deletions
Original file line numberDiff line numberDiff line change
@@ -137,6 +137,7 @@ REAL FUNCTION CLANHB( NORM, UPLO, N, K, AB, LDAB,
137137
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
138138
* December 2016
139139
*
140+
IMPLICIT NONE
140141
* .. Scalar Arguments ..
141142
CHARACTER NORM, UPLO
142143
INTEGER K, LDAB, N
@@ -154,14 +155,17 @@ REAL FUNCTION CLANHB( NORM, UPLO, N, K, AB, LDAB,
154155
* ..
155156
* .. Local Scalars ..
156157
INTEGER I, J, L
157-
REAL ABSA, SCALE, SUM, VALUE
158+
REAL ABSA, SUM, VALUE
159+
* ..
160+
* .. Local Arrays ..
161+
REAL SSQ( 2 ), COLSSQ( 2 )
158162
* ..
159163
* .. External Functions ..
160164
LOGICAL LSAME, SISNAN
161165
EXTERNAL LSAME, SISNAN
162166
* ..
163167
* .. External Subroutines ..
164-
EXTERNAL CLASSQ
168+
EXTERNAL CLASSQ, SCOMBSSQ
165169
* ..
166170
* .. Intrinsic Functions ..
167171
INTRINSIC ABS, MAX, MIN, REAL, SQRT
@@ -233,39 +237,57 @@ REAL FUNCTION CLANHB( NORM, UPLO, N, K, AB, LDAB,
233237
ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
234238
*
235239
* Find normF(A).
240+
* SSQ(1) is scale
241+
* SSQ(2) is sum-of-squares
242+
* For better accuracy, sum each column separately.
243+
*
244+
SSQ( 1 ) = ZERO
245+
SSQ( 2 ) = ONE
246+
*
247+
* Sum off-diagonals
236248
*
237-
SCALE = ZERO
238-
SUM = ONE
239249
IF( K.GT.0 ) THEN
240250
IF( LSAME( UPLO, 'U' ) ) THEN
241251
DO 110 J = 2, N
252+
COLSSQ( 1 ) = ZERO
253+
COLSSQ( 2 ) = ONE
242254
CALL CLASSQ( MIN( J-1, K ), AB( MAX( K+2-J, 1 ), J ),
243-
$ 1, SCALE, SUM )
255+
$ 1, COLSSQ( 1 ), COLSSQ( 2 ) )
256+
CALL SCOMBSSQ( SSQ, COLSSQ )
244257
110 CONTINUE
245258
L = K + 1
246259
ELSE
247260
DO 120 J = 1, N - 1
248-
CALL CLASSQ( MIN( N-J, K ), AB( 2, J ), 1, SCALE,
249-
$ SUM )
261+
COLSSQ( 1 ) = ZERO
262+
COLSSQ( 2 ) = ONE
263+
CALL CLASSQ( MIN( N-J, K ), AB( 2, J ), 1,
264+
$ COLSSQ( 1 ), COLSSQ( 2 ) )
265+
CALL SCOMBSSQ( SSQ, COLSSQ )
250266
120 CONTINUE
251267
L = 1
252268
END IF
253-
SUM = 2*SUM
269+
SSQ( 2 ) = 2*SSQ( 2 )
254270
ELSE
255271
L = 1
256272
END IF
273+
*
274+
* Sum diagonal
275+
*
276+
COLSSQ( 1 ) = ZERO
277+
COLSSQ( 2 ) = ONE
257278
DO 130 J = 1, N
258279
IF( REAL( AB( L, J ) ).NE.ZERO ) THEN
259280
ABSA = ABS( REAL( AB( L, J ) ) )
260-
IF( SCALE.LT.ABSA ) THEN
261-
SUM = ONE + SUM*( SCALE / ABSA )**2
262-
SCALE = ABSA
281+
IF( COLSSQ( 1 ).LT.ABSA ) THEN
282+
COLSSQ( 2 ) = ONE + COLSSQ(2)*( COLSSQ(1) / ABSA )**2
283+
COLSSQ( 1 ) = ABSA
263284
ELSE
264-
SUM = SUM + ( ABSA / SCALE )**2
285+
COLSSQ( 2 ) = COLSSQ( 2 ) + ( ABSA / COLSSQ( 1 ) )**2
265286
END IF
266287
END IF
267288
130 CONTINUE
268-
VALUE = SCALE*SQRT( SUM )
289+
CALL SCOMBSSQ( SSQ, COLSSQ )
290+
VALUE = SSQ( 1 )*SQRT( SSQ( 2 ) )
269291
END IF
270292
*
271293
CLANHB = VALUE

‎SRC/clanhe.f

Lines changed: 33 additions & 12 deletions
Original file line numberDiff line numberDiff line change
@@ -129,6 +129,7 @@ REAL FUNCTION CLANHE( NORM, UPLO, N, A, LDA, WORK )
129129
* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
130130
* December 2016
131131
*
132+
IMPLICIT NONE
132133
* .. Scalar Arguments ..
133134
CHARACTER NORM, UPLO
134135
INTEGER LDA, N
@@ -146,14 +147,17 @@ REAL FUNCTION CLANHE( NORM, UPLO, N, A, LDA, WORK )
146147
* ..
147148
* .. Local Scalars ..
148149
INTEGER I, J
149-
REAL ABSA, SCALE, SUM, VALUE
150+
REAL ABSA, SUM, VALUE
151+
* ..
152+
* .. Local Arrays ..
153+
REAL SSQ( 2 ), COLSSQ( 2 )
150154
* ..
151155
* .. External Functions ..
152156
LOGICAL LSAME, SISNAN
153157
EXTERNAL LSAME, SISNAN
154158
* ..
155159
* .. External Subroutines ..
156-
EXTERNAL CLASSQ
160+
EXTERNAL CLASSQ, SCOMBSSQ
157161
* ..
158162
* .. Intrinsic Functions ..
159163
INTRINSIC ABS, REAL, SQRT
@@ -223,31 +227,48 @@ REAL FUNCTION CLANHE( NORM, UPLO, N, A, LDA, WORK )
223227
ELSE IF( ( LSAME( NORM, 'F' ) ) .OR. ( LSAME( NORM, 'E' ) ) ) THEN
224228
*
225229
* Find normF(A).
230+
* SSQ(1) is scale
231+
* SSQ(2) is sum-of-squares
232+
* For better accuracy, sum each column separately.
233+
*
234+
SSQ( 1 ) = ZERO
235+
SSQ( 2 ) = ONE
236+
*
237+
* Sum off-diagonals
226238
*
227-
SCALE = ZERO
228-
SUM = ONE
229239
IF( LSAME( UPLO, 'U' ) ) THEN
230240
DO 110 J = 2, N
231-
CALL CLASSQ( J-1, A( 1, J ), 1, SCALE, SUM )
241+
COLSSQ( 1 ) = ZERO
242+
COLSSQ( 2 ) = ONE
243+
CALL CLASSQ( J-1, A( 1, J ), 1,
244+
$ COLSSQ( 1 ), COLSSQ( 2 ) )
245+
CALL SCOMBSSQ( SSQ, COLSSQ )
232246
110 CONTINUE
233247
ELSE
234248
DO 120 J = 1, N - 1
235-
CALL CLASSQ( N-J, A( J+1, J ), 1, SCALE, SUM )
249+
COLSSQ( 1 ) = ZERO
250+
COLSSQ( 2 ) = ONE
251+
CALL CLASSQ( N-J, A( J+1, J ), 1,
252+
$ COLSSQ( 1 ), COLSSQ( 2 ) )
253+
CALL SCOMBSSQ( SSQ, COLSSQ )
236254
120 CONTINUE
237255
END IF
238-
SUM = 2*SUM
256+
SSQ( 2 ) = 2*SSQ( 2 )
257+
*
258+
* Sum diagonal
259+
*
239260
DO 130 I = 1, N
240261
IF( REAL( A( I, I ) ).NE.ZERO ) THEN
241262
ABSA = ABS( REAL( A( I, I ) ) )
242-
IF( SCALE.LT.ABSA ) THEN
243-
SUM = ONE + SUM*( SCALE / ABSA )**2
244-
SCALE = ABSA
263+
IF( SSQ( 1 ).LT.ABSA ) THEN
264+
SSQ( 2 ) = ONE + SSQ( 2 )*( SSQ( 1 ) / ABSA )**2
265+
SSQ( 1 ) = ABSA
245266
ELSE
246-
SUM = SUM + ( ABSA / SCALE )**2
267+
SSQ( 2 ) = SSQ( 2 ) + ( ABSA / SSQ( 1 ) )**2
247268
END IF
248269
END IF
249270
130 CONTINUE
250-
VALUE = SCALE*SQRT( SUM )
271+
VALUE = SSQ( 1 )*SQRT( SSQ( 2 ) )
251272
END IF
252273
*
253274
CLANHE = VALUE

0 commit comments

Comments
 (0)