MINICERN/mathlib/gen/g/vavzro.F

   1 *
   2 * $Id$
   3 *
   4 * $Log$
   5 * Revision 1.1.1.1  1996/04/01 15:02:48  mclareni
   6 * Mathlib gen
   7 *
   8 *
   9 #include "gen/pilot.h"
  10       SUBROUTINE VAVZRO(A,B,X,RKA,BE2,LU)
  11 C
  12 C     RZERO SEARCHES FOR THE ROOT OF THE EQUATION FCN=0 IN THE INTERVAL
  13 C
  14       COMMON /FORFCN/ SS,LFCN
  15       DATA E,EPSI,MAXFUN/1E-9,1E-5,100/
  16       MC=0
  17 C
  18       XA=MIN(A,B)
  19       XB=MAX(A,B)
  20       FA=VAVFCN(A,RKA,BE2)
  21       MC=MC+1
  22       FB=VAVFCN(B,RKA,BE2)
  23       IF(FA*FB .GT. 0.0) GO TO 16
  24       MC=MC+1
  25 C
  26     4 X=0.5*(XA+XB)
  27       R=X-XA
  28       EE=ABS(X)+E
  29       IF(R .LE. EE*EPSI) GO TO 18
  30       F1=FA
  31       X1=XA
  32       F2=FB
  33       X2=XB
  34     1 MC=MC+1
  35       G=VAVFCN(X,RKA,BE2)
  36       IF(MC .GT. MAXFUN) GO TO 17
  37       FX=G
  38 C
  39       IF(FX*FA .GT. 0.0) GO TO 2
  40       FB=FX
  41       XB=X
  42       GO TO 3
  43     2 XA=X
  44       FA=FX
  45 C
  46 C     PARABOLA ITERATION
  47 C
  48     3 IF((X1-X2)*(X2-X)*(X1-X) .EQ. 0.0) GO TO 4
  49       F3=FX
  50       X3=X
  51       U1=(F1-F2)/(X1-X2)
  52       U2=(F2-FX)/(X2-X)
  53       CA=U1-U2
  54       CB=(X1+X2)*U2-(X2+X)*U1
  55       CC=(X1-X)*F1-X1*(CA*X1+CB)
  56       IF(CA .EQ. 0.0) GO TO 8
  57       U3=0.5*CB/CA
  58       U4=U3**2-CC/CA
  59       IF(U4 .LT. 0.0) GO TO 4
  60       U5=SQRT(U4)
  61       IF(X .GE. -U3) GO TO 10
  62       X=-U3-U5
  63       GO TO 9
  64    10 X=-U3+U5
  65       GO TO 9
  66     8 X=-CC/CB
  67     9 IF(X .LT. XA) GO TO 4
  68       IF(X .GT. XB) GO TO 4
  69 C
  70 C     TEST FOR OUTPUT
  71 C
  72       R=ABS(X-X3)
  73       R1=ABS(X-X2)
  74       IF(R .GT. R1) R=R1
  75       EE=ABS(X)+E
  76       IF(R/EE .GT. EPSI) GO TO 5
  77       MC=MC+1
  78       G=VAVFCN(X,RKA,BE2)
  79       IF(MC .GT. MAXFUN) GO TO 17
  80       FX=G
  81       IF(FX .EQ. 0.0) GO TO 18
  82       IF(FX*FA .LT. 0.0) GO TO 7
  83       XX=X+EPSI*EE
  84       IF(XX .GE. XB) GO TO 18
  85       MC=MC+1
  86       G=VAVFCN(X,RKA,BE2)
  87       IF(MC .GT. MAXFUN) GO TO 17
  88       FF=VAVFCN(XX,RKA,BE2)
  89       FA=FF
  90       XA=XX
  91       GO TO 6
  92     7 XX=X-EPSI*EE
  93       IF(XX .LE. XA) GO TO 18
  94       MC=MC+1
  95       FX=G
  96       IF(MC .GT. MAXFUN) GO TO 17
  97       FF=VAVFCN(XX,RKA,BE2)
  98       FB=FF
  99       XB=XX
 100     6 IF(FX*FF .GT. 0.0) GO TO 14
 101    18 R=EPSI*EE
 102       FF=VAVFCN(X,RKA,BE2)
 103       RETURN
 104    14 F1=F2
 105       X1=X2
 106       F2=FX
 107       X2=X
 108       X=XX
 109       FX=FF
 110       GO TO 3
 111 C
 112     5 F1=F2
 113       X1=X2
 114       F2=F3
 115       X2=X3
 116       GO TO 1
 117 C
 118    16 WRITE(LU,301)
 119       RETURN
 120 C
 121    17 WRITE(LU,300) X,G,LFCN
 122   301 FORMAT(/10X,' RZERO   FCN(A) AND FCN(B) HAVE SAME SIGN'/)
 123   300 FORMAT(/10X,' RZERO   NUMBER OF ITERATIONS EXCEEDED'/
 124      1  10X,' X=',E15.5,2X,' FCN(X)=',E15.5,2X,' LFCN=',I2/)
 125       RETURN
 126 C
 127       END