LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ sb1nrm2()

 subroutine sb1nrm2 ( integer n, integer incx, real thresh )

Definition at line 1088 of file sblat1.f.

1089* Compare NRM2 with a reference computation using combinations
1090* of the following values:
1091*
1092* 0, very small, small, ulp, 1, 1/ulp, big, very big, infinity, NaN
1093*
1094* one of these values is used to initialize x(1) and x(2:N) is
1095* filled with random values from [-1,1] scaled by another of
1096* these values.
1097*
1098* This routine is adapted from the test suite provided by
1099* Anderson E. (2017)
1100* Algorithm 978: Safe Scaling in the Level 1 BLAS
1101* ACM Trans Math Softw 44:1--28
1102* https://doi.org/10.1145/3061665
1103*
1104 IMPLICIT NONE
1105* .. Scalar Arguments ..
1106 INTEGER INCX, N
1107 REAL THRESH
1108*
1109* =====================================================================
1110* .. Parameters ..
1111 INTEGER NMAX, NOUT, NV
1112 parameter(nmax=20, nout=6, nv=10)
1113 REAL HALF, ONE, TWO, ZERO
1114 parameter(half=0.5e+0, one=1.0e+0, two= 2.0e+0,
1115 & zero=0.0e+0)
1116* .. External Functions ..
1117 REAL SNRM2
1118 EXTERNAL snrm2
1119* .. Intrinsic Functions ..
1120 INTRINSIC abs, max, min, real, sqrt
1121* .. Model parameters ..
1122 REAL BIGNUM, SAFMAX, SAFMIN, SMLNUM, ULP
1123 parameter(bignum=0.1014120480e+32,
1124 & safmax=0.8507059173e+38,
1125 & safmin=0.1175494351e-37,
1126 & smlnum=0.9860761315e-31,
1127 & ulp=0.1192092896e-06)
1128* .. Local Scalars ..
1129 REAL ROGUE, SNRM, TRAT, V0, V1, WORKSSQ, Y1, Y2,
1130 & YMAX, YMIN, YNRM, ZNRM
1131 INTEGER I, IV, IW, IX
1132 LOGICAL FIRST
1133* .. Local Arrays ..
1134 REAL VALUES(NV), WORK(NMAX), X(NMAX), Z(NMAX)
1135* .. Executable Statements ..
1136 values(1) = zero
1137 values(2) = two*safmin
1138 values(3) = smlnum
1139 values(4) = ulp
1140 values(5) = one
1141 values(6) = one / ulp
1142 values(7) = bignum
1143 values(8) = safmax
1144 values(9) = sxvals(v0,2)
1145 values(10) = sxvals(v0,3)
1146 rogue = -1234.5678e+0
1147 first = .true.
1148*
1149* Check that the arrays are large enough
1150*
1151 IF (n*abs(incx).GT.nmax) THEN
1152 WRITE (nout,99) "SNRM2", nmax, incx, n, n*abs(incx)
1153 RETURN
1154 END IF
1155*
1156* Zero-sized inputs are tested in STEST1.
1157 IF (n.LE.0) THEN
1158 RETURN
1159 END IF
1160*
1161* Generate (N-1) values in (-1,1).
1162*
1163 DO i = 2, n
1164 CALL random_number(work(i))
1165 work(i) = one - two*work(i)
1166 END DO
1167*
1168* Compute the sum of squares of the random values
1169* by an unscaled algorithm.
1170*
1171 workssq = zero
1172 DO i = 2, n
1173 workssq = workssq + work(i)*work(i)
1174 END DO
1175*
1176* Construct the test vector with one known value
1177* and the rest from the random work array multiplied
1178* by a scaling factor.
1179*
1180 DO iv = 1, nv
1181 v0 = values(iv)
1182 IF (abs(v0).GT.one) THEN
1183 v0 = v0*half
1184 END IF
1185 z(1) = v0
1186 DO iw = 1, nv
1187 v1 = values(iw)
1188 IF (abs(v1).GT.one) THEN
1189 v1 = (v1*half) / sqrt(real(n))
1190 END IF
1191 DO i = 2, n
1192 z(i) = v1*work(i)
1193 END DO
1194*
1195* Compute the expected value of the 2-norm
1196*
1197 y1 = abs(v0)
1198 IF (n.GT.1) THEN
1199 y2 = abs(v1)*sqrt(workssq)
1200 ELSE
1201 y2 = zero
1202 END IF
1203 ymin = min(y1, y2)
1204 ymax = max(y1, y2)
1205*
1206* Expected value is NaN if either is NaN. The test
1207* for YMIN == YMAX avoids further computation if both
1208* are infinity.
1209*
1210 IF ((y1.NE.y1).OR.(y2.NE.y2)) THEN
1212 ynrm = y1 + y2
1213 ELSE IF (ymin == ymax) THEN
1214 ynrm = sqrt(two)*ymax
1215 ELSE IF (ymax == zero) THEN
1216 ynrm = zero
1217 ELSE
1218 ynrm = ymax*sqrt(one + (ymin / ymax)**2)
1219 END IF
1220*
1221* Fill the input array to SNRM2 with steps of incx
1222*
1223 DO i = 1, n
1224 x(i) = rogue
1225 END DO
1226 ix = 1
1227 IF (incx.LT.0) ix = 1 - (n-1)*incx
1228 DO i = 1, n
1229 x(ix) = z(i)
1230 ix = ix + incx
1231 END DO
1232*
1233* Call SNRM2 to compute the 2-norm
1234*
1235 snrm = snrm2(n,x,incx)
1236*
1237* Compare SNRM and ZNRM. Roundoff error grows like O(n)
1238* in this implementation so we scale the test ratio accordingly.
1239*
1240 IF (incx.EQ.0) THEN
1241 znrm = sqrt(real(n))*abs(x(1))
1242 ELSE
1243 znrm = ynrm
1244 END IF
1245*
1246* The tests for NaN rely on the compiler not being overly
1247* aggressive and removing the statements altogether.
1248 IF ((snrm.NE.snrm).OR.(znrm.NE.znrm)) THEN
1249 IF ((snrm.NE.snrm).NEQV.(znrm.NE.znrm)) THEN
1250 trat = one / ulp
1251 ELSE
1252 trat = zero
1253 END IF
1254 ELSE IF (snrm == znrm) THEN
1255 trat = zero
1256 ELSE IF (znrm == zero) THEN
1257 trat = snrm / ulp
1258 ELSE
1259 trat = (abs(snrm-znrm) / znrm) / (real(n)*ulp)
1260 END IF
1261 IF ((trat.NE.trat).OR.(trat.GE.thresh)) THEN
1262 IF (first) THEN
1263 first = .false.
1264 WRITE(nout,99999)
1265 END IF
1266 WRITE (nout,98) "SNRM2", n, incx, iv, iw, trat
1267 END IF
1268 END DO
1269 END DO
127099999 FORMAT (' FAIL')
1271 99 FORMAT ( ' Not enough space to test ', a6, ': NMAX = ',i6,
1272 + ', INCX = ',i6,/,' N = ',i6,', must be at least ',i6 )
1273 98 FORMAT( 1x, a6, ': N=', i6,', INCX=', i4, ', IV=', i2, ', IW=',
1274 + i2, ', test=', e15.8 )
1275 RETURN
1276 CONTAINS
1277 REAL FUNCTION SXVALS(XX,K)
1278* .. Scalar Arguments ..
1279 REAL XX
1280 INTEGER K
1281* .. Local Scalars ..
1282 REAL X, Y, YY, Z
1283* .. Intrinsic Functions ..
1284 INTRINSIC huge
1285* .. Executable Statements ..
1286 y = huge(xx)
1287 z = yy
1288 IF (k.EQ.1) THEN
1289 x = -z
1290 ELSE IF (k.EQ.2) THEN
1291 x = z
1292 ELSE IF (k.EQ.3) THEN
1293 x = z / z
1294 END IF
1295 sxvals = x
1296 RETURN
1297 END
real(wp) function snrm2(n, x, incx)
SNRM2
Definition snrm2.f90:89
Here is the caller graph for this function: