319
320
321
322
323
324
325
326 LOGICAL RND
327 INTEGER BETA, EMAX, EMIN, T
328 REAL EPS, RMAX, RMIN
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384 LOGICAL FIRST, IEEE, IWARN, LIEEE1, LRND
385 INTEGER GNMIN, GPMIN, I, LBETA, LEMAX, LEMIN, LT,
386 $ NGNMIN, NGPMIN
387 REAL A, B, C, HALF, LEPS, LRMAX, LRMIN, ONE, RBASE,
388 $ SIXTH, SMALL, THIRD, TWO, ZERO
389
390
391 REAL SLAMC3
393
394
396
397
399
400
401 SAVE first, iwarn, lbeta, lemax, lemin, leps, lrmax,
402 $ lrmin, lt
403
404
405 DATA first / .true. / , iwarn / .false. /
406
407
408
409 IF( first ) THEN
410 first = .false.
411 zero = 0
412 one = 1
413 two = 2
414
415
416
417
418
419
420
421
422
423
424 CALL slamc1( lbeta, lt, lrnd, lieee1 )
425
426
427
428 b = lbeta
429 a = b**( -lt )
430 leps = a
431
432
433
434 b = two / 3
435 half = one / 2
436 sixth =
slamc3( b, -half )
437 third =
slamc3( sixth, sixth )
438 b =
slamc3( third, -half )
440 b = abs( b )
441 IF( b.LT.leps )
442 $ b = leps
443
444 leps = 1
445
446
447 10 CONTINUE
448 IF( ( leps.GT.b ) .AND. ( b.GT.zero ) ) THEN
449 leps = b
450 c =
slamc3( half*leps, ( two**5 )*( leps**2 ) )
455 GO TO 10
456 END IF
457
458
459 IF( a.LT.leps )
460 $ leps = a
461
462
463
464
465
466
467
468 rbase = one / lbeta
469 small = one
470 DO 20 i = 1, 3
471 small =
slamc3( small*rbase, zero )
472 20 CONTINUE
474 CALL slamc4( ngpmin, one, lbeta )
475 CALL slamc4( ngnmin, -one, lbeta )
476 CALL slamc4( gpmin, a, lbeta )
477 CALL slamc4( gnmin, -a, lbeta )
478 ieee = .false.
479
480 IF( ( ngpmin.EQ.ngnmin ) .AND. ( gpmin.EQ.gnmin ) ) THEN
481 IF( ngpmin.EQ.gpmin ) THEN
482 lemin = ngpmin
483
484
485 ELSE IF( ( gpmin-ngpmin ).EQ.3 ) THEN
486 lemin = ngpmin - 1 + lt
487 ieee = .true.
488
489
490 ELSE
491 lemin =
min( ngpmin, gpmin )
492
493 iwarn = .true.
494 END IF
495
496 ELSE IF( ( ngpmin.EQ.gpmin ) .AND. ( ngnmin.EQ.gnmin ) ) THEN
497 IF( abs( ngpmin-ngnmin ).EQ.1 ) THEN
498 lemin =
max( ngpmin, ngnmin )
499
500
501 ELSE
502 lemin =
min( ngpmin, ngnmin )
503
504 iwarn = .true.
505 END IF
506
507 ELSE IF( ( abs( ngpmin-ngnmin ).EQ.1 ) .AND.
508 $ ( gpmin.EQ.gnmin ) ) THEN
509 IF( ( gpmin-
min( ngpmin, ngnmin ) ).EQ.3 )
THEN
510 lemin =
max( ngpmin, ngnmin ) - 1 + lt
511
512
513 ELSE
514 lemin =
min( ngpmin, ngnmin )
515
516 iwarn = .true.
517 END IF
518
519 ELSE
520 lemin =
min( ngpmin, ngnmin, gpmin, gnmin )
521
522 iwarn = .true.
523 END IF
524
525
526 IF( iwarn ) THEN
527 first = .true.
528 WRITE( 6, fmt = 9999 )lemin
529 END IF
530
531
532
533
534
535
536
537 ieee = ieee .OR. lieee1
538
539
540
541
542
543 lrmin = 1
544 DO 30 i = 1, 1 - lemin
545 lrmin =
slamc3( lrmin*rbase, zero )
546 30 CONTINUE
547
548
549
550 CALL slamc5( lbeta, lt, lemin, ieee, lemax, lrmax )
551 END IF
552
553 beta = lbeta
554 t = lt
555 rnd = lrnd
556 eps = leps
557 emin = lemin
558 rmin = lrmin
559 emax = lemax
560 rmax = lrmax
561
562 RETURN
563
564 9999 FORMAT( / / ' WARNING. The value EMIN may be incorrect:-',
565 $ ' EMIN = ', i8, /
566 $ ' If, after inspection, the value EMIN looks',
567 $ ' acceptable please comment out ',
568 $ / ' the IF block as marked within the code of routine',
569 $ ' SLAMC2,', / ' otherwise supply EMIN explicitly.', / )
570
571
572