*
*     F02WGF Example Program Text
*     Mark 22 Release. NAG Copyright 2008.
*     .. Parameters ..
      INTEGER          NIN, NOUT
      PARAMETER        (NIN=5,NOUT=6)
      INTEGER          MAXM, MAXN, MAXNCV
      PARAMETER        (MAXM=500,MAXN=500,MAXNCV=25)
      INTEGER          LDU, LDV
      PARAMETER        (LDU=MAXM,LDV=MAXN)
*     .. Local Scalars ..
      INTEGER          I, IFAIL, K, M, N, NCONV, NCV
*     .. External Subroutines ..
      EXTERNAL         AV, F02WGF
*     .. Local Arrays ..
      DOUBLE PRECISION RESID(MAXNCV), RUSER(1), SIGMA(MAXNCV),
     +                 U(LDU,MAXNCV), V(LDV,MAXNCV)
      INTEGER          IUSER(1)
*     .. Executable Statements ..
      WRITE (NOUT,*) 'F02WGF Example Program Results'
      WRITE (NOUT,*)
*     Skip heading in data file
      READ (NIN,*)
      READ (NIN,*) M, N, K, NCV
      IF (N.LT.1 .OR. N.GT.MAXN) THEN
         WRITE (NOUT,99999) 'N is out of range: N = ', N
      ELSE IF (M.LT.1 .OR. M.GT.MAXM) THEN
         WRITE (NOUT,99999) 'M is out of range: M = ', M
      ELSE IF (NCV.GT.MAXNCV) THEN
         WRITE (NOUT,99999) 'NCV is out of range: NCV = ', NCV
      ELSE
*        Initialize for problem.
         IFAIL = 1
         CALL F02WGF(M,N,K,NCV,AV,NCONV,SIGMA,U,LDU,V,LDV,RESID,IUSER,
     +               RUSER,IFAIL)
         IF (IFAIL.EQ.0) THEN
*
*           Print computed residuals
*
            WRITE (NOUT,*) '  Singular Value    Residual'
            WRITE (NOUT,99998) (SIGMA(I),RESID(I),I=1,NCONV)
         ELSE
            WRITE (NOUT,99997) IFAIL
         END IF
      END IF
*
99999 FORMAT (1X,A,I5)
99998 FORMAT (1X,F10.5,8X,G10.2)
99997 FORMAT (1X,' ** F02WGF Returned with IFAIL = ',I5)
      END
*
*     Matrix vector subroutines
*
      SUBROUTINE AV(IFLAG,M,N,X,AX,IUSER,RUSER)
*
*     Computes  w <- A*x or w <- Trans(A)*x.
*
*     .. Parameters ..
      DOUBLE PRECISION ONE, ZERO
      PARAMETER        (ONE=1.0D+0,ZERO=0.0D+0)
*     .. Scalar Arguments ..
      INTEGER          IFLAG, M, N
*     .. Array Arguments ..
      DOUBLE PRECISION AX(*), RUSER(1), X(*)
      INTEGER          IUSER(1)
*     .. Local Scalars ..
      DOUBLE PRECISION H, K, S, T
      INTEGER          I, J
*     .. Intrinsic Functions ..
      INTRINSIC        DBLE, MIN
*     .. Executable Statements ..
      CONTINUE
      H = ONE/DBLE(M+1)
      K = ONE/DBLE(N+1)
      IF (IFLAG.EQ.1) THEN
         DO 20 I = 1, M
            AX(I) = ZERO
   20    CONTINUE
         T = ZERO
*
         DO 80 J = 1, N
            T = T + K
            S = ZERO
            DO 40 I = 1, MIN(J,M)
               S = S + H
               AX(I) = AX(I) + K*S*(T-ONE)*X(J)
   40       CONTINUE
            DO 60 I = J + 1, M
               S = S + H
               AX(I) = AX(I) + K*T*(S-ONE)*X(J)
   60       CONTINUE
   80    CONTINUE
      ELSE
         DO 100 I = 1, N
            AX(I) = ZERO
  100    CONTINUE
         T = ZERO
*
         DO 160 J = 1, N
            T = T + K
            S = ZERO
            DO 120 I = 1, MIN(J,M)
               S = S + H
               AX(J) = AX(J) + K*S*(T-ONE)*X(I)
  120       CONTINUE
            DO 140 I = J + 1, M
               S = S + H
               AX(J) = AX(J) + K*T*(S-ONE)*X(I)
  140       CONTINUE
  160    CONTINUE
      END IF
*
      RETURN
      END