c
c ---------------------------------------------------------------
c
      subroutine add(mx,my,g,mmx,mmy,gg)
c
c prolongation step - add correction
c high-order interpolation
c
      implicit   real*8 (a-h,o-z) 
      parameter  (inode=257,jnode=129)
      dimension  g(mx,my),gg(mmx,mmy)
      common/bi/ mode,itmax,ires,jres,ierr,jerr,icor,jcor
      common/br/ resmax,errmax,cormax
      common/g/  cor(inode,jnode)
c
      nx         = mx  -1
      ny         = my  -1
c
      do 1 i=1,mx
      do 1 j=1,my
    1 cor(i,j)   = 0.
c
c calculate correction at common points
c
      i          = 1
      do 10 ii=2,mmx-1
      i          = i  +2
      j          = 1
      do 10 jj=2,mmy-1
      j          = j  +2
   10 cor(i,j)   = gg(ii,jj)  -g(i,j)
c
c interpolate data along common i-lines
c
      do 200 i=3,mx-2,2
      j          = 2
      cor(i,j)   = .5*(cor(i,j-1)  +cor(i,j+1))
      if (ny .eq. 2) go to 200
      j          = 2
      cor(i,j)   = .375*cor(i,j-1)  +.75*cor(i,j+1)
     .             -.125*cor(i,j+3)
      j          = my  -1
      cor(i,j)   = .375*cor(i,j+1)  +.75*cor(i,j-1)
     .             -.125*cor(i,j-3)
      if (ny .eq. 4) go to 200
      do 210 j=4,my-3,2
  210 cor(i,j)   = .5625*(cor(i,j-1)  +cor(i,j+1))
     .            -.0625*(cor(i,j-3)  +cor(i,j+3))
  200 continue
c
c interpolate remaining points
c
      i          = 2
      do 220 j=2,my-1
  220 cor(i,j)   = .5*(cor(i+1,j)  +cor(i-1,j))
      if (nx .eq. 2) go to 260
      i          = 2
      do 230 j=2,my-1
  230 cor(i,j)   = .375*cor(i-1,j)  +.75*cor(i+1,j)
     .             -.125*cor(i+3,j)
      i          = mx  -1
      do 240 j=2,my-1
  240 cor(i,j)   = .375*cor(i+1,j)  +.75*cor(i-1,j)
     .             -.125*cor(i-3,j)
      if (nx .eq. 4) go to 260
      do 250 i=4,mx-3,2
      do 250 j=2,my-1
  250 cor(i,j)   = .5625*(cor(i+1,j)  +cor(i-1,j))
     .            -.0625*(cor(i+3,j)  +cor(i-3,j))
  260 continue
c
c correct fine mesh solution
c
      do 310 i=2,mx-1
      do 310 j=2,my-1
  310 g(i,j)   = g(i,j)  +cor(i,j)  
c
      return
      end
c
c ---------------------------------------------------------------
c
      subroutine add1(mx,my,g,mmx,mmy,gg)
c
c prolongation step - add correction
c linear interpolation
c
      implicit   real*8 (a-h,o-z)
      parameter  (inode=257,jnode=129)
      dimension  g(mx,my),gg(mmx,mmy)
      common/bi/ mode,itmax,ires,jres,ierr,jerr,icor,jcor
      common/br/ resmax,errmax,cormax
      common/g/  cor(inode,jnode)
c
      nx         = mx  -1
      ny         = my  -1
c
      do 1 j=1,my
      do 1 i=1,mx
    1 cor(i,j)   = 0.
c
c calculate correction at common points
c
      j          = 1
      do 10 jj=2,mmy-1
      j          = j  +2
      i          = 1
      do 10 ii=2,mmx-1
      i          = i  +2
   10 cor(i,j)   = gg(ii,jj)  -g(i,j)
c
c interpolate data along common i-lines
c
      do 100 j=2,my-1,2
      do 100 i=3,mx-2,2
  100 cor(i,j)   = (cor(i,j+1)  +cor(i,j-1))/2.
c
c interpolate data along non-common i-lines
c
      do 110 j=1,my
      do 110 i=2,mx-1,2
  110 cor(i,j)   = (cor(i+1,j)  +cor(i-1,j))/2.
c
c correct fine mesh solution
c
      do 300 j=2,my-1
      do 300 i=2,mx-1
  300 g(i,j)   = g(i,j)  +cor(i,j)
c
      return
      end
c
c ---------------------------------------------------------------
c
      subroutine inject(mx,my,f,mmx,mmy,ff)
c
c injection step - moving vector from fine mesh to coarse mesh
c
      implicit   real*8 (a-h,o-z)
      dimension  f(mx,my),ff(mmx,mmy)
c
c simple injection
c
      j          = -1
      do 10 jj=1,mmy
      j          = j  +2
      i          = -1
      do 10 ii=1,mmx
      i          = i  +2
   10 ff(ii,jj)  = f(i,j)
c
      return
      end
