SCALAPACK 2.2.2
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ slamc2()

subroutine slamc2 ( integer  beta,
integer  t,
logical  rnd,
real  eps,
integer  emin,
real  rmin,
integer  emax,
real  rmax 
)

Definition at line 1183 of file tools.f.

1184*
1185* -- LAPACK auxiliary routine (version 2.0) --
1186* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd.,
1187* Courant Institute, Argonne National Lab, and Rice University
1188* October 31, 1992
1189*
1190* .. Scalar Arguments ..
1191 LOGICAL RND
1192 INTEGER BETA, EMAX, EMIN, T
1193 REAL EPS, RMAX, RMIN
1194* ..
1195*
1196* Purpose
1197* =======
1198*
1199* SLAMC2 determines the machine parameters specified in its argument
1200* list.
1201*
1202* Arguments
1203* =========
1204*
1205* BETA (output) INTEGER
1206* The base of the machine.
1207*
1208* T (output) INTEGER
1209* The number of ( BETA ) digits in the mantissa.
1210*
1211* RND (output) LOGICAL
1212* Specifies whether proper rounding ( RND = .TRUE. ) or
1213* chopping ( RND = .FALSE. ) occurs in addition. This may not
1214* be a reliable guide to the way in which the machine performs
1215* its arithmetic.
1216*
1217* EPS (output) REAL
1218* The smallest positive number such that
1219*
1220* fl( 1.0 - EPS ) .LT. 1.0,
1221*
1222* where fl denotes the computed value.
1223*
1224* EMIN (output) INTEGER
1225* The minimum exponent before (gradual) underflow occurs.
1226*
1227* RMIN (output) REAL
1228* The smallest normalized number for the machine, given by
1229* BASE**( EMIN - 1 ), where BASE is the floating point value
1230* of BETA.
1231*
1232* EMAX (output) INTEGER
1233* The maximum exponent before overflow occurs.
1234*
1235* RMAX (output) REAL
1236* The largest positive number for the machine, given by
1237* BASE**EMAX * ( 1 - EPS ), where BASE is the floating point
1238* value of BETA.
1239*
1240* Further Details
1241* ===============
1242*
1243* The computation of EPS is based on a routine PARANOIA by
1244* W. Kahan of the University of California at Berkeley.
1245*
1246* =====================================================================
1247*
1248* .. Local Scalars ..
1249 LOGICAL FIRST, IEEE, IWARN, LIEEE1, LRND
1250 INTEGER GNMIN, GPMIN, I, LBETA, LEMAX, LEMIN, LT,
1251 $ NGNMIN, NGPMIN
1252 REAL A, B, C, HALF, LEPS, LRMAX, LRMIN, ONE, RBASE,
1253 $ SIXTH, SMALL, THIRD, TWO, ZERO
1254* ..
1255* .. External Functions ..
1256 REAL SLAMC3
1257 EXTERNAL slamc3
1258* ..
1259* .. External Subroutines ..
1260 EXTERNAL slamc1, slamc4, slamc5
1261* ..
1262* .. Intrinsic Functions ..
1263 INTRINSIC abs, max, min
1264* ..
1265* .. Save statement ..
1266 SAVE first, iwarn, lbeta, lemax, lemin, leps, lrmax,
1267 $ lrmin, lt
1268* ..
1269* .. Data statements ..
1270 DATA first / .true. / , iwarn / .false. /
1271* ..
1272* .. Executable Statements ..
1273*
1274 IF( first ) THEN
1275 first = .false.
1276 zero = 0
1277 one = 1
1278 two = 2
1279*
1280* LBETA, LT, LRND, LEPS, LEMIN and LRMIN are the local values of
1281* BETA, T, RND, EPS, EMIN and RMIN.
1282*
1283* Throughout this routine we use the function SLAMC3 to ensure
1284* that relevant values are stored and not held in registers, or
1285* are not affected by optimizers.
1286*
1287* SLAMC1 returns the parameters LBETA, LT, LRND and LIEEE1.
1288*
1289 CALL slamc1( lbeta, lt, lrnd, lieee1 )
1290*
1291* Start to find EPS.
1292*
1293 b = lbeta
1294 a = b**( -lt )
1295 leps = a
1296*
1297* Try some tricks to see whether or not this is the correct EPS.
1298*
1299 b = two / 3
1300 half = one / 2
1301 sixth = slamc3( b, -half )
1302 third = slamc3( sixth, sixth )
1303 b = slamc3( third, -half )
1304 b = slamc3( b, sixth )
1305 b = abs( b )
1306 IF( b.LT.leps )
1307 $ b = leps
1308*
1309 leps = 1
1310*
1311*+ WHILE( ( LEPS.GT.B ).AND.( B.GT.ZERO ) )LOOP
1312 10 CONTINUE
1313 IF( ( leps.GT.b ) .AND. ( b.GT.zero ) ) THEN
1314 leps = b
1315 c = slamc3( half*leps, ( two**5 )*( leps**2 ) )
1316 c = slamc3( half, -c )
1317 b = slamc3( half, c )
1318 c = slamc3( half, -b )
1319 b = slamc3( half, c )
1320 GO TO 10
1321 END IF
1322*+ END WHILE
1323*
1324 IF( a.LT.leps )
1325 $ leps = a
1326*
1327* Computation of EPS complete.
1328*
1329* Now find EMIN. Let A = + or - 1, and + or - (1 + BASE**(-3)).
1330* Keep dividing A by BETA until (gradual) underflow occurs. This
1331* is detected when we cannot recover the previous A.
1332*
1333 rbase = one / lbeta
1334 small = one
1335 DO 20 i = 1, 3
1336 small = slamc3( small*rbase, zero )
1337 20 CONTINUE
1338 a = slamc3( one, small )
1339 CALL slamc4( ngpmin, one, lbeta )
1340 CALL slamc4( ngnmin, -one, lbeta )
1341 CALL slamc4( gpmin, a, lbeta )
1342 CALL slamc4( gnmin, -a, lbeta )
1343 ieee = .false.
1344*
1345 IF( ( ngpmin.EQ.ngnmin ) .AND. ( gpmin.EQ.gnmin ) ) THEN
1346 IF( ngpmin.EQ.gpmin ) THEN
1347 lemin = ngpmin
1348* ( Non twos-complement machines, no gradual underflow;
1349* e.g., VAX )
1350 ELSE IF( ( gpmin-ngpmin ).EQ.3 ) THEN
1351 lemin = ngpmin - 1 + lt
1352 ieee = .true.
1353* ( Non twos-complement machines, with gradual underflow;
1354* e.g., IEEE standard followers )
1355 ELSE
1356 lemin = min( ngpmin, gpmin )
1357* ( A guess; no known machine )
1358 iwarn = .true.
1359 END IF
1360*
1361 ELSE IF( ( ngpmin.EQ.gpmin ) .AND. ( ngnmin.EQ.gnmin ) ) THEN
1362 IF( abs( ngpmin-ngnmin ).EQ.1 ) THEN
1363 lemin = max( ngpmin, ngnmin )
1364* ( Twos-complement machines, no gradual underflow;
1365* e.g., CYBER 205 )
1366 ELSE
1367 lemin = min( ngpmin, ngnmin )
1368* ( A guess; no known machine )
1369 iwarn = .true.
1370 END IF
1371*
1372 ELSE IF( ( abs( ngpmin-ngnmin ).EQ.1 ) .AND.
1373 $ ( gpmin.EQ.gnmin ) ) THEN
1374 IF( ( gpmin-min( ngpmin, ngnmin ) ).EQ.3 ) THEN
1375 lemin = max( ngpmin, ngnmin ) - 1 + lt
1376* ( Twos-complement machines with gradual underflow;
1377* no known machine )
1378 ELSE
1379 lemin = min( ngpmin, ngnmin )
1380* ( A guess; no known machine )
1381 iwarn = .true.
1382 END IF
1383*
1384 ELSE
1385 lemin = min( ngpmin, ngnmin, gpmin, gnmin )
1386* ( A guess; no known machine )
1387 iwarn = .true.
1388 END IF
1389***
1390* Comment out this if block if EMIN is ok
1391 IF( iwarn ) THEN
1392 first = .true.
1393 WRITE( 6, fmt = 9999 )lemin
1394 END IF
1395***
1396*
1397* Assume IEEE arithmetic if we found denormalised numbers above,
1398* or if arithmetic seems to round in the IEEE style, determined
1399* in routine SLAMC1. A true IEEE machine should have both things
1400* true; however, faulty machines may have one or the other.
1401*
1402 ieee = ieee .OR. lieee1
1403*
1404* Compute RMIN by successive division by BETA. We could compute
1405* RMIN as BASE**( EMIN - 1 ), but some machines underflow during
1406* this computation.
1407*
1408 lrmin = 1
1409 DO 30 i = 1, 1 - lemin
1410 lrmin = slamc3( lrmin*rbase, zero )
1411 30 CONTINUE
1412*
1413* Finally, call SLAMC5 to compute EMAX and RMAX.
1414*
1415 CALL slamc5( lbeta, lt, lemin, ieee, lemax, lrmax )
1416 END IF
1417*
1418 beta = lbeta
1419 t = lt
1420 rnd = lrnd
1421 eps = leps
1422 emin = lemin
1423 rmin = lrmin
1424 emax = lemax
1425 rmax = lrmax
1426*
1427 RETURN
1428*
1429 9999 FORMAT( / / ' WARNING. The value EMIN may be incorrect:-',
1430 $ ' EMIN = ', i8, /
1431 $ ' If, after inspection, the value EMIN looks',
1432 $ ' acceptable please comment out ',
1433 $ / ' the IF block as marked within the code of routine',
1434 $ ' SLAMC2,', / ' otherwise supply EMIN explicitly.', / )
1435*
1436* End of SLAMC2
1437*
#define max(A, B)
Definition pcgemr.c:180
#define min(A, B)
Definition pcgemr.c:181
subroutine slamc5(beta, p, emin, ieee, emax, rmax)
Definition tools.f:1565
subroutine slamc1(beta, t, rnd, ieee1)
Definition tools.f:997
subroutine slamc4(emin, start, base)
Definition tools.f:1481
real function slamc3(a, b)
Definition tools.f:1443
Here is the call graph for this function:
Here is the caller graph for this function: