dorth.f

Go to the documentation of this file.

C Work performed under the auspices of the U.S. Department of Energy
C by Lawrence Livermore National Laboratory under contract number
C W-7405-Eng-48.
C
      SUBROUTINE dorth (VNEW, V, HES, N, LL, LDHES, KMP, SNORMW)
C
C***BEGIN PROLOGUE  DORTH
C***DATE WRITTEN   890101   (YYMMDD)
C***REVISION DATE  900926   (YYMMDD)
C
C
C-----------------------------------------------------------------------
C***DESCRIPTION
C
C This routine orthogonalizes the vector VNEW against the previous
C KMP vectors in the V array.  It uses a modified Gram-Schmidt
C orthogonalization procedure with conditional reorthogonalization.
C
C      On entry
C
C         VNEW = The vector of length N containing a scaled product
C                OF The Jacobian and the vector V(*,LL).
C
C         V    = The N x LL array containing the previous LL
C                orthogonal vectors V(*,1) to V(*,LL).
C
C         HES  = An LL x LL upper Hessenberg matrix containing,
C                in HES(I,K), K.LT.LL, scaled inner products of
C                A*V(*,K) and V(*,I).
C
C        LDHES = The leading dimension of the HES array.
C
C         N    = The order of the matrix A, and the length of VNEW.
C
C         LL   = The current order of the matrix HES.
C
C          KMP = The number of previous vectors the new vector VNEW
C                must be made orthogonal to (KMP .LE. MAXL).
C
C
C      On return
C
C         VNEW = The new vector orthogonal to V(*,I0),
C                where I0 = MAX(1, LL-KMP+1).
C
C         HES  = Upper Hessenberg matrix with column LL filled in with
C                scaled inner products of A*V(*,LL) and V(*,I).
C
C       SNORMW = L-2 norm of VNEW.
C
C-----------------------------------------------------------------------
C***ROUTINES CALLED
C   DDOT, DNRM2, DAXPY
C
C***END PROLOGUE  DORTH
C
      INTEGER N, LL, LDHES, KMP
      DOUBLE PRECISION VNEW, V, HES, SNORMW
      dimension vnew(*), v(n,*), hes(ldhes,*)
      INTEGER I, I0
      DOUBLE PRECISION ARG, DDOT, DNRM2, SUMDSQ, TEM, VNRM
C
C-----------------------------------------------------------------------
C Get norm of unaltered VNEW for later use.
C-----------------------------------------------------------------------
      vnrm = dnrm2(n, vnew, 1)
C-----------------------------------------------------------------------
C Do Modified Gram-Schmidt on VNEW = A*V(LL).
C Scaled inner products give new column of HES.
C Projections of earlier vectors are subtracted from VNEW.
C-----------------------------------------------------------------------
      i0 = max0(1,ll-kmp+1)
      DO 10 i = i0,ll
        hes(i,ll) = ddot(n, v(1,i), 1, vnew, 1)
        tem = -hes(i,ll)
        CALL daxpy(n, tem, v(1,i), 1, vnew, 1)
 10     CONTINUE
C-----------------------------------------------------------------------
C Compute SNORMW = norm of VNEW.
C If VNEW is small compared to its input value (in norm), then
C Reorthogonalize VNEW to V(*,1) through V(*,LL).
C Correct if relative correction exceeds 1000*(unit roundoff).
C Finally, correct SNORMW using the dot products involved.
C-----------------------------------------------------------------------
      snormw = dnrm2(n, vnew, 1)
      IF (vnrm + 0.001d0*snormw .NE. vnrm) RETURN
      sumdsq = 0.0d0
      DO 30 i = i0,ll
        tem = -ddot(n, v(1,i), 1, vnew, 1)
        IF (hes(i,ll) + 0.001d0*tem .EQ. hes(i,ll)) go to 30
        hes(i,ll) = hes(i,ll) - tem
        CALL daxpy(n, tem, v(1,i), 1, vnew, 1)
        sumdsq = sumdsq + tem**2
 30     CONTINUE
      IF (sumdsq .EQ. 0.0d0) RETURN
      arg = max(0.0d0,snormw**2 - sumdsq)
      snormw = sqrt(arg)
      RETURN
C
C------END OF SUBROUTINE DORTH------------------------------------------
      END