beskes.f 2.7 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677
  1. *DECK BESKES
  2. SUBROUTINE BESKES (XNU, X, NIN, BKE)
  3. C***BEGIN PROLOGUE BESKES
  4. C***PURPOSE Compute a sequence of exponentially scaled modified Bessel
  5. C functions of the third kind of fractional order.
  6. C***LIBRARY SLATEC (FNLIB)
  7. C***CATEGORY C10B3
  8. C***TYPE SINGLE PRECISION (BESKES-S, DBSKES-D)
  9. C***KEYWORDS EXPONENTIALLY SCALED, FNLIB, FRACTIONAL ORDER,
  10. C MODIFIED BESSEL FUNCTION, SEQUENCE OF BESSEL FUNCTIONS,
  11. C SPECIAL FUNCTIONS, THIRD KIND
  12. C***AUTHOR Fullerton, W., (LANL)
  13. C***DESCRIPTION
  14. C
  15. C BESKES computes a sequence of exponentially scaled
  16. C (i.e., multipled by EXP(X)) modified Bessel
  17. C functions of the third kind of order XNU + I at X, where X .GT. 0,
  18. C XNU lies in (-1,1), and I = 0, 1, ... , NIN - 1, if NIN is positive
  19. C and I = 0, -1, ... , NIN + 1, if NIN is negative. On return, the
  20. C vector BKE(.) contains the results at X for order starting at XNU.
  21. C
  22. C***REFERENCES (NONE)
  23. C***ROUTINES CALLED R1MACH, R9KNUS, XERMSG
  24. C***REVISION HISTORY (YYMMDD)
  25. C 770601 DATE WRITTEN
  26. C 890531 Changed all specific intrinsics to generic. (WRB)
  27. C 890911 Removed unnecessary intrinsics. (WRB)
  28. C 890911 REVISION DATE from Version 3.2
  29. C 891214 Prologue converted to Version 4.0 format. (BAB)
  30. C 900315 CALLs to XERROR changed to CALLs to XERMSG. (THJ)
  31. C 900326 Removed duplicate information from DESCRIPTION section.
  32. C (WRB)
  33. C***END PROLOGUE BESKES
  34. DIMENSION BKE(*)
  35. SAVE ALNBIG
  36. DATA ALNBIG / 0. /
  37. C***FIRST EXECUTABLE STATEMENT BESKES
  38. IF (ALNBIG.EQ.0.) ALNBIG = LOG (R1MACH(2))
  39. C
  40. V = ABS(XNU)
  41. N = ABS(NIN)
  42. C
  43. IF (V .GE. 1.) CALL XERMSG ('SLATEC', 'BESKES',
  44. + 'ABS(XNU) MUST BE LT 1', 2, 2)
  45. IF (X .LE. 0.) CALL XERMSG ('SLATEC', 'BESKES', 'X IS LE 0', 3,
  46. + 2)
  47. IF (N .EQ. 0) CALL XERMSG ('SLATEC', 'BESKES',
  48. + 'N THE NUMBER IN THE SEQUENCE IS 0', 4, 2)
  49. C
  50. CALL R9KNUS (V, X, BKE(1), BKNU1, ISWTCH)
  51. IF (N.EQ.1) RETURN
  52. C
  53. VINCR = SIGN (1.0, REAL(NIN))
  54. DIRECT = VINCR
  55. IF (XNU.NE.0.) DIRECT = VINCR*SIGN(1.0,XNU)
  56. IF (ISWTCH .EQ. 1 .AND. DIRECT .GT. 0.) CALL XERMSG ('SLATEC',
  57. + 'BESKES', 'X SO SMALL BESSEL K-SUB-XNU+1 OVERFLOWS', 5, 2)
  58. BKE(2) = BKNU1
  59. C
  60. IF (DIRECT.LT.0.) CALL R9KNUS (ABS(XNU+VINCR), X, BKE(2), BKNU1,
  61. 1 ISWTCH)
  62. IF (N.EQ.2) RETURN
  63. C
  64. VEND = ABS(XNU+NIN) - 1.0
  65. IF ((VEND-0.5)*LOG(VEND)+0.27-VEND*(LOG(X)-.694) .GT. ALNBIG)
  66. 1CALL XERMSG ( 'SLATEC', 'BESKES',
  67. 2'X SO SMALL OR ABS(NU) SO BIG THAT BESSEL K-SUB-NU OVERFLOWS',
  68. 35, 2)
  69. C
  70. V = XNU
  71. DO 10 I=3,N
  72. V = V + VINCR
  73. BKE(I) = 2.0*V*BKE(I-1)/X + BKE(I-2)
  74. 10 CONTINUE
  75. C
  76. RETURN
  77. END