SUBROUTINE SCHDD(R,LDR,P,X,Z,LDZ,NZ,Y,RHO,C,S,INFO)
      INTEGER LDR,P,LDZ,NZ,INFO
      REAL R(LDR,1),X(1),Z(LDZ,1),Y(1),S(1)
      REAL RHO(1),C(1)
C
C     SCHDD DOWNDATES AN AUGMENTED CHOLESKY DECOMPOSITION OR THE
C     TRIANGULAR FACTOR OF AN AUGMENTED QR DECOMPOSITION.
C     SPECIFICALLY, GIVEN AN UPPER TRIANGULAR MATRIX R OF ORDER P,  A
C     ROW VECTOR X, A COLUMN VECTOR Z, AND A SCALAR Y, SCHDD
C     DETERMINEDS A ORTHOGONAL MATRIX U AND A SCALAR ZETA SUCH THAT
C
C                        (R   Z )     (RR  ZZ)
C                    U * (      )  =  (      ) ,
C                        (0 ZETA)     ( X   Y)
C
C     WHERE RR IS UPPER TRIANGULAR.  IF R AND Z HAVE BEEN OBTAINED
C     FROM THE FACTORIZATION OF A LEAST SQUARES PROBLEM, THEN
C     RR AND ZZ ARE THE FACTORS CORRESPONDING TO THE PROBLEM
C     WITH THE OBSERVATION (X,Y) REMOVED.  IN THIS CASE, IF RHO
C     IS THE NORM OF THE RESIDUAL VECTOR, THEN THE NORM OF
C     THE RESIDUAL VECTOR OF THE DOWNDATED PROBLEM IS
C     SQRT(RHO**2 - ZETA**2). SCHDD WILL SIMULTANEOUSLY DOWNDATE
C     SEVERAL TRIPLETS (Z,Y,RHO) ALONG WITH R.
C     FOR A LESS TERSE DESCRIPTION OF WHAT SCHDD DOES AND HOW
C     IT MAY BE APPLIED, SEE THE LINPACK GUIDE.
C
C     THE MATRIX U IS DETERMINED AS THE PRODUCT U(1)*...*U(P)
C     WHERE U(I) IS A ROTATION IN THE (P+1,I)-PLANE OF THE
C     FORM
C
C                       ( C(I)     -S(I)     )
C                       (                    ) .
C                       ( S(I)       C(I)    )
C
C     THE ROTATIONS ARE CHOSEN SO THAT C(I) IS REAL.
C
C     THE USER IS WARNED THAT A GIVEN DOWNDATING PROBLEM MAY
C     BE IMPOSSIBLE TO ACCOMPLISH OR MAY PRODUCE
C     INACCURATE RESULTS.  FOR EXAMPLE, THIS CAN HAPPEN
C     IF X IS NEAR A VECTOR WHOSE REMOVAL WILL REDUCE THE
C     RANK OF R.  BEWARE.
C
C     ON ENTRY
C
C         R      REAL(LDR,P), WHERE LDR .GE. P.
C                R CONTAINS THE UPPER TRIANGULAR MATRIX
C                THAT IS TO BE DOWNDATED.  THE PART OF  R
C                BELOW THE DIAGONAL IS NOT REFERENCED.
C
C         LDR    INTEGER.
C                LDR IS THE LEADING DIMENSION FO THE ARRAY R.
C
C         P      INTEGER.
C                P IS THE ORDER OF THE MATRIX R.
C
C         X      REAL(P).
C                X CONTAINS THE ROW VECTOR THAT IS TO
C                BE REMOVED FROM R.  X IS NOT ALTERED BY SCHDD.
C
C         Z      REAL(LDZ,NZ), WHERE LDZ .GE. P.
C                Z IS AN ARRAY OF NZ P-VECTORS WHICH
C                ARE TO BE DOWNDATED ALONG WITH R.
C
C         LDZ    INTEGER.
C                LDZ IS THE LEADING DIMENSION OF THE ARRAY Z.
C
C         NZ     INTEGER.
C                NZ IS THE NUMBER OF VECTORS TO BE DOWNDATED
C                NZ MAY BE ZERO, IN WHICH CASE Z, Y, AND RHO
C                ARE NOT REFERENCED.
C
C         Y      REAL(NZ).
C                Y CONTAINS THE SCALARS FOR THE DOWNDATING
C                OF THE VECTORS Z.  Y IS NOT ALTERED BY SCHDD.
C
C         RHO    REAL(NZ).
C                RHO CONTAINS THE NORMS OF THE RESIDUAL
C                VECTORS THAT ARE TO BE DOWNDATED.
C
C     ON RETURN
C
C         R
C         Z      CONTAIN THE DOWNDATED QUANTITIES.
C         RHO
C
C         C      REAL(P).
C                C CONTAINS THE COSINES OF THE TRANSFORMING
C                ROTATIONS.
C
C         S      REAL(P).
C                S CONTAINS THE SINES OF THE TRANSFORMING
C                ROTATIONS.
C
C         INFO   INTEGER.
C                INFO IS SET AS FOLLOWS.
C
C                   INFO = 0  IF THE ENTIRE DOWNDATING
C                             WAS SUCCESSFUL.
C
C                   INFO =-1  IF R COULD NOT BE DOWNDATED.
C                             IN THIS CASE, ALL QUANTITIES
C                             ARE LEFT UNALTERED.
C
C                   INFO = 1  IF SOME RHO COULD NOT BE
C                             DOWNDATED.  THE OFFENDING RHOS ARE
C                             SET TO -1.
C
C     LINPACK. THIS VERSION DATED 08/14/78 .
C     G.W. STEWART, UNIVERSITY OF MARYLAND, ARGONNE NATIONAL LAB.
C
C     SCHDD USES THE FOLLOWING FUNCTIONS AND SUBPROGRAMS.
C
C     FORTRAN ABS
C     BLAS SDOT, SNRM2
C
      INTEGER I,II,J
      REAL A,ALPHA,AZETA,NORM,SNRM2
      REAL SDOT,T,ZETA,B,XX
C
C     SOLVE THE SYSTEM TRANS(R)*A = X, PLACING THE RESULT
C     IN THE ARRAY S.
C
      INFO = 0
      S(1) = X(1)/R(1,1)
      IF (P .LT. 2) GO TO 20
      DO 10 J = 2, P
         S(J) = X(J) - SDOT(J-1,R(1,J),1,S,1)
         S(J) = S(J)/R(J,J)
   10 CONTINUE
   20 CONTINUE
      NORM = SNRM2(P,S,1)
      IF (NORM .LT. 1.0E0) GO TO 30
         INFO = -1
      GO TO 120
   30 CONTINUE
         ALPHA = SQRT(1.0E0-NORM**2)
C
C        DETERMINE THE TRANSFORMATIONS.
C
         DO 40 II = 1, P
            I = P - II + 1
            SCALE = ALPHA + ABS(S(I))
            A = ALPHA/SCALE
            B = S(I)/SCALE
            NORM = SQRT(A**2+B**2+0.0E0**2)
            C(I) = A/NORM
            S(I) = B/NORM
            ALPHA = SCALE*NORM
   40    CONTINUE
C
C        APPLY THE TRANSFORMATIONS TO R.
C
         DO 60 J = 1, P
            XX = 0.0E0
            DO 50 II = 1, J
               I = J - II + 1
               T = C(I)*XX + S(I)*R(I,J)
               R(I,J) = C(I)*R(I,J) - S(I)*XX
               XX = T
   50       CONTINUE
   60    CONTINUE
C
C        IF REQUIRED, DOWNDATE Z AND RHO.
C
         IF (NZ .LT. 1) GO TO 110
         DO 100 J = 1, NZ
            ZETA = Y(J)
            DO 70 I = 1, P
               Z(I,J) = (Z(I,J) - S(I)*ZETA)/C(I)
               ZETA = C(I)*ZETA - S(I)*Z(I,J)
   70       CONTINUE
            AZETA = ABS(ZETA)
            IF (AZETA .LE. RHO(J)) GO TO 80
               INFO = 1
               RHO(J) = -1.0E0
            GO TO 90
   80       CONTINUE
               RHO(J) = RHO(J)*SQRT(1.0E0-(AZETA/RHO(J))**2)
   90       CONTINUE
  100    CONTINUE
  110    CONTINUE
  120 CONTINUE
      RETURN
      END