src/numerical/slatec/fortran/dyairy.f

   1 *DECK DYAIRY
   2       SUBROUTINE DYAIRY (X, RX, C, BI, DBI)
   3 C***BEGIN PROLOGUE  DYAIRY
   4 C***SUBSIDIARY
   5 C***PURPOSE  Subsidiary to DBESJ and DBESY
   6 C***LIBRARY   SLATEC
   7 C***TYPE      DOUBLE PRECISION (YAIRY-S, DYAIRY-D)
   8 C***AUTHOR  Amos, D. E., (SNLA)
   9 C           Daniel, S. L., (SNLA)
  10 C***DESCRIPTION
  11 C
  12 C                  DYAIRY computes the Airy function BI(X)
  13 C                   and its derivative DBI(X) for DASYJY
  14 C
  15 C                                     INPUT
  16 C
  17 C         X  - Argument, computed by DASYJY, X unrestricted
  18 C        RX  - RX=SQRT(ABS(X)), computed by DASYJY
  19 C         C  - C=2.*(ABS(X)**1.5)/3., computed by DASYJY
  20 C
  21 C                                    OUTPUT
  22 C        BI  - Value of function BI(X)
  23 C       DBI  - Value of the derivative DBI(X)
  24 C
  25 C***SEE ALSO  DBESJ, DBESY
  26 C***ROUTINES CALLED  (NONE)
  27 C***REVISION HISTORY  (YYMMDD)
  28 C   750101  DATE WRITTEN
  29 C   890531  Changed all specific intrinsics to generic.  (WRB)
  30 C   891214  Prologue converted to Version 4.0 format.  (BAB)
  31 C   900328  Added TYPE section.  (WRB)
  32 C   910408  Updated the AUTHOR section.  (WRB)
  33 C***END PROLOGUE  DYAIRY
  34 C
  35       INTEGER I, J, M1, M1D, M2, M2D, M3, M3D, M4D, N1, N1D, N2, N2D,
  36      1 N3, N3D, N4D
  37       DOUBLE PRECISION AA,AX,BB,BI,BJN,BJP,BK1,BK2,BK3,BK4,C,CON1,CON2,
  38      1 CON3, CV, DAA, DBB, DBI, DBJN, DBJP, DBK1, DBK2, DBK3, DBK4, D1,
  39      2 D2, EX, E1, E2, FPI12, F1, F2, RTRX, RX, SPI12, S1, S2, T, TC,
  40      3 TEMP1, TEMP2, TT, X
  41       DIMENSION BK1(20), BK2(20), BK3(20), BK4(14)
  42       DIMENSION BJP(19), BJN(19), AA(14), BB(14)
  43       DIMENSION DBK1(21), DBK2(20), DBK3(20), DBK4(14)
  44       DIMENSION DBJP(19), DBJN(19), DAA(14), DBB(14)
  45       SAVE N1, N2, N3, M1, M2, M3, N1D, N2D, N3D, N4D,
  46      1 M1D, M2D, M3D, M4D, FPI12, SPI12, CON1, CON2, CON3,
  47      2 BK1, BK2, BK3, BK4, BJN, BJP, AA, BB, DBK1, DBK2, DBK3, DBK4,
  48      3 DBJP, DBJN, DAA, DBB
  49       DATA N1,N2,N3/20,19,14/
  50       DATA M1,M2,M3/18,17,12/
  51       DATA N1D,N2D,N3D,N4D/21,20,19,14/
  52       DATA M1D,M2D,M3D,M4D/19,18,17,12/
  53       DATA FPI12,SPI12,CON1,CON2,CON3/
  54      1 1.30899693899575D+00, 1.83259571459405D+00, 6.66666666666667D-01,
  55      2 7.74148278841779D+00, 3.64766105490356D-01/
  56       DATA BK1(1),  BK1(2),  BK1(3),  BK1(4),  BK1(5),  BK1(6),
  57      1     BK1(7),  BK1(8),  BK1(9),  BK1(10), BK1(11), BK1(12),
  58      2     BK1(13), BK1(14), BK1(15), BK1(16), BK1(17), BK1(18),
  59      3     BK1(19), BK1(20)/ 2.43202846447449D+00, 2.57132009754685D+00,
  60      4 1.02802341258616D+00, 3.41958178205872D-01, 8.41978629889284D-02,
  61      5 1.93877282587962D-02, 3.92687837130335D-03, 6.83302689948043D-04,
  62      6 1.14611403991141D-04, 1.74195138337086D-05, 2.41223620956355D-06,
  63      7 3.24525591983273D-07, 4.03509798540183D-08, 4.70875059642296D-09,
  64      8 5.35367432585889D-10, 5.70606721846334D-11, 5.80526363709933D-12,
  65      9 5.76338988616388D-13, 5.42103834518071D-14, 4.91857330301677D-15/
  66       DATA BK2(1),  BK2(2),  BK2(3),  BK2(4),  BK2(5),  BK2(6),
  67      1     BK2(7),  BK2(8),  BK2(9),  BK2(10), BK2(11), BK2(12),
  68      2     BK2(13), BK2(14), BK2(15), BK2(16), BK2(17), BK2(18),
  69      3     BK2(19), BK2(20)/ 5.74830555784088D-01,-6.91648648376891D-03,
  70      4 1.97460263052093D-03,-5.24043043868823D-04, 1.22965147239661D-04,
  71      5-2.27059514462173D-05, 2.23575555008526D-06, 4.15174955023899D-07,
  72      6-2.84985752198231D-07, 8.50187174775435D-08,-1.70400826891326D-08,
  73      7 2.25479746746889D-09,-1.09524166577443D-10,-3.41063845099711D-11,
  74      8 1.11262893886662D-11,-1.75542944241734D-12, 1.36298600401767D-13,
  75      9 8.76342105755664D-15,-4.64063099157041D-15, 7.78772758732960D-16/
  76       DATA BK3(1),  BK3(2),  BK3(3),  BK3(4),  BK3(5),  BK3(6),
  77      1     BK3(7),  BK3(8),  BK3(9),  BK3(10), BK3(11), BK3(12),
  78      2     BK3(13), BK3(14), BK3(15), BK3(16), BK3(17), BK3(18),
  79      3     BK3(19), BK3(20)/ 5.66777053506912D-01, 2.63672828349579D-03,
  80      4 5.12303351473130D-05, 2.10229231564492D-06, 1.42217095113890D-07,
  81      5 1.28534295891264D-08, 7.28556219407507D-10,-3.45236157301011D-10,
  82      6-2.11919115912724D-10,-6.56803892922376D-11,-8.14873160315074D-12,
  83      7 3.03177845632183D-12, 1.73447220554115D-12, 1.67935548701554D-13,
  84      8-1.49622868806719D-13,-5.15470458953407D-14, 8.75741841857830D-15,
  85      9 7.96735553525720D-15,-1.29566137861742D-16,-1.11878794417520D-15/
  86       DATA BK4(1),  BK4(2),  BK4(3),  BK4(4),  BK4(5),  BK4(6),
  87      1     BK4(7),  BK4(8),  BK4(9),  BK4(10), BK4(11), BK4(12),
  88      2     BK4(13), BK4(14)/ 4.85444386705114D-01,-3.08525088408463D-03,
  89      3 6.98748404837928D-05,-2.82757234179768D-06, 1.59553313064138D-07,
  90      4-1.12980692144601D-08, 9.47671515498754D-10,-9.08301736026423D-11,
  91      5 9.70776206450724D-12,-1.13687527254574D-12, 1.43982917533415D-13,
  92      6-1.95211019558815D-14, 2.81056379909357D-15,-4.26916444775176D-16/
  93       DATA BJP(1),  BJP(2),  BJP(3),  BJP(4),  BJP(5),  BJP(6),
  94      1     BJP(7),  BJP(8),  BJP(9),  BJP(10), BJP(11), BJP(12),
  95      2     BJP(13), BJP(14), BJP(15), BJP(16), BJP(17), BJP(18),
  96      3     BJP(19)         / 1.34918611457638D-01,-3.19314588205813D-01,
  97      4 5.22061946276114D-02, 5.28869112170312D-02,-8.58100756077350D-03,
  98      5-2.99211002025555D-03, 4.21126741969759D-04, 8.73931830369273D-05,
  99      6-1.06749163477533D-05,-1.56575097259349D-06, 1.68051151983999D-07,
 100      7 1.89901103638691D-08,-1.81374004961922D-09,-1.66339134593739D-10,
 101      8 1.42956335780810D-11, 1.10179811626595D-12,-8.60187724192263D-14,
 102      9-5.71248177285064D-15, 4.08414552853803D-16/
 103       DATA BJN(1),  BJN(2),  BJN(3),  BJN(4),  BJN(5),  BJN(6),
 104      1     BJN(7),  BJN(8),  BJN(9),  BJN(10), BJN(11), BJN(12),
 105      2     BJN(13), BJN(14), BJN(15), BJN(16), BJN(17), BJN(18),
 106      3     BJN(19)         / 6.59041673525697D-02,-4.24905910566004D-01,
 107      4 2.87209745195830D-01, 1.29787771099606D-01,-4.56354317590358D-02,
 108      5-1.02630175982540D-02, 2.50704671521101D-03, 3.78127183743483D-04,
 109      6-7.11287583284084D-05,-8.08651210688923D-06, 1.23879531273285D-06,
 110      7 1.13096815867279D-07,-1.46234283176310D-08,-1.11576315688077D-09,
 111      8 1.24846618243897D-10, 8.18334132555274D-12,-8.07174877048484D-13,
 112      9-4.63778618766425D-14, 4.09043399081631D-15/
 113       DATA AA(1),   AA(2),   AA(3),   AA(4),   AA(5),   AA(6),
 114      1     AA(7),   AA(8),   AA(9),   AA(10),  AA(11),  AA(12),
 115      2     AA(13),  AA(14) /-2.78593552803079D-01, 3.52915691882584D-03,
 116      3 2.31149677384994D-05,-4.71317842263560D-06, 1.12415907931333D-07,
 117      4 2.00100301184339D-08,-2.60948075302193D-09, 3.55098136101216D-11,
 118      5 3.50849978423875D-11,-5.83007187954202D-12, 2.04644828753326D-13,
 119      6 1.10529179476742D-13,-2.87724778038775D-14, 2.88205111009939D-15/
 120       DATA BB(1),   BB(2),   BB(3),   BB(4),   BB(5),   BB(6),
 121      1     BB(7),   BB(8),   BB(9),   BB(10),  BB(11),  BB(12),
 122      2     BB(13),  BB(14) /-4.90275424742791D-01,-1.57647277946204D-03,
 123      3 9.66195963140306D-05,-1.35916080268815D-07,-2.98157342654859D-07,
 124      4 1.86824767559979D-08, 1.03685737667141D-09,-3.28660818434328D-10,
 125      5 2.57091410632780D-11, 2.32357655300677D-12,-9.57523279048255D-13,
 126      6 1.20340828049719D-13, 2.90907716770715D-15,-4.55656454580149D-15/
 127       DATA DBK1(1), DBK1(2), DBK1(3), DBK1(4), DBK1(5), DBK1(6),
 128      1     DBK1(7), DBK1(8), DBK1(9), DBK1(10),DBK1(11),DBK1(12),
 129      2     DBK1(13),DBK1(14),DBK1(15),DBK1(16),DBK1(17),DBK1(18),
 130      3     DBK1(19),DBK1(20),
 131      4     DBK1(21)        / 2.95926143981893D+00, 3.86774568440103D+00,
 132      5 1.80441072356289D+00, 5.78070764125328D-01, 1.63011468174708D-01,
 133      6 3.92044409961855D-02, 7.90964210433812D-03, 1.50640863167338D-03,
 134      7 2.56651976920042D-04, 3.93826605867715D-05, 5.81097771463818D-06,
 135      8 7.86881233754659D-07, 9.93272957325739D-08, 1.21424205575107D-08,
 136      9 1.38528332697707D-09, 1.50190067586758D-10, 1.58271945457594D-11,
 137      1 1.57531847699042D-12, 1.50774055398181D-13, 1.40594335806564D-14,
 138      2 1.24942698777218D-15/
 139       DATA DBK2(1), DBK2(2), DBK2(3), DBK2(4), DBK2(5), DBK2(6),
 140      1     DBK2(7), DBK2(8), DBK2(9), DBK2(10),DBK2(11),DBK2(12),
 141      2     DBK2(13),DBK2(14),DBK2(15),DBK2(16),DBK2(17),DBK2(18),
 142      3    DBK2(19),DBK2(20)/ 5.49756809432471D-01, 9.13556983276901D-03,
 143      4-2.53635048605507D-03, 6.60423795342054D-04,-1.55217243135416D-04,
 144      5 3.00090325448633D-05,-3.76454339467348D-06,-1.33291331611616D-07,
 145      6 2.42587371049013D-07,-8.07861075240228D-08, 1.71092818861193D-08,
 146      7-2.41087357570599D-09, 1.53910848162371D-10, 2.56465373190630D-11,
 147      8-9.88581911653212D-12, 1.60877986412631D-12,-1.20952524741739D-13,
 148      9-1.06978278410820D-14, 5.02478557067561D-15,-8.68986130935886D-16/
 149       DATA DBK3(1), DBK3(2), DBK3(3), DBK3(4), DBK3(5), DBK3(6),
 150      1     DBK3(7), DBK3(8), DBK3(9), DBK3(10),DBK3(11),DBK3(12),
 151      2     DBK3(13),DBK3(14),DBK3(15),DBK3(16),DBK3(17),DBK3(18),
 152      3    DBK3(19),DBK3(20)/ 5.60598509354302D-01,-3.64870013248135D-03,
 153      4-5.98147152307417D-05,-2.33611595253625D-06,-1.64571516521436D-07,
 154      5-2.06333012920569D-08,-4.27745431573110D-09,-1.08494137799276D-09,
 155      6-2.37207188872763D-10,-2.22132920864966D-11, 1.07238008032138D-11,
 156      7 5.71954845245808D-12, 7.51102737777835D-13,-3.81912369483793D-13,
 157      8-1.75870057119257D-13, 6.69641694419084D-15, 2.26866724792055D-14,
 158      9 2.69898141356743D-15,-2.67133612397359D-15,-6.54121403165269D-16/
 159       DATA DBK4(1), DBK4(2), DBK4(3), DBK4(4), DBK4(5), DBK4(6),
 160      1     DBK4(7), DBK4(8), DBK4(9), DBK4(10),DBK4(11),DBK4(12),
 161      2    DBK4(13),DBK4(14)/ 4.93072999188036D-01, 4.38335419803815D-03,
 162      3-8.37413882246205D-05, 3.20268810484632D-06,-1.75661979548270D-07,
 163      4 1.22269906524508D-08,-1.01381314366052D-09, 9.63639784237475D-11,
 164      5-1.02344993379648D-11, 1.19264576554355D-12,-1.50443899103287D-13,
 165      6 2.03299052379349D-14,-2.91890652008292D-15, 4.42322081975475D-16/
 166       DATA DBJP(1), DBJP(2), DBJP(3), DBJP(4), DBJP(5), DBJP(6),
 167      1     DBJP(7), DBJP(8), DBJP(9), DBJP(10),DBJP(11),DBJP(12),
 168      2     DBJP(13),DBJP(14),DBJP(15),DBJP(16),DBJP(17),DBJP(18),
 169      3     DBJP(19)        / 1.13140872390745D-01,-2.08301511416328D-01,
 170      4 1.69396341953138D-02, 2.90895212478621D-02,-3.41467131311549D-03,
 171      5-1.46455339197417D-03, 1.63313272898517D-04, 3.91145328922162D-05,
 172      6-3.96757190808119D-06,-6.51846913772395D-07, 5.98707495269280D-08,
 173      7 7.44108654536549D-09,-6.21241056522632D-10,-6.18768017313526D-11,
 174      8 4.72323484752324D-12, 3.91652459802532D-13,-2.74985937845226D-14,
 175      9-1.95036497762750D-15, 1.26669643809444D-16/
 176       DATA DBJN(1), DBJN(2), DBJN(3), DBJN(4), DBJN(5), DBJN(6),
 177      1     DBJN(7), DBJN(8), DBJN(9), DBJN(10),DBJN(11),DBJN(12),
 178      2     DBJN(13),DBJN(14),DBJN(15),DBJN(16),DBJN(17),DBJN(18),
 179      3     DBJN(19)        /-1.88091260068850D-02,-1.47798180826140D-01,
 180      4 5.46075900433171D-01, 1.52146932663116D-01,-9.58260412266886D-02,
 181      5-1.63102731696130D-02, 5.75364806680105D-03, 7.12145408252655D-04,
 182      6-1.75452116846724D-04,-1.71063171685128D-05, 3.24435580631680D-06,
 183      7 2.61190663932884D-07,-4.03026865912779D-08,-2.76435165853895D-09,
 184      8 3.59687929062312D-10, 2.14953308456051D-11,-2.41849311903901D-12,
 185      9-1.28068004920751D-13, 1.26939834401773D-14/
 186       DATA DAA(1),  DAA(2),  DAA(3),  DAA(4),  DAA(5),  DAA(6),
 187      1     DAA(7),  DAA(8),  DAA(9),  DAA(10), DAA(11), DAA(12),
 188      2     DAA(13), DAA(14)/ 2.77571356944231D-01,-4.44212833419920D-03,
 189      3 8.42328522190089D-05, 2.58040318418710D-06,-3.42389720217621D-07,
 190      4 6.24286894709776D-09, 2.36377836844577D-09,-3.16991042656673D-10,
 191      5 4.40995691658191D-12, 5.18674221093575D-12,-9.64874015137022D-13,
 192      6 4.90190576608710D-14, 1.77253430678112D-14,-5.55950610442662D-15/
 193       DATA DBB(1),  DBB(2),  DBB(3),  DBB(4),  DBB(5),  DBB(6),
 194      1     DBB(7),  DBB(8),  DBB(9),  DBB(10), DBB(11), DBB(12),
 195      2     DBB(13), DBB(14)/ 4.91627321104601D-01, 3.11164930427489D-03,
 196      3 8.23140762854081D-05,-4.61769776172142D-06,-6.13158880534626D-08,
 197      4 2.87295804656520D-08,-1.81959715372117D-09,-1.44752826642035D-10,
 198      5 4.53724043420422D-11,-3.99655065847223D-12,-3.24089119830323D-13,
 199      6 1.62098952568741D-13,-2.40765247974057D-14, 1.69384811284491D-16/
 200 C***FIRST EXECUTABLE STATEMENT  DYAIRY
 201       AX = ABS(X)
 202       RX = SQRT(AX)
 203       C = CON1*AX*RX
 204       IF (X.LT.0.0D0) GO TO 120
 205       IF (C.GT.8.0D0) GO TO 60
 206       IF (X.GT.2.5D0) GO TO 30
 207       T = (X+X-2.5D0)*0.4D0
 208       TT = T + T
 209       J = N1
 210       F1 = BK1(J)
 211       F2 = 0.0D0
 212       DO 10 I=1,M1
 213         J = J - 1
 214         TEMP1 = F1
 215         F1 = TT*F1 - F2 + BK1(J)
 216         F2 = TEMP1
 217    10 CONTINUE
 218       BI = T*F1 - F2 + BK1(1)
 219       J = N1D
 220       F1 = DBK1(J)
 221       F2 = 0.0D0
 222       DO 20 I=1,M1D
 223         J = J - 1
 224         TEMP1 = F1
 225         F1 = TT*F1 - F2 + DBK1(J)
 226         F2 = TEMP1
 227    20 CONTINUE
 228       DBI = T*F1 - F2 + DBK1(1)
 229       RETURN
 230    30 CONTINUE
 231       RTRX = SQRT(RX)
 232       T = (X+X-CON2)*CON3
 233       TT = T + T
 234       J = N1
 235       F1 = BK2(J)
 236       F2 = 0.0D0
 237       DO 40 I=1,M1
 238         J = J - 1
 239         TEMP1 = F1
 240         F1 = TT*F1 - F2 + BK2(J)
 241         F2 = TEMP1
 242    40 CONTINUE
 243       BI = (T*F1-F2+BK2(1))/RTRX
 244       EX = EXP(C)
 245       BI = BI*EX
 246       J = N2D
 247       F1 = DBK2(J)
 248       F2 = 0.0D0
 249       DO 50 I=1,M2D
 250         J = J - 1
 251         TEMP1 = F1
 252         F1 = TT*F1 - F2 + DBK2(J)
 253         F2 = TEMP1
 254    50 CONTINUE
 255       DBI = (T*F1-F2+DBK2(1))*RTRX
 256       DBI = DBI*EX
 257       RETURN
 258 C
 259    60 CONTINUE
 260       RTRX = SQRT(RX)
 261       T = 16.0D0/C - 1.0D0
 262       TT = T + T
 263       J = N1
 264       F1 = BK3(J)
 265       F2 = 0.0D0
 266       DO 70 I=1,M1
 267         J = J - 1
 268         TEMP1 = F1
 269         F1 = TT*F1 - F2 + BK3(J)
 270         F2 = TEMP1
 271    70 CONTINUE
 272       S1 = T*F1 - F2 + BK3(1)
 273       J = N2D
 274       F1 = DBK3(J)
 275       F2 = 0.0D0
 276       DO 80 I=1,M2D
 277         J = J - 1
 278         TEMP1 = F1
 279         F1 = TT*F1 - F2 + DBK3(J)
 280         F2 = TEMP1
 281    80 CONTINUE
 282       D1 = T*F1 - F2 + DBK3(1)
 283       TC = C + C
 284       EX = EXP(C)
 285       IF (TC.GT.35.0D0) GO TO 110
 286       T = 10.0D0/C - 1.0D0
 287       TT = T + T
 288       J = N3
 289       F1 = BK4(J)
 290       F2 = 0.0D0
 291       DO 90 I=1,M3
 292         J = J - 1
 293         TEMP1 = F1
 294         F1 = TT*F1 - F2 + BK4(J)
 295         F2 = TEMP1
 296    90 CONTINUE
 297       S2 = T*F1 - F2 + BK4(1)
 298       BI = (S1+EXP(-TC)*S2)/RTRX
 299       BI = BI*EX
 300       J = N4D
 301       F1 = DBK4(J)
 302       F2 = 0.0D0
 303       DO 100 I=1,M4D
 304         J = J - 1
 305         TEMP1 = F1
 306         F1 = TT*F1 - F2 + DBK4(J)
 307         F2 = TEMP1
 308   100 CONTINUE
 309       D2 = T*F1 - F2 + DBK4(1)
 310       DBI = RTRX*(D1+EXP(-TC)*D2)
 311       DBI = DBI*EX
 312       RETURN
 313   110 BI = EX*S1/RTRX
 314       DBI = EX*RTRX*D1
 315       RETURN
 316 C
 317   120 CONTINUE
 318       IF (C.GT.5.0D0) GO TO 150
 319       T = 0.4D0*C - 1.0D0
 320       TT = T + T
 321       J = N2
 322       F1 = BJP(J)
 323       E1 = BJN(J)
 324       F2 = 0.0D0
 325       E2 = 0.0D0
 326       DO 130 I=1,M2
 327         J = J - 1
 328         TEMP1 = F1
 329         TEMP2 = E1
 330         F1 = TT*F1 - F2 + BJP(J)
 331         E1 = TT*E1 - E2 + BJN(J)
 332         F2 = TEMP1
 333         E2 = TEMP2
 334   130 CONTINUE
 335       BI = (T*E1-E2+BJN(1)) - AX*(T*F1-F2+BJP(1))
 336       J = N3D
 337       F1 = DBJP(J)
 338       E1 = DBJN(J)
 339       F2 = 0.0D0
 340       E2 = 0.0D0
 341       DO 140 I=1,M3D
 342         J = J - 1
 343         TEMP1 = F1
 344         TEMP2 = E1
 345         F1 = TT*F1 - F2 + DBJP(J)
 346         E1 = TT*E1 - E2 + DBJN(J)
 347         F2 = TEMP1
 348         E2 = TEMP2
 349   140 CONTINUE
 350       DBI = X*X*(T*F1-F2+DBJP(1)) + (T*E1-E2+DBJN(1))
 351       RETURN
 352 C
 353   150 CONTINUE
 354       RTRX = SQRT(RX)
 355       T = 10.0D0/C - 1.0D0
 356       TT = T + T
 357       J = N3
 358       F1 = AA(J)
 359       E1 = BB(J)
 360       F2 = 0.0D0
 361       E2 = 0.0D0
 362       DO 160 I=1,M3
 363         J = J - 1
 364         TEMP1 = F1
 365         TEMP2 = E1
 366         F1 = TT*F1 - F2 + AA(J)
 367         E1 = TT*E1 - E2 + BB(J)
 368         F2 = TEMP1
 369         E2 = TEMP2
 370   160 CONTINUE
 371       TEMP1 = T*F1 - F2 + AA(1)
 372       TEMP2 = T*E1 - E2 + BB(1)
 373       CV = C - FPI12
 374       BI = (TEMP1*COS(CV)+TEMP2*SIN(CV))/RTRX
 375       J = N4D
 376       F1 = DAA(J)
 377       E1 = DBB(J)
 378       F2 = 0.0D0
 379       E2 = 0.0D0
 380       DO 170 I=1,M4D
 381         J = J - 1
 382         TEMP1 = F1
 383         TEMP2 = E1
 384         F1 = TT*F1 - F2 + DAA(J)
 385         E1 = TT*E1 - E2 + DBB(J)
 386         F2 = TEMP1
 387         E2 = TEMP2
 388   170 CONTINUE
 389       TEMP1 = T*F1 - F2 + DAA(1)
 390       TEMP2 = T*E1 - E2 + DBB(1)
 391       CV = C - SPI12
 392       DBI = (TEMP1*COS(CV)-TEMP2*SIN(CV))*RTRX
 393       RETURN
 394       END