SUBROUTINE DW7ZBF (L, N, S, W, Y, Z)
C
C  ***  COMPUTE  Y  AND  Z  FOR  DL7UPD  CORRESPONDING TO BFGS UPDATE.
C
      INTEGER N
      DOUBLE PRECISION L(1), S(N), W(N), Y(N), Z(N)
C     DIMENSION L(N*(N+1)/2)
C
C--------------------------  PARAMETER USAGE  --------------------------
C
C L (I/O) CHOLESKY FACTOR OF HESSIAN, A LOWER TRIANG. MATRIX STORED
C             COMPACTLY BY ROWS.
C N (INPUT) ORDER OF  L  AND LENGTH OF  S,  W,  Y,  Z.
C S (INPUT) THE STEP JUST TAKEN.
C W (OUTPUT) RIGHT SINGULAR VECTOR OF RANK 1 CORRECTION TO L.
C Y (INPUT) CHANGE IN GRADIENTS CORRESPONDING TO S.
C Z (OUTPUT) LEFT SINGULAR VECTOR OF RANK 1 CORRECTION TO L.
C
C-------------------------------  NOTES  -------------------------------
C
C  ***  ALGORITHM NOTES  ***
C
C        WHEN  S  IS COMPUTED IN CERTAIN WAYS, E.G. BY  GQTSTP  OR
C     DBLDOG,  IT IS POSSIBLE TO SAVE N**2/2 OPERATIONS SINCE  (L**T)*S
C     OR  L*(L**T)*S IS THEN KNOWN.
C        IF THE BFGS UPDATE TO L*(L**T) WOULD REDUCE ITS DETERMINANT TO
C     LESS THAN EPS TIMES ITS OLD VALUE, THEN THIS ROUTINE IN EFFECT
C     REPLACES  Y  BY  THETA*Y + (1 - THETA)*L*(L**T)*S,  WHERE  THETA
C     (BETWEEN 0 AND 1) IS CHOSEN TO MAKE THE REDUCTION FACTOR = EPS.
C
C  ***  GENERAL  ***
C
C     CODED BY DAVID M. GAY (FALL 1979).
C     THIS SUBROUTINE WAS WRITTEN IN CONNECTION WITH RESEARCH SUPPORTED
C     BY THE NATIONAL SCIENCE FOUNDATION UNDER GRANTS MCS-7600324 AND
C     MCS-7906671.
C
C------------------------  EXTERNAL QUANTITIES  ------------------------
C
C  ***  FUNCTIONS AND SUBROUTINES CALLED  ***
C
      DOUBLE PRECISION DD7TPR
      EXTERNAL DD7TPR, DL7IVM, DL7TVM
C DD7TPR RETURNS INNER PRODUCT OF TWO VECTORS.
C DL7IVM MULTIPLIES L**-1 TIMES A VECTOR.
C DL7TVM MULTIPLIES L**T TIMES A VECTOR.
C
C  ***  INTRINSIC FUNCTIONS  ***
C/+
      DOUBLE PRECISION DSQRT
C/
C--------------------------  LOCAL VARIABLES  --------------------------
C
      INTEGER I
      DOUBLE PRECISION CS, CY, EPS, EPSRT, ONE, SHS, YS, THETA
C
C  ***  DATA INITIALIZATIONS  ***
C
C/6
C     DATA EPS/0.1D+0/, ONE/1.D+0/
C/7
      PARAMETER (EPS=0.1D+0, ONE=1.D+0)
C/
C
C+++++++++++++++++++++++++++++++  BODY  ++++++++++++++++++++++++++++++++
C
      CALL DL7TVM(N, W, L, S)
      SHS = DD7TPR(N, W, W)
      YS = DD7TPR(N, Y, S)
      IF (YS .GE. EPS*SHS) GO TO 10
         THETA = (ONE - EPS) * SHS / (SHS - YS)
         EPSRT = DSQRT(EPS)
         CY = THETA / (SHS * EPSRT)
         CS = (ONE + (THETA-ONE)/EPSRT) / SHS
         GO TO 20
 10   CY = ONE / (DSQRT(YS) * DSQRT(SHS))
      CS = ONE / SHS
 20   CALL DL7IVM(N, Z, L, Y)
      DO 30 I = 1, N
 30      Z(I) = CY * Z(I)  -  CS * W(I)
C
 999  RETURN
C  ***  LAST CARD OF DW7ZBF FOLLOWS  ***
      END