C      _______________________________________________
C     |                                               |
C     |    SOLVE  -X''(T) = 1, X(0) = X(1) = 0        |
C     |          BY THE MULTIGRID METHOD              |
C     |_______________________________________________|
      REAL R(518),X(518),D0,D1,S,T,U
      INTEGER I,IT,J,K,L,LL,M,N,NL
      K = 8
      N = 2**K
      IT = 4
      T = .0001
      U = .7
C     --------------------------
C     |*** INPUT RIGHT SIDE ***|
C     --------------------------
      DO 10 I = 1,N
10         R(I) = 1.
C     --------------------
C     |*** INITIALIZE ***|
C     --------------------
      S = (1./N)**2
      DO 20 I = 1,N
20         R(I) = S*R(I)
      M = N
      L = 1
      LL = N
      NL = N + N + K - 2
      DO 30 I = 1,NL
30         X(I) = 0.
      D1 = 0.
40    J = 0
C     --------------------------------
C     |*** GAUSS-SEIDEL ITERATION ***|
C     --------------------------------
50    D0 = D1
      D1 = 0.
      J = J + 1
      I = L
60    I = I + 1
           S = .5*(X(I-1)+X(I+1)+R(I))
           D1 = D1 + ABS(S-X(I))
           X(I) = S
      IF ( I .LT. LL ) GOTO 60
      WRITE(6,70) D1
70    FORMAT(' DIF:',F20.10)
      IF ( J .LT. IT ) GOTO 50
      IF ( D1  .LT.  T ) GOTO 100
      IF ( D1/D0 .LT. U ) GOTO 50
      IF ( N .EQ. 2 ) GOTO 50
C     -----------------------------------------
C     |*** COARSER MESH (SLOW CONVERGENCE) ***|
C     -----------------------------------------
      I = LL + 2
80    L = L + 2
      I = I + 1
      IF ( L .GT. LL ) GOTO 90
           X(I) = 0.
           R(I) = 4*(R(L)+X(L-1)-2*X(L)+X(L+1))
      GOTO 80
90    N = N/2
      WRITE(6,*) 'NUMBER OF MESH INTERVALS:',N
      LL = LL + N + 1
      L = L + 1
      GOTO 40
C     ---------------------------------------
C     |*** FINER MESH (FAST CONVERGENCE) ***|
C     ---------------------------------------
100   IF ( N .EQ. M ) GOTO 120
      I = L - 3
      J = LL
110        X(I) = X(I) + X(J)
           X(I+1) = X(I+1) + .5*(X(J)+X(J+1))
           I = I - 2
           J = J - 1
      IF ( J .GT. L ) GOTO 110
      N = N + N
      WRITE(6,*) 'NUMBER OF MESH INTERVALS:',N
      LL = L - 2
      L = I
      X(I+1) = X(I+1) + .5*X(J+1)
      GOTO 40
C     ----------------------
C     |*** PRINT ANSWER ***|
C     ----------------------
120   M = N + 1
      DO 130 I = 1,M
           J = I - 1
           S = J/FLOAT(N)
130        WRITE(6,140) J,S,X(I)
140   FORMAT(I5,F10.5,F15.6)
      END
