[CONTACT]

[ABOUT]

[POLICY]

[ADVERTISE]

SUBROUTINE M,X,F,WORK,

Found at: ftp.icm.edu.pl:70/packages/netlib/vfnlib/dvi0.f

      SUBROUTINE DVI0 (M, X, F, WORK, IWORK, INFO)
C***BEGIN PROLOGUE  DVI0
C***PURPOSE  Computes the hyperbolic Bessel function of the first kind
C            of order zero (I0) for a vector of real arguments
C***LIBRARY   VFNLIB
C***CATEGORY  C10B1
C***TYPE      DOUBLE PRECISION (VI0-S, DVI0-D)
C***KEYWORDS  BESSEL FUNCTION,FIRST KIND,HYPERBOLIC BESSEL FUNCTION,
C             MODIFIED BESSEL FUNCTION, ORDER ZERO, VECTORIZED
C***AUTHOR  SAUNDERS, B. V., (NIST)
C             COMPUTING AND APPLIED MATHEMATICS LABORATORY
C             NATIONAL INSTITUTE OF STANDARDS AND TECHNOLOGY
C             GAITHERSBURG, MD 20899
C           BOISVERT, R. F., (NIST)
C             COMPUTING AND APPLIED MATHEMATICS LABORATORY
C             NATIONAL INSTITUTE OF STANDARDS AND TECHNOLOGY
C             GAITHERSBURG, MD 20899
C***DESCRIPTION
C
C   DVI0 computes the modified (hyperbolic) Bessel function of the
C   first kind of order zero (I0) for real arguments using uniform
C   approximation by Chebyshev polynomials.
C
C
C   P A R A M E T E R S
C
C   M      (Input) Integer (M .GT. 0)
C          The number of arguments at which the function is to be
C          evaluated.
C
C   X      (Input) Double precision array of length M
C          The arguments at which the function is to be evaluated are
C          stored in X(1) to X(M) in any order.
C
C   F      (Output) Double precision array of length M
C          F(i) contains the value of the function at X(i), i=1,..,M.
C
C   WORK   (Work) Double precision vector of length 7*M 
C
C   IWORK  (Work) Integer vector of length M 
C
C   INFO   (Output) Integer
C          Indicates status of computed result. The following table
C          lists possible values and their meanings.  If OK=Yes then
C          all F(i) have been set, otherwise none have been set.
C
C          INFO  OK            Description
C          ------------------------------------------------------------
C            0   Yes  Successfull execution.
C            1   No   Error: M .LE. 0
C            2   No   Error: Some abs(X(i)) so big I0 overflows.
C                     The index of the first offending argument is 
C                     returned in IWORK(1).
C
C
C *********************************************************************
C   This routine is a modification of the function DBESI0 developed by
C   W. Fullerton of LANL.
C *********************************************************************
C
C***REFERENCES  Ronald F. Boisvert and Bonita V. Saunders, Portable
C               Vectorized Software for Bessel Function Evaluation,
C               ACM Transactions on Mathematical Software 18 (1992),
C               pp 456-469.
C***ROUTINES CALLED  DWI0 
C***REVISION HISTORY  (YYMMDD)
C   910226  DATE WRITTEN 
C   930203  Minor modifications to prologue.
C***END PROLOGUE  DVI0
C
C     ----------
C     PARAMETERS
C     ----------
C
      INTEGER INFO, IWORK, M
      DOUBLE PRECISION F, X, WORK
C
      DIMENSION X(M), F(M), WORK(7*M), IWORK(M)
C
C     ---------------
C     LOCAL VARIABLES
C     ---------------
C
      INTEGER IWB0, IWB1, IWB2, IWY, IWTC, IWYC, IWZC, JIN
C
C***FIRST EXECUTABLE STATEMENT  DVI0
C
C     ... PARTITION WORK ARRAYS
C
      IWY   = 1
      IWTC  = IWY  + M
      IWYC  = IWTC + M
      IWZC  = IWYC + M
      IWB0  = IWZC + M
      IWB1  = IWB0 + M
      IWB2  = IWB1 + M
C     Total = IWB2 + M
C
      JIN   = 1
C     Total = JIN  + M
C
C     ... DWI0 DOES ALL THE WORK
C
      CALL DWI0(M,X,F,WORK(IWY),WORK(IWTC),WORK(IWYC),WORK(IWZC),
     +         IWORK(JIN),WORK(IWB0),WORK(IWB1),WORK(IWB2),INFO)
C
      RETURN
      END
      SUBROUTINE DWI0 (M, X, F, Y, TCMP, YCMP, ZCMP, INDX, B0, B1, B2, 
     +   INFO)
C***BEGIN PROLOGUE  DWI0
C***SUBSIDIARY
C***PURPOSE  Computes the hyperbolic Bessel function of the first kind
C            of order zero (I0) for a vector of arguments
C***LIBRARY   VFNLIB
C***AUTHOR  SAUNDERS, B. V., (NIST)
C             COMPUTING AND APPLIED MATHEMATICS LABORATORY
C             NATIONAL INSTITUTE OF STANDARDS AND TECHNOLOGY
C             GAITHERSBURG, MD 20899
C           BOISVERT, R. F., (NIST)
C             COMPUTING AND APPLIED MATHEMATICS LABORATORY
C             NATIONAL INSTITUTE OF STANDARDS AND TECHNOLOGY
C             GAITHERSBURG, MD 20899
C***ROUTINES CALLED  IDWCS, D1MACH, DWNGT, DWGTHR, DWGTLE, DWGT, DWLE,
C                    DWSCTR, DWCS
C***REVISION HISTORY  (YYMMDD)
C   910226  DATE WRITTEN 
C***END PROLOGUE  DWI0
C
C  ----------
C  PARAMETERS
C  ----------
C
      INTEGER INFO, INDX, M
      DOUBLE PRECISION B0, B1, B2, F, X, Y, TCMP, YCMP, ZCMP
C
      DIMENSION B0(M), B1(M), B2(M), F(M), INDX(M), X(M), Y(M), 
     +          TCMP(M), YCMP(M), ZCMP(M)
C
C  ---------------
C  LOCAL VARIABLES
C  ---------------
C
      INTEGER LAI0, LAI02, LBI0
      PARAMETER ( LAI0=46, LAI02=69, LBI0=18 )
C
      INTEGER I, IDWCS, J, KEY, N, NTI0, NTAI0, NTAI02
      DOUBLE PRECISION AI0CS, AI02CS, BI0CS, EPMACH, EPS, D1MACH, XSML, 
     +        XMAX
C
      DIMENSION AI0CS(LAI0), AI02CS(LAI02), BI0CS(LBI0)
C
      SAVE AI0CS, AI02CS, BI0CS, N, NTAI0, NTAI02, NTI0, XSML, XMAX
C
C----------------------------------------------------------------------
C
C Series for BI0        on the interval  0.          to  9.00000D+00
C                                        with weighted error   9.51E-34
C                                         log weighted error  33.02
C                               significant figures required  33.31
C                                    decimal places required  33.65
C
      DATA BI0 CS(  1) / -.7660547252 8391449510 8189497624 3285 D-1   /
      DATA BI0 CS(  2) / +.1927337953 9938082699 5240875088 1196 D+1   /
      DATA BI0 CS(  3) / +.2282644586 9203013389 3702929233 0415 D+0   /
      DATA BI0 CS(  4) / +.1304891466 7072904280 7933421069 1888 D-1   /
      DATA BI0 CS(  5) / +.4344270900 8164874513 7868268102 6107 D-3   /
      DATA BI0 CS(  6) / +.9422657686 0019346639 2317174411 8766 D-5   /
      DATA BI0 CS(  7) / +.1434006289 5106910799 6209187817 9957 D-6   /
      DATA BI0 CS(  8) / +.1613849069 6617490699 1541971999 4611 D-8   /
      DATA BI0 CS(  9) / +.1396650044 5356696994 9509270814 2522 D-10  /
      DATA BI0 CS( 10) / +.9579451725 5054453446 2752317189 3333 D-13  /
      DATA BI0 CS( 11) / +.5333981859 8625021310 1510774400 0000 D-15  /
      DATA BI0 CS( 12) / +.2458716088 4374707746 9678591999 9999 D-17  /
      DATA BI0 CS( 13) / +.9535680890 2487700269 4434133333 3333 D-20  /
      DATA BI0 CS( 14) / +.3154382039 7214273367 8933333333 3333 D-22  /
      DATA BI0 CS( 15) / +.9004564101 0946374314 6666666666 6666 D-25  /
      DATA BI0 CS( 16) / +.2240647369 1236700160 0000000000 0000 D-27  /
      DATA BI0 CS( 17) / +.4903034603 2428373333 3333333333 3333 D-30  /
      DATA BI0 CS( 18) / +.9508172606 1226666666 6666666666 6666 D-33  /

		
C
C-------------------------------------------------------------------
C
C Series for AI0        on the interval  1.25000D-01 to  3.33333D-01
C                                        with weighted error   2.74E-32
C                                         log weighted error  31.56
C                               significant figures required  30.15
C                                    decimal places required  32.39
C
      DATA AI0 CS(  1) / +.7575994494 0237959427 2987203743 8 D-1      /
      DATA AI0 CS(  2) / +.7591380810 8233455072 9297873320 4 D-2      /
      DATA AI0 CS(  3) / +.4153131338 9237505018 6319749138 2 D-3      /
      DATA AI0 CS(  4) / +.1070076463 4390730735 8242970217 0 D-4      /
      DATA AI0 CS(  5) / -.7901179979 2128946607 5031948573 0 D-5      /
      DATA AI0 CS(  6) / -.7826143501 4387522697 8898980690 9 D-6      /
      DATA AI0 CS(  7) / +.2783849942 9488708063 8118538985 7 D-6      /
      DATA AI0 CS(  8) / +.8252472600 6120271919 6682913319 8 D-8      /
      DATA AI0 CS(  9) / -.1204463945 5201991790 5496089110 3 D-7      /
      DATA AI0 CS( 10) / +.1559648598 5060764436 1228752792 8 D-8      /
      DATA AI0 CS( 11) / +.2292556367 1033165434 7725480285 7 D-9      /
      DATA AI0 CS( 12) / -.1191622884 2790646036 7777423447 8 D-9      /
      DATA AI0 CS( 13) / +.1757854916 0324098302 1833124774 3 D-10     /
      DATA AI0 CS( 14) / +.1128224463 2189005171 4441135682 4 D-11     /
      DATA AI0 CS( 15) / -.1146848625 9272988777 2963387698 2 D-11     /
      DATA AI0 CS( 16) / +.2715592054 8036628726 4365192160 6 D-12     /
      DATA AI0 CS( 17) / -.2415874666 5626878384 4247572028 1 D-13     /
      DATA AI0 CS( 18) / -.6084469888 2551250646 0609963922 4 D-14     /
      DATA AI0 CS( 19) / +.3145705077 1754772937 0836026730 3 D-14     /
      DATA AI0 CS( 20) / -.7172212924 8711877179 6217505917 6 D-15     /
      DATA AI0 CS( 21) / +.7874493403 4541033960 8390960332 7 D-16     /
      DATA AI0 CS( 22) / +.1004802753 0094624023 4524457183 9 D-16     /
      DATA AI0 CS( 23) / -.7566895365 3505348534 2843588881 0 D-17     /
      DATA AI0 CS( 24) / +.2150380106 8761198878 1205128784 5 D-17     /
      DATA AI0 CS( 25) / -.3754858341 8308744291 5158445260 8 D-18     /
      DATA AI0 CS( 26) / +.2354065842 2269925769 0075710532 2 D-19     /
      DATA AI0 CS( 27) / +.1114667612 0479285302 2637335511 0 D-19     /
      DATA AI0 CS( 28) / -.5398891884 3969903786 9677932270 9 D-20     /
      DATA AI0 CS( 29) / +.1439598792 2407526770 4285840452 2 D-20     /
      DATA AI0 CS( 30) / -.2591916360 1110934064 6081840196 2 D-21     /
      DATA AI0 CS( 31) / +.2238133183 9985839074 3409229824 0 D-22     /
      DATA AI0 CS( 32) / +.5250672575 3647711727 7221683199 9 D-23     /
      DATA AI0 CS( 33) / -.3249904138 5332307841 7343228586 6 D-23     /
      DATA AI0 CS( 34) / +.9924214103 2050379278 5728471040 0 D-24     /
      DATA AI0 CS( 35) / -.2164992254 2446695231 4655429973 3 D-24     /
      DATA AI0 CS( 36) / +.3233609471 9435940839 7333299199 9 D-25     /
      DATA AI0 CS( 37) / -.1184620207 3967424898 2473386666 6 D-26     /
      DATA AI0 CS( 38) / -.1281671853 9504986505 4833868799 9 D-26     /
      DATA AI0 CS( 39) / +.5827015182 2793905116 0556885333 3 D-27     /
      DATA AI0 CS( 40) / -.1668222326 0261097193 6450150399 9 D-27     /
      DATA AI0 CS( 41) / +.3625309510 5415699757 0068480000 0 D-28     /
      DATA AI0 CS( 42) / -.5733627999 0557135899 4595839999 9 D-29     /
      DATA AI0 CS( 43) / +.3736796722 0630982296 4258133333 3 D-30     /
      DATA AI0 CS( 44) / +.1602073983 1568519633 6551253333 3 D-30     /
      DATA AI0 CS( 45) / -.8700424864 0572298845 2249599999 9 D-31     /
      DATA AI0 CS( 46) / +.2741320937 9374811456 0341333333 3 D-31     /
C
C-------------------------------------------------------------------
C
C Series for AI02       on the interval  0.          to  1.25000D-01
C                                        with weighted error   1.97E-32
C                                         log weighted error  31.71
C                               significant figures required  30.15
C                                    decimal places required  32.63
C
      DATA AI02CS(  1) / +.5449041101 4108831607 8960962268 0 D-1      /
      DATA AI02CS(  2) / +.3369116478 2556940898 9785662979 9 D-2      /
      DATA AI02CS(  3) / +.6889758346 9168239842 6263914301 1 D-4      /
      DATA AI02CS(  4) / +.2891370520 8347564829 6692402323 2 D-5      /
      DATA AI02CS(  5) / +.2048918589 4690637418 2760534093 1 D-6      /
      DATA AI02CS(  6) / +.2266668990 4981780645 9327743136 1 D-7      /
      DATA AI02CS(  7) / +.3396232025 7083863451 5084396952 3 D-8      /
      DATA AI02CS(  8) / +.4940602388 2249695891 0482449783 5 D-9      /
      DATA AI02CS(  9) / +.1188914710 7846438342 4084525196 3 D-10     /
      DATA AI02CS( 10) / -.3149916527 9632413645 3864862961 9 D-10     /
      DATA AI02CS( 11) / -.1321581184 0447713118 7540739926 7 D-10     /
      DATA AI02CS( 12) / -.1794178531 5068061177 7943574026 9 D-11     /
      DATA AI02CS( 13) / +.7180124451 3836662336 7106429346 9 D-12     /
      DATA AI02CS( 14) / +.3852778382 7421427011 4089801777 6 D-12     /
      DATA AI02CS( 15) / +.1540086217 5214098269 1325823339 7 D-13     /
      DATA AI02CS( 16) / -.4150569347 2872220866 2689972015 6 D-13     /
      DATA AI02CS( 17) / -.9554846698 8283076487 0214494312 5 D-14     /
      DATA AI02CS( 18) / +.3811680669 3526224207 4605535511 8 D-14     /
      DATA AI02CS( 19) / +.1772560133 0565263836 0493266675 8 D-14     /
      DATA AI02CS( 20) / -.3425485619 6772191346 1924790328 2 D-15     /
      DATA AI02CS( 21) / -.2827623980 5165834849 4205593759 4 D-15     /
      DATA AI02CS( 22) / +.3461222867 6974610930 9706250813 4 D-16     /
      DATA AI02CS( 23) / +.4465621420 2967599990 1042054284 3 D-16     /
      DATA AI02CS( 24) / -.4830504485 9441820712 5525403795 4 D-17     /
      DATA AI02CS( 25) / -.7233180487 8747539545 6227240924 5 D-17     /
      DATA AI02CS( 26) / +.9921475412 1736985988 8046093981 0 D-18     /
      DATA AI02CS( 27) / +.1193650890 8459820855 0439949924 2 D-17     /
      DATA AI02CS( 28) / -.2488709837 1508072357 2054491660 2 D-18     /
      DATA AI02CS( 29) / -.1938426454 1609059289 8469781132 6 D-18     /
      DATA AI02CS( 30) / +.6444656697 3734438687 8301949394 9 D-19     /
      DATA AI02CS( 31) / +.2886051596 2892243264 8171383073 4 D-19     /
      DATA AI02CS( 32) / -.1601954907 1749718070 6167156200 7 D-19     /
      DATA AI02CS( 33) / -.3270815010 5923147208 9193567485 9 D-20     /
      DATA AI02CS( 34) / +.3686932283 8264091811 4600723939 3 D-20     /
      DATA AI02CS( 35) / +.1268297648 0309501530 1359529710 9 D-22     /
      DATA AI02CS( 36) / -.7549825019 3772739076 9636664410 1 D-21     /
      DATA AI02CS( 37) / +.1502133571 3778353496 3712789053 4 D-21     /
      DATA AI02CS( 38) / +.1265195883 5096485349 3208799248 3 D-21     /
      DATA AI02CS( 39) / -.6100998370 0836807086 2940891600 2 D-22     /
      DATA AI02CS( 40) / -.1268809629 2601282643 6872095924 2 D-22     /
      DATA AI02CS( 41) / +.1661016099 8907414578 4038487490 5 D-22     /
      DATA AI02CS( 42) / -.1585194335 7658855793 7970504881 4 D-23     /
      DATA AI02CS( 43) / -.3302645405 9682178009 5381766755 6 D-23     /
      DATA AI02CS( 44) / +.1313580902 8392397817 4039623117 4 D-23     /
      DATA AI02CS( 45) / +.3689040246 6711567933 1425637280 4 D-24     /
      DATA AI02CS( 46) / -.4210141910 4616891492 1978247249 9 D-24     /
      DATA AI02CS( 47) / +.4791954591 0828657806 3171401373 0 D-25     /
      DATA AI02CS( 48) / +.8459470390 2218217952 9971707412 4 D-25     /
      DATA AI02CS( 49) / -.4039800940 8728324931 4607937181 0 D-25     /
      DATA AI02CS( 50) / -.6434714653 6504313473 0100850469 5 D-26     /
      DATA AI02CS( 51) / +.1225743398 8756659903 4464736990 5 D-25     /
      DATA AI02CS( 52) / -.2934391316 0257089231 9879821175 4 D-26     /
      DATA AI02CS( 53) / -.1961311309 1949829262 0371205728 9 D-26     /
      DATA AI02CS( 54) / +.1503520374 8221934241 6229900309 8 D-26     /
      DATA AI02CS( 55) / -.9588720515 7448265520 3386388206 9 D-28     /
      DATA AI02CS( 56) / -.3483339380 8170454863 9441108511 4 D-27     /
      DATA AI02CS( 57) / +.1690903610 2630436730 6244960725 6 D-27     /
      DATA AI02CS( 58) / +.1982866538 7356030438 9400115718 8 D-28     /
      DATA AI02CS( 59) / -.5317498081 4918162145 7583002528 4 D-28     /
      DATA AI02CS( 60) / +.1803306629 8883929462 3501450390 1 D-28     /
      DATA AI02CS( 61) / +.6213093341 4548931758 8405311242 2 D-29     /
      DATA AI02CS( 62) / -.7692189292 7721618632 0072806673 0 D-29     /
      DATA AI02CS( 63) / +.1858252826 1117025426 2556016596 3 D-29     /
      DATA AI02CS( 64) / +.1237585142 2813957248 9927154554 1 D-29     /
      DATA AI02CS( 65) / -.1102259120 4092238032 1779478779 2 D-29     /
      DATA AI02CS( 66) / +.1886287118 0397044900 7787447943 1 D-30     /
      DATA AI02CS( 67) / +.2160196872 2436589131 4903141406 0 D-30     /
      DATA AI02CS( 68) / -.1605454124 9197432005 8446594965 5 D-30     /
      DATA AI02CS( 69) / +.1965352984 5942906039 3884807331 8 D-31     /
C
C-------------------------------------------------------------------
C
      DATA NTI0 / 0 /
C
C***FIRST EXECUTABLE STATEMENT  DWI0
C
      IF (M .LE. 0)  GO TO 910
C
      IF (NTI0.EQ.0) THEN
         EPMACH = D1MACH(3)
         EPS = 0.10D0*EPMACH
         NTI0 = IDWCS(BI0CS, LBI0, EPS)
         NTAI0 = IDWCS(AI0CS, LAI0, EPS)
         NTAI02 = IDWCS(AI02CS, LAI02, EPS)
         XSML = SQRT(4.0D0*EPMACH)
         XMAX = LOG(D1MACH(2))
      ENDIF
C
      DO 10 I=1,M
         Y(I) = ABS(X(I))
  10  CONTINUE
C
      CALL DWNGT(M,Y,XMAX,KEY)
      IF (KEY .NE. 0)  GO TO 920
C
C  ----------------
C  CASE Y .LE. XSML
C  ----------------
C
      DO 15 I=1,M
         F(I) = 1.0D0
  15  CONTINUE 
C
C  --------------------------
C  CASE  XSML .LT. Y .LE. 3.0
C  --------------------------
C
      CALL DWGTLE(M,Y,XSML,3.0D0,N,INDX)
      IF (N .GT. 0) THEN
         CALL DWGTHR(N,Y,INDX,YCMP)
         DO 20 J=1,N
            TCMP(J) = YCMP(J)**2/4.50D0 - 1.0D0
  20     CONTINUE
         CALL DWCS(N,TCMP,BI0CS,NTI0,ZCMP,B0,B1,B2)
         DO 30 J=1,N
            ZCMP(J) = 2.750D0 + ZCMP(J)
  30     CONTINUE
         CALL DWSCTR(N,ZCMP,INDX,F)
      ENDIF
C
C  -------------------------
C  CASE  3.0 .LT. Y .LE. 8.0
C  -------------------------
C
      CALL DWGTLE(M,Y,3.0D0,8.0D0,N,INDX)
      IF (N .GT. 0) THEN
         CALL DWGTHR(N,Y,INDX,YCMP)
         DO 50 J=1,N
            TCMP(J) = (48.0D0/YCMP(J) - 11.0D0)/5.0D0
  50     CONTINUE
         CALL DWCS(N,TCMP,AI0CS,NTAI0,ZCMP,B0,B1,B2)
         DO 60 J=1,N
            ZCMP(J) = EXP(YCMP(J))*(0.3750D0+ZCMP(J))/SQRT(YCMP(J))
  60     CONTINUE
         CALL DWSCTR(N,ZCMP,INDX,F)
      ENDIF
C
C  ----------------
C  CASE  Y .GT. 8.0
C  ----------------
C
      CALL DWGT(M,Y,8.0D0,N,INDX)
      IF (N .GT. 0) THEN
         CALL DWGTHR(N,Y,INDX,YCMP)
        DO 80 J=1,N
            TCMP(J) = 16.0D0/YCMP(J) - 1.0D0
  80     CONTINUE
         CALL DWCS(N,TCMP,AI02CS,NTAI02,ZCMP,B0,B1,B2)
         DO 90 J=1,N
            ZCMP(J) = EXP(YCMP(J))*(0.3750D0+ZCMP(J))/SQRT(YCMP(J))
  90     CONTINUE
         CALL DWSCTR(N,ZCMP,INDX,F)
      ENDIF
C
C  -----
C  EXITS
C  -----
C
C     ... NORMAL
C
      INFO = 0
      GO TO 999
C
C     ... M .LE. 0
C
  910 CONTINUE
      INFO = 1
      GO TO 999
C
C     ... ABS(X) SO LARGE I0 OVERFLOWS
C
  920 CONTINUE
      INFO = 2
      INDX(1) = KEY
      GO TO 999
C
  999 CONTINUE
      RETURN
      END

		
.


AD:

NEW PAGES:

[ODDNUGGET]

[GOPHER]