home *** CD-ROM | disk | FTP | other *** search
/ Power-Programmierung / CD1.mdf / fortran / library / ssp / sequence / teas.for < prev    next >
Encoding:
Text File  |  1985-11-29  |  4.8 KB  |  205 lines

  1. C
  2. C     ..................................................................
  3. C
  4. C        SUBROUTINE TEAS
  5. C
  6. C        PURPOSE
  7. C           CALCULATE THE LIMIT OF A GIVEN SEQUENCE BY MEANS OF THE
  8. C           EPSILON-ALGORITHM.
  9. C
  10. C        USAGE
  11. C           CALL TEAS(X,N,FIN,EPS,IER)
  12. C
  13. C        DESCRIPTION OF PARAMETERS
  14. C           X      - VECTOR WHOSE COMPONENTS ARE TERMS OF THE GIVEN
  15. C                    SEQUENCE. ON RETURN THE COMPONENTS OF VECTOR X
  16. C                    ARE DESTROYED.
  17. C           N      - DIMENSION OF INPUT VECTOR X.
  18. C           FIN    - RESULTANT SCALAR CONTAINING ON RETURN THE LIMIT
  19. C                    OF THE GIVEN SEQUENCE.
  20. C           EPS    - AN INPUT VALUE, WHICH SPECIFIES THE UPPER BOUND
  21. C                    OF THE RELATIVE (ABSOLUTE) ERROR IF THE COMPONENTS
  22. C                    OF X ARE ABSOLUTELY GREATER (LESS) THAN ONE.
  23. C                    CALCULATION IS TERMINATED AS SOON AS THREE TIMES IN
  24. C                    SUCCESSION THE RELATIVE (ABSOLUTE) DIFFERENCE
  25. C                    BETWEEN NEIGHBOURING TERMS IS NOT GREATER THAN EPS.
  26. C           IER    - RESULTANT ERROR PARAMETER CODED IN THE FOLLOWING
  27. C                    FORM
  28. C                     IER=0  - NO ERROR
  29. C                     IER=1  - REQUIRED ACCURACY NOT REACHED WITH
  30. C                              MAXIMAL NUMBER OF ITERATIONS
  31. C                     IER=-1 - INTEGER N IS LESS THAN TEN.
  32. C
  33. C        REMARKS
  34. C           NO ACTION BESIDES ERROR MESSAGE IN CASE N LESS THAN TEN.
  35. C           THE CHARACTER OF THE GIVEN INFINITE SEQUENCE MUST BE
  36. C           RECOGNIZABLE BY THOSE N COMPONENTS OF THE INPUT VECTOR X.
  37. C
  38. C        SUBROUTINES AND FUNCTION SUBPROGRAMS REQUIRED
  39. C           NONE
  40. C
  41. C        METHOD
  42. C           THE CONVERGENCE OF THE GIVEN SEQUENCE IS ACCELERATED BY
  43. C           MEANS OF THE E(2)-TRANSFORMATION, USED IN AN ITERATIVE WAY.
  44. C           FOR REFERENCE, SEE
  45. C           ALGORITHM 215,SHANKS, CACM 1963, NO. 11, PP. 662. AND
  46. C           P. WYNN, SINGULAR RULES FOR CERTAIN NON-LINEAR ALGORITHMS
  47. C           BIT VOL. 3, 1963, PP. 175-195.
  48. C
  49. C     ..................................................................
  50. C
  51.       SUBROUTINE TEAS(X,N,FIN,EPS,IER)
  52. C
  53.       DIMENSION X(1)
  54. C
  55. C        TEST ON WRONG INPUT PARAMETER N
  56. C
  57.       NEW=N
  58.       IF(NEW-10)1,2,2
  59.     1 IER=-1
  60.       RETURN
  61. C
  62. C        CALCULATE INITIAL VALUES FOR THE EPSILON ARRAY
  63. C
  64.     2 ISW1=0
  65.       ISW2=0
  66.       W1=1.E38
  67.       W7=X(4)-X(3)
  68.       IF(W7)3,4,3
  69.     3 W1=1./W7
  70. C
  71.     4 W5=1.E38
  72.       W7=X(2)-X(1)
  73.       IF(W7)5,6,5
  74.     5 W5=1./W7
  75. C
  76.     6 W4=X(3)-X(2)
  77.       IF(W4)9,7,9
  78.     7 W4=1.E38
  79.       T=X(2)
  80.       W2=X(3)
  81.     8 W3=1.E38
  82.       GO TO 17
  83. C
  84.     9 W4=1./W4
  85. C
  86.       T=1.E38
  87.       W7=W4-W5
  88.       IF(W7)10,11,10
  89.    10 T=X(2)+1./W7
  90. C
  91.    11 W2=W1-W4
  92.       IF(W2)15,12,15
  93.    12 W2=1.E38
  94.       IF(T-1.E38)13,14,14
  95.    13 ISW2=1
  96.    14 W3=W4
  97.       GO TO 17
  98. C
  99.    15 W2=X(3)+1./W2
  100.       W7=W2-T
  101.       IF(W7)16,8,16
  102.    16 W3=W4+1./W7
  103. C
  104.    17 ISW1=ISW2
  105.       ISW2=0
  106.       IMIN=4
  107. C
  108. C        CALCULATE DIAGONALS OF THE EPSILON ARRAY IN A DO-LOOP
  109. C
  110.       DO 40 I=5,NEW
  111.       IAUS=I-IMIN
  112.       W4=1.E38
  113.       W5=X(I-1)
  114.       W7=X(I)-X(I-1)
  115.       IF(W7)18,24,18
  116.    18 W4=1./W7
  117. C
  118.       IF(W1-1.E38)19,25,25
  119.    19 W6=W4-W1
  120. C
  121. C        TEST FOR NECESSITY OF A SINGULAR RULE
  122. C
  123.       IF(ABS(W6)-ABS(W4)*1.E-4)20,20,22
  124.    20 ISW2=1
  125.       IF(W6)22,21,22
  126.    21 W5=1.E38
  127.       W6=W1
  128.       IF(W2-1.E38)28,26,26
  129.    22 W5=X(I-1)+1./W6
  130. C
  131. C        FIRST TEST FOR LOSS OF SIGNIFICANCE
  132. C
  133.       IF(ABS(W5)-ABS(X(I-1))*1.E-5)23,24,24
  134.    23 IF(W5)36,24,36
  135. C
  136.    24 W7=W5-W2
  137.       IF(W7)27,25,27
  138.    25 W6=1.E38
  139.    26 ISW2=0
  140.       X(IAUS)=W2
  141.       GO TO 37
  142.    27 W6=W1+1./W7
  143.    28 IF(ISW1-1)33,29,29
  144. C
  145. C        CALCULATE X(IAUS) WITH HELP OF SINGULAR RULE
  146. C
  147.    29 IF(W2-1.E38)30,32,32
  148.    30 W7=W5/(W2-W5)+T/(W2-T)+X(I-2)/(X(I-2)-W2)
  149.       IF(1.+W7)31,38,31
  150.    31 X(IAUS)=W7*W2/(1.+W7)
  151.       GO TO 39
  152. C
  153.    32 X(IAUS)=W5+T-X(I-2)
  154.       GO TO 39
  155. C
  156.    33 W7=W6-W3
  157.       IF(W7)34,38,34
  158.    34 X(IAUS)=W2+1./W7
  159. C
  160. C        SECOND TEST FOR LOSS OF SIGNIFICANCE
  161. C
  162.       IF(ABS(X(IAUS))-ABS(W2)*1.E-5)35,37,37
  163.    35 IF(X(IAUS))36,37,36
  164. C
  165.    36 NEW=IAUS-1
  166.       ISW2=0
  167.       GO TO 41
  168. C
  169.    37 IF(W2-1.E38)39,38,38
  170.    38 X(IAUS)=1.E38
  171.       IMIN=I
  172. C
  173.    39 W1=W4
  174.       T=W2
  175.       W2=W5
  176.       W3=W6
  177.       ISW1=ISW2
  178.    40 ISW2=0
  179. C
  180.       NEW=NEW-IMIN
  181. C
  182. C        TEST FOR ACCURACY
  183. C
  184.    41 IEND=NEW-1
  185.       DO 47 I=1,IEND
  186.       W1=ABS(X(I)-X(I+1))
  187.       W2=ABS(X(I+1))
  188.       IF(W1-EPS)44,44,42
  189.    42 IF(W2-1.)46,46,43
  190.    43 IF(W1-EPS*W2)44,44,46
  191.    44 ISW2=ISW2+1
  192.       IF(3-ISW2)45,45,47
  193.    45 FIN=X(I)
  194.       IER=0
  195.       RETURN
  196. C
  197.    46 ISW2=0
  198.    47 CONTINUE
  199. C
  200.       IF(NEW-6)48,2,2
  201.    48 FIN=X(NEW)
  202.       IER=1
  203.       RETURN
  204.       END
  205.