SUBROUTINE ddamp
  !
  ! $Author: CPTEC $
  ! $Date: 2000/02/24 13:15:43 $
  ! $Revision: 1.0 $
  !
  ! subprogram:    ddamp       divergence damping
  !   prgrmmr: janjic          org: w/np2      date: 94-03-08       
  !     
  ! abstract:
  !     ddamp modifies the wind components so as to reduce the
  !     horizontal divergence.  a switch provides the option of
  !     also modifying the temperature from an energy viewpoint.
  !     
  ! program history log:
  !   87-08-??  janjic     - originator
  !   95-03-25  black      - conversion from 1-d to 2-d in horizontal
  !   95-03-28  black      - added external edge
  !     
  ! usage: call ddamp from main program ebu
  !
  !   input argument list:
  !       none     
  !  
  !   output argument list: 
  !     none
  !     
  !   output files:
  !     none
  !     
  !   subprograms called:
  !  
  !     unique: none
  !  
  !     library: none
  !  
  !   common blocks: ctlblk
  !                  masks
  !                  dynam
  !                  scrch
  !                  vrbls
  !                  indx
  !   
  IMPLICIT NONE

  !*** problem size

  INCLUDE "parmeta.f90"
  INTEGER, PARAMETER :: jam=6+2*(jm-10)
  INTEGER, PARAMETER :: lp1=lm+1

  !*** global variables

  INCLUDE "masks.f90"  ! hbm2(in); vtm(in)
  INCLUDE "dynam.f90"  ! ddmpu(in); ddmpv(in); deta(in)
  INCLUDE "vrbls.f90"  ! t(inout); u(inout); v(inout)
  INCLUDE "contin.f90" ! div(inout); pdsl(in)
  INCLUDE "indx.f90"   ! ihe(in); ihw(in); ive(in); ivw(in)

  !*** local variables

  REAL,    PARAMETER                :: d00=.0e0
  REAL,    PARAMETER                :: d50=.5e0
  REAL,    PARAMETER                :: rfcp=.25e0/1004.6e0
  LOGICAL, PARAMETER                :: heat=.FALSE.
  INTEGER                           :: i, j, l
  REAL,    DIMENSION(im,jm)         :: rdpdx
  REAL,    DIMENSION(im,jm)         :: rdpdy
  REAL,    DIMENSION(im,jm)         :: ut
  REAL,    DIMENSION(im,jm)         :: vt
  REAL,    DIMENSION(0:im+1,0:jm+1) :: cke
  REAL,    DIMENSION(0:im+1,0:jm+1) :: dpde

  DO j=2,jm-1
     DO i=1,im
        cke(i,j)=d00
     END DO
  END DO
  DO l=1,lm
     CALL zero2(cke)
     CALL zero2(dpde)
     DO j=1,jm
        DO i=1,im
           dpde(i,j)=deta(l)*pdsl(i,j)
           div(i,j,l)=div(i,j,l)*hbm2(i,j)
        END DO
     END DO

     DO j=3,jm-2
        DO i=1,im
           rdpdx(i,j)=vtm(i,j,l)/(dpde(i+ivw(j),j)+dpde(i+ive(j),j))
           rdpdy(i,j)=vtm(i,j,l)/(dpde(i,j-1)+dpde(i,j+1))
        END DO
     END DO

     DO j=3,jm-2
        DO i=1,im
           ut(i,j)=u(i,j,l)
           vt(i,j)=v(i,j,l)
           u(i,j,l)=u(i,j,l)+(div(i+ive(j),j,l)-div(i+ivw(j),j,l)) &
                *rdpdx(i,j)*ddmpu(i,j)
           v(i,j,l)=v(i,j,l) &
                +(div(i,j+1,l)-div(i,j-1,l))*rdpdy(i,j)*ddmpv(i,j)
           cke(i,j)=d50*(u(i,j,l)*u(i,j,l)-ut(i,j)*ut(i,j) &
                +v(i,j,l)*v(i,j,l)-vt(i,j)*vt(i,j))
        END DO
     END DO

     IF(heat)THEN
        DO j=3,jm-2
           DO i=1,im
              t(i,j,l)=t(i,j,l)-rfcp*(cke(i+ihe(j),j)+cke(i,j+1) &
                   +cke(i+ihw(j),j)+cke(i,j-1))*hbm2(i,j)
           END DO
        END DO
     ENDIF
  END DO
END SUBROUTINE ddamp
