SUBROUTINE vdifv (lmv, ktm, dtq2, uz0, vz0, &
     akmsvwork, u, v, akmwork, zcolwork)
  !
  ! $Author: CPTEC $
  ! $Date: 2000/02/24 13:15:43 $
  ! $Revision: 1.0 $
  !
  !
  ! vertical diffusion
  !
  IMPLICIT NONE
  INCLUDE "parmeta.f90"
  INTEGER, PARAMETER                                     :: lp1=lm+1
  INTEGER, PARAMETER                                     :: lm1=lm-1
  INTEGER, DIMENSION(  im  ,   jm       ), INTENT(IN   ) :: lmv
  INTEGER,                                 INTENT(IN   ) :: ktm
  REAL   ,                                 INTENT(IN   ) :: dtq2
  REAL   , DIMENSION(0:im+1, 0:jm+1     ), INTENT(IN   ) :: uz0
  REAL   , DIMENSION(0:im+1, 0:jm+1     ), INTENT(IN   ) :: vz0
  REAL   , DIMENSION(  im  ,   jm       ), INTENT(IN   ) :: akmsvwork
  REAL   , DIMENSION(0:im+1, 0:jm+1, lm ), INTENT(INOUT) :: u, v
  REAL   , DIMENSION(  im  ,   jm  , lm1), INTENT(IN   ) :: akmwork
  REAL   , DIMENSION(  im  ,   jm  , lp1), INTENT(IN   ) :: zcolwork

  INTEGER :: i, j, l, kt, lmax
  REAL    :: dtdif, dtozs, dtozak, akmh, rcmvb, rcml, cf
  REAL   , DIMENSION(  im  ,   jm  , lm1) :: dtoz
  REAL   , DIMENSION(  im  ,   jm  , lm1) :: cm
  REAL   , DIMENSION(  im  ,   jm  , lm1) :: cr
  REAL   , DIMENSION(  im  ,   jm  , lm1) :: rsu
  REAL   , DIMENSION(  im  ,   jm  , lm1) :: rsv


  dtdif=dtq2/float(ktm)

  lmax = MAXVAL(lmv(1:im, 3:jm-2))

  DO l = 1, lmax -1
     DO j = 3, jm-2
        DO i = 1, im
           IF (l <= lmv(i,j)-1) THEN
              dtoz(i,j,l)=dtdif/(zcolwork(i,j,l)-zcolwork(i,j,l+1))
              cr(i,j,l)=-dtoz(i,j,l)*akmwork(i,j,l)
           END IF
        END DO
     END DO
  END DO

  DO kt = 1, ktm

     DO j = 3, jm-2
        DO i = 1, im
           cm(i,j,1)=dtoz(i,j,1)*akmwork(i,j,1)+1.
           rsu(i,j,1)=u(i,j,1)
           rsv(i,j,1)=v(i,j,1)
        END DO
     END DO

     DO l = 2, lmax-1
        DO j = 3, jm-2
           DO i = 1, im
              IF (l <= lmv(i,j)-1) THEN
                 cf=-dtoz(i,j,l)*akmwork(i,j,l-1)/cm(i,j,l-1)
                 cm(i,j,l)=-cr(i,j,l-1)*cf &
                      +(akmwork(i,j,l-1)+akmwork(i,j,l))*dtoz(i,j,l)+1.
                 rsu(i,j,l)=-rsu(i,j,l-1)*cf+u(i,j,l)
                 rsv(i,j,l)=-rsv(i,j,l-1)*cf+v(i,j,l)
              END IF
           END DO
        END DO
     END DO

     DO j = 3, jm-2
        DO i = 1, im
           dtozs=dtdif/(zcolwork(i,j,lmv(i,j))-zcolwork(i,j,lmv(i,j)+1))
           akmh=akmwork(i,j,lmv(i,j)-1)
           cf=-dtozs*akmh/cm(i,j,lmv(i,j)-1)
           rcmvb=1./((akmh+akmsvwork(i,j))*dtozs-cr(i,j,lmv(i,j)-1)*cf+1.)
           dtozak=dtozs*akmsvwork(i,j)
           u(i,j,lmv(i,j))=(dtozak*uz0(i,j)-rsu(i,j,lmv(i,j)-1)*cf &
                +u(i,j,lmv(i,j)))*rcmvb
           v(i,j,lmv(i,j))=(dtozak*vz0(i,j)-rsv(i,j,lmv(i,j)-1)*cf &
                +v(i,j,lmv(i,j)))*rcmvb
        END DO
     END DO

     DO l = lmax-1, 1, -1
        DO j = 3, jm-2
           DO i = 1, im
              IF (l <= lmv(i,j) - 1) THEN
                 rcml=1./cm(i,j,l)
                 u(i,j,l)=(-cr(i,j,l)*u(i,j,l+1)+rsu(i,j,l))*rcml
                 v(i,j,l)=(-cr(i,j,l)*v(i,j,l+1)+rsv(i,j,l))*rcml
              END IF
           END DO
        END DO
     END DO
  END DO
END SUBROUTINE vdifv
