C
B = UPPER LIMIT
C
C = ARRAY OF FUNCTION DEFINITION PARAMETERS IF REQUIRED
C
ERR = THE DESIRED ACCURACY
C
RES = THE RESULTING VALUE FOR THE INTEGRAL
EXTERNAL FUN
C Z IS THE ARRAY OF APPROXIMATIONS
DIMENSION Z(10,10)
C INITIALIZE AND COMPUTE THE FIRST APPROXIMATION.
I=1
DEL=B-A
Z(1,1)=.5*DEL*(FUN(A,C)+FUN(B,C))
C THE MAIN LOOP. THE FIRST PART COMPUTES THE INTEGRAL USING A
C 2J+1 POINT TRAPEZOID RULE. THE METHODS MAKES MAXIMAL USE
OF THE
C VALUES ALREADY COMPUTED.
10 J=2**(I-1)
DEL=DEL/2.
I=I+1
Z(I,1)=.5*Z(I-l,l)
DO 1 K=1,J
X=A+(2.*K-1)*DEL
Z(I,l)=Z(I,l)+DEL*FUN(X,C)
1 CONTINUE
C NOW WE DO THE RICHARDSION EXTRAPOLATION
DO 2 K=2,I
Z(I,K)=(4.**(K-1)*Z(I,K-l)-Z(I-l,K-1))/(4.**(K-1)-1.)
2 CONTINUE
C ERROR CONTROL
DIFF=ABS(Z(I,I)-Z(I,I-1))
PRINT 108, DIFF,Z(I,I)
108 FORMAT(5X,2E15.5)
IF(DIFF.LT.ERR) GO TO 20
C THE MAXIMUM NUMBER OF ITERATIONS ALLOWED IS 10.
IF(I.LT.10) GO TO 10
PRINT 100
100 FORMAT(` MORE THAN 10 ITERATIONS REQUIRED, CHECK
PARAMETERS.')
STOP
20 RES=Z(I,I)
RETURN
END
Function FUN(T,C)
FUNCTION FUN(T,C)
C THIS FUNCTION CONTAINS THE INTEGRAND OF THE I2 PARAMETER.
C It has been modified to include the effect of the consumers model
C using the GMTD model. See EQ. 3.25 for explanation of symbols
C used in flow rate equation.
89
Previous Page
