1184
1185
1186
1187
1188
1189
1190
1191 LOGICAL RND
1192 INTEGER BETA, EMAX, EMIN, T
1193 REAL EPS, RMAX, RMIN
1194
1195
1196
1197
1198
1199
1200
1201
1202
1203
1204
1205
1206
1207
1208
1209
1210
1211
1212
1213
1214
1215
1216
1217
1218
1219
1220
1221
1222
1223
1224
1225
1226
1227
1228
1229
1230
1231
1232
1233
1234
1235
1236
1237
1238
1239
1240
1241
1242
1243
1244
1245
1246
1247
1248
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
1256 REAL SLAMC3
1258
1259
1261
1262
1264
1265
1266 SAVE first, iwarn, lbeta, lemax, lemin, leps, lrmax,
1267 $ lrmin, lt
1268
1269
1270 DATA first / .true. / , iwarn / .false. /
1271
1272
1273
1274 IF( first ) THEN
1275 first = .false.
1276 zero = 0
1277 one = 1
1278 two = 2
1279
1280
1281
1282
1283
1284
1285
1286
1287
1288
1289 CALL slamc1( lbeta, lt, lrnd, lieee1 )
1290
1291
1292
1293 b = lbeta
1294 a = b**( -lt )
1295 leps = a
1296
1297
1298
1299 b = two / 3
1300 half = one / 2
1301 sixth =
slamc3( b, -half )
1302 third =
slamc3( sixth, sixth )
1303 b =
slamc3( third, -half )
1305 b = abs( b )
1306 IF( b.LT.leps )
1307 $ b = leps
1308
1309 leps = 1
1310
1311
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 ) )
1320 GO TO 10
1321 END IF
1322
1323
1324 IF( a.LT.leps )
1325 $ leps = a
1326
1327
1328
1329
1330
1331
1332
1333 rbase = one / lbeta
1334 small = one
1335 DO 20 i = 1, 3
1336 small =
slamc3( small*rbase, zero )
1337 20 CONTINUE
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
1349
1350 ELSE IF( ( gpmin-ngpmin ).EQ.3 ) THEN
1351 lemin = ngpmin - 1 + lt
1352 ieee = .true.
1353
1354
1355 ELSE
1356 lemin =
min( ngpmin, gpmin )
1357
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
1365
1366 ELSE
1367 lemin =
min( ngpmin, ngnmin )
1368
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
1377
1378 ELSE
1379 lemin =
min( ngpmin, ngnmin )
1380
1381 iwarn = .true.
1382 END IF
1383
1384 ELSE
1385 lemin =
min( ngpmin, ngnmin, gpmin, gnmin )
1386
1387 iwarn = .true.
1388 END IF
1389
1390
1391 IF( iwarn ) THEN
1392 first = .true.
1393 WRITE( 6, fmt = 9999 )lemin
1394 END IF
1395
1396
1397
1398
1399
1400
1401
1402 ieee = ieee .OR. lieee1
1403
1404
1405
1406
1407
1408 lrmin = 1
1409 DO 30 i = 1, 1 - lemin
1410 lrmin =
slamc3( lrmin*rbase, zero )
1411 30 CONTINUE
1412
1413
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
1437