Trailing-Edge
-
PDP-10 Archives
-
decuslib20-01
-
decus/20-0025/abeta.for
There is 1 other file named abeta.for in the archive. Click here to see a list.
00010 C ABETA.SRC HFW JULY 15, 1969 75 LINES
00020 C EVALUATE THE INCOMPLETE BETA FUNCTION (NEEDS GAMXX*)
00030 FUNCTION ABETA(XX,A,B)
00040 F=0.
00050 X=XX
00060 IF(X)70,70,2
00070 2 C=1.-X
00080 IF(C)60,60,4
00090 4 IF(A)70,70,6
00100 6 IF(B)70,70,8
00110 8 R=A+B-2.
00120 IF(R)9,10,9
00130 9 IF(X-(A-1.)/R)11,11,12
00140 10 IF(X-.5)11,11,12
00150 11 LOOP=1
00160 TA=A
00170 TB=B
00180 GO TO 14
00190 12 LOOP=2
00200 TA=B
00210 TB=A
00220 X=C
00230 C=1.-X
00240 14 AB=A+B
00250 BL=TB-1.
00260 A1=TA-1.
00270 A2=TA-2.
00280 HA=TA
00290 HB=TB
00300 F=(X**TA)*(C**BL)/TA
00310 18 IF(TA-2.)22,22,20
00320 20 TA=TA-1.
00330 F=F/TA
00340 22 IF(TB-2.)26,26,24
00350 24 TB=TB-1.
00360 F=F/TB
00370 26 IF(AB-2.)30,30,28
00380 28 AB=AB-1.
00390 F=F*AB
00400 GO TO 18
00410 30 GA=GAM(TA)
00420 GB=GAM(TB)
00430 GAB=GAM(AB)
00440 F=F*GAB/(GA*GB)
00450 R=X/C
00460 ROLD=0.
00470 AL=0.
00480 AC=1.
00490 BL=1.
00500 BC=1.
00510 H=0.
00520 DO 45 K=1,20
00530 DO 40 M=1,5
00540 H=H+1.
00550 HH=H+H
00560 E=(A1+H)*(H-HB)*R/((A2+HH)*(A1+HH))
00570 EP=H*(A1+HB+H)*R/((A1+HH)*(HA+HH))
00580 AL=AL*E+AC
00590 BL=BL*E+BC
00600 AC=AC*EP+AL
00610 BC=BC*EP+BL
00620 RNEW=AC/BC
00630 IF(ABS((ROLD/RNEW)-1.)-.0000007)50,50,40
00640 40 ROLD=RNEW
00650 BIG=AMAX1(ABS(AL),ABS(BL),ABS(AC),ABS(BC))
00660 AL=AL/BIG
00670 BL=BL/BIG
00680 AC=AC/BIG
00690 45 BC=BC/BIG
00700 50 F=F*RNEW
00710 GO TO (70,60),LOOP
00720 60 F=1.-F
00730 70 ABETA=F
00740 RETURN
00750 END