MINICERN/mathlib/gen/e/e406cod.inc

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1996/04/01 15:02:28  mclareni
   6 * Mathlib gen
   7 *
   8 *
   9 *
  10 * e406cod.inc
  11 *
  12       DIMENSION C(0:128),CC(0:128),AC(0:64),BC(0:64)
  13
  14       PARAMETER (Z1 = 1, Z2 = 2, HF = Z1/2)
  15
  16       ALFA=HF*(B-A)
  17       BETA=HF*(B+A)
  18       C1=F(A)
  19       C2=F(B)
  20       AC(0)=C2+C1
  21       AC(1)=C2-C1
  22
  23       DO 1 I = 1,7
  24       I1=2**(I-1)
  25       I2=I1-1
  26       I3=2*I1
  27       C1=Z2/I1
  28       C2=PI/I1
  29
  30       DO 2 J = 0,I2
  31     2 C(J)=F(ALFA*COS((J+HF)*C2)+BETA)
  32
  33 C    COMPUTE B-COEFFICIENTS
  34
  35       DO 3 J = 0,I2
  36       F1=J*C2
  37       F2=-HF*F1
  38
  39       C3=2*COS(F1)
  40       A2=0
  41       A1=0
  42       A0=C(I2)
  43       DO 4 K = I2-1,0,-1
  44       A2=A1
  45       A1=A0
  46     4 A0=C(K)+C3*A1-A2
  47     3 BC(J)=C1*(A0*COS(F1+F2)-A1*COS(F2))
  48       BC(I1)=0
  49
  50 C    COMPUTE NEW C-COEFFICIENTS
  51
  52       DO 5 J = 0,I1
  53       C(J)=HF*(AC(J)+BC(J))
  54     5 C(I1+J)=HF*(AC(I1-J)-BC(I1-J))
  55
  56 C    CHECK IF CALCULATION IS FINISHED
  57
  58       CMX=0
  59       DO 6 J = 0,I3
  60       CC(J)=ABS(C(J))
  61     6 CMX=MAX(CMX,CC(J))
  62
  63       IF(CMX .GT. 0) THEN
  64        CMX=1/CMX
  65        CC(I3)=HF*CC(I3)
  66        A0=CC(I2)*CMX
  67        A1=CC(I1)*CMX
  68        DO 10 J = I1+2,I3
  69        A2=CC(J)*CMX
  70        IF(A0 .LE. EPS .AND. A1. LE. EPS .AND. A2 .LE. EPS) GO TO 9
  71
  72        A0=A1
  73    10  A1=A2
  74       ENDIF
  75
  76 C    DOUBLE THE NUMBER OF COEFFICIENTS.
  77
  78       IF(I .EQ. 7) GO TO 1
  79
  80       DO 12 J = 0,I3
  81    12 AC(J)=C(J)
  82
  83     1 CONTINUE
  84
  85 C    REQUIRED ACCURACY NOT OBTAINED
  86
  87       NC=64
  88       DELTA=0
  89       DO 13 J = 60,NC
  90    13 DELTA=DELTA+ABS(C(J))
  91       C(0)=HF*C(0)
  92       CALL MTLPRT(NAME,'E406.1','REQUIRED ACCURACY NOT OBTAINED')
  93       RETURN
  94
  95 C    REQUIRED ACCURACY OBTAINED
  96 C    SUM NEGLECTED TERMS IN EXPANSION
  97
  98     9 DELTA=0
  99       DO 7 K = J,I3
 100     7 DELTA=DELTA+CC(K)
 101
 102 C    CHECK IF FURTHER REDUCTION OF COEFFICIENTS IS POSSIBLE.
 103
 104       NC=J-1
 105       REST=EPS-DELTA
 106       IF(REST .GT. 0) THEN
 107     8  IF(CC(NC) .GE. REST) GO TO 14
 108        DELTA=DELTA+CC(NC)
 109        REST=REST-CC(NC)
 110        NC=NC-1
 111        GO TO 8
 112       ENDIF
 113    14 C(0)=HF*C(0)
 114       RETURN
 115       END