52 SAVE alneps, sqeps, bot, first
56 alneps = -
log(r1mach(3))
57 sqeps =
sqrt(r1mach(4))
62 IF (x .LT. 0.0) CALL
xermsg(
'SLATEC',
'GAMIT',
'X IS NEGATIVE',
65 IF (x.NE.0.0) alx =
log(x)
67 IF (a.NE.0.0) sga = sign(1.0, a)
68 ainta = aint(a+0.5*sga)
71 IF (x.GT.0.0) go to 20
73 IF (ainta.GT.0.0 .OR. aeps.NE.0.0)
gamit =
gamr(a+1.0)
76 20
IF (x.GT.1.0) go to 40
77 IF (a.GE.(-0.5) .OR. aeps.NE.0.0) CALL
algams(a+1.0, algap1,
82 40
IF (a.LT.x) go to 50
88 50 alng =
r9lgic(a, x, alx)
93 IF (aeps.EQ.0.0 .AND. ainta.LE.0.0) go to 60
94 CALL
algams(a+1.0, algap1, sgngam)
95 t =
log(
abs(a)) + alng - algap1
96 IF (t.GT.alneps) go to 70
97 IF (t.GT.(-alneps)) h = 1.0 - sga*sgngam*
exp(t)
98 IF (
abs(h).GT.sqeps) go to 60
100 CALL
xermsg(
'SLATEC',
'GAMIT',
'RESULT LT HALF PRECISION', 1, 1)
102 60 t = -a*alx +
log(
abs(h))
function r9lgit(A, X, ALGAP1)
function r9gmit(A, X, ALGAP1, SGNGAM, ALX)
subroutine algams(X, ALGAM, SGNGAM)
octave_value log(void) const
subroutine xermsg(LIBRAR, SUBROU, MESSG, NERR, LEVEL)
function r9lgic(A, X, ALX)
real function gamit(A, X)
octave_value sqrt(void) const