home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD1.mdf / fortran / library / ssp / spclfnct / gmmma.for < prev    next >
Encoding:
Text File  |  1986-01-03  |  1.7 KB  |  70 lines

  1. C
  2. C     ..................................................................
  3. C
  4. C        SUBROUTINE GMMMA
  5. C
  6. C        PURPOSE
  7. C           COMPUTES THE GAMMA FUNCTION FOR A GIVEN ARGUMENT
  8. C
  9. C        USAGE
  10. C           CALL GMMMA(XX,GX,IER)
  11. C
  12. C        DESCRIPTION OF PARAMETERS
  13. C           XX -THE ARGUMENT FOR THE GAMMA FUNCTION
  14. C           GX -THE RESULTANT GAMMA FUNCTION VALUE
  15. C           IER-RESULTANT ERROR CODE WHERE
  16. C               IER=0  NO ERROR
  17. C               IER=1  XX IS WITHIN .000001 OF BEING A NEGATIVE INTEGER
  18. C               IER=2  XX GT 57, OVERFLOW, GX SET TO 1.0E75
  19. C
  20. C        REMARKS
  21. C           NONE
  22. C
  23. C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
  24. C           NONE
  25. C
  26. C        METHOD
  27. C           THE RECURSION RELATION AND POLYNOMIAL APPROXIMATION
  28. C           BY C.HASTINGS,JR., 'APPROXIMATIONS FOR DIGITAL COMPUTERS',
  29. C           PRINCETON UNIVERSITY PRESS, 1955
  30. C
  31. C     ..................................................................
  32. C
  33.       SUBROUTINE GMMMA(XX,GX,IER)
  34.       IF(XX-57.)6,6,4
  35.     4 IER=2
  36.       GX=1.E38
  37.       RETURN
  38.     6 X=XX
  39.       ERR=1.0E-6
  40.       IER=0
  41.       GX=1.0
  42.       IF(X-2.0)50,50,15
  43.    10 IF(X-2.0)110,110,15
  44.    15 X=X-1.0
  45.       GX=GX*X
  46.       GO TO 10
  47.    50 IF(X-1.0)60,120,110
  48. C
  49. C        SEE IF X IS NEAR NEGATIVE INTEGER OR ZERO
  50. C
  51.    60 IF(X-ERR)62,62,80
  52.    62 Y=FLOAT(INT(X))-X
  53.       IF(ABS(Y)-ERR)130,130,64
  54.    64 IF(1.0-Y-ERR)130,130,70
  55. C
  56. C        X NOT NEAR A NEGATIVE INTEGER OR ZERO
  57. C
  58.    70 IF(X-1.0)80,80,110
  59.    80 GX=GX/X
  60.       X=X+1.0
  61.       GO TO 70
  62.   110 Y=X-1.0
  63.       GY=1.0+Y*(-0.5771017+Y*(+0.9858540+Y*(-0.8764218+Y*(+0.8328212+
  64.      1Y*(-0.5684729+Y*(+0.2548205+Y*(-0.05149930)))))))
  65.       GX=GX*GY
  66.   120 RETURN
  67.   130 IER=1
  68.       RETURN
  69.       END
  70.