!+ horizontal diffusion.
      !+ $Id: hdiff.f90,v 1.20 1999/10/18 13:18:45 m214003 Exp $
      !OCL NOALIAS
      
   5: SUBROUTINE hdiff
      
        ! Description:
        !
        ! Horizontal diffusion.
  10:   !
        ! Method:
        !
        ! This subroutine performs horizontal diffusion.
        !
  15:   ! *hdiff* is called from *stepon*
        !
        ! Authors:
        !
        ! U. Schlese, DKRZ, June 1994, original source
  20:   ! L. Kornblueh, MPI, May 1998, f90 rewrite
        ! U. Schulzweida, MPI, May 1998, f90 rewrite
        ! I. Kirchner, MPI, August 1998, tendency diagnostics
        ! T. Diehl, DKRZ, July 1999, parallel version 
        ! 
  25:   ! for more details see file AUTHORS
        !
      
        USE mo_decomposition, ONLY: lc => local_decomposition
        USE mo_memory_sp,     ONLY: sd, stp, svo
  30:   USE mo_control,       ONLY: dtime, lmidatm, ltdiag, nk, nkp1, nlev, nlevp1, &
                                    twodt
        USE mo_truncation,    ONLY: ntrn
        USE mo_semi_impl,     ONLY: hdamp, vcrit, vmax
        USE mo_hdiff,         ONLY: difd, dift, diftcor, difvo, enstdif, ldiahdf, &
  35:                               nlvstd1, nlvstd2
        USE mo_constants,     ONLY: a
        USE mo_start_dataset, ONLY: ndiahdf, nstart, nstep
        USE mo_diff,          ONLY: iq, ncdif
        USE mo_midatm,        ONLY: shigh, damhih
  40:   USE mo_diag_tendency, ONLY: pdvor, pdtem, pddiv
      
        IMPLICIT NONE
      
        !  Local scalars: 
  45:   REAL :: zaa, zcons, zcor, zdamp, zdifd, zdift, zdifvo, znn, ztwodt, zzcor, &
      &      zztemp, zzzn
        INTEGER :: ilev2, in, is, jk, jlev, jn, jndif, jr, ic, i, snsp, nns
      
        !  Local arrays: 
  50:   REAL :: zcor1(nkp1,nlev), zcor2(nkp1,nlev), znd(nkp1,2*nlev), &
      &      zndifd(nkp1,nlev), zndift(nkp1,nlev), zndifvo(nkp1,nlev), &
      &      znt(nkp1,nlev), znvo(nkp1,2*nlev), zzus(nkp1,nlev)
      
        INTEGER :: np1(lc%snsp), nindex(lc%nns)
  55: 
      !DIR$ NOBOUNDS sd,svo
      
        !  Executable statements
      
  60:   snsp   = lc%snsp
        nns    = lc%nns
        np1    = lc%np1
        nindex = lc%nindex
      
  65:   ! Each PE loops only over n's owned by itself 
      
        DO i=1,nns
           jn=nindex(i)
           zzus(jn,:) = 0.
  70:   ENDDO
      
      !-- 1. Compute diffusion correction factor
      
        ztwodt = twodt
  75:   IF (nstep==nstart) ztwodt = ztwodt*.5
      
        ilev2 = 2*nlev
        
        zaa = 1./(a*a)
  80:   zcons = 2.5*dtime/a
      
        DO i=1,nns
           jn=nindex(i)
           IF (jn >= 2) THEN
  85:         in = jn - 1
              zdamp = 1.
              IF (lmidatm) THEN
                 DO jk = nlev, 1, -1
                    IF (in>ntrn(jk)) zdamp = zdamp*hdamp
  90:               IF (jk<=nlvstd2 .AND. jk>=nlvstd1) zdamp = zdamp*enstdif
                    zzus(jn,jk) = zdamp
                 END DO
              ELSE
                 IF (jn<=nk/3) THEN
  95:               DO jk = nlev, 1, -1
                       IF (in>ntrn(jk)) zdamp = zdamp*hdamp
                       IF (jk<=nlvstd2 .AND. jk>=nlvstd1) zdamp = zdamp*enstdif
                       zzus(jn,jk) = zdamp
                    END DO
 100:            ELSE
                    DO jk = 1, nlev
                       zzus(jn,jk) = 1.
                    END DO
                 END IF
 105:         END IF
           ENDIF
        END DO
      
        IF (lmidatm) THEN
 110:      DO i=1,nns
              jn=nindex(i)
              DO jk = 1,nlev
                 shigh(jn,jk)=1.
              END DO
 115:      END DO
      
           IF (.NOT.ldiahdf) THEN
              DO i=1,nns
                 jn=nindex(i)
 120:            IF (jn >= 3) THEN
                    znn = jn-1
                    DO jk = 1,nlev
                       zcor = (1.+(znn*vmax(jk)-vcrit)) 
                       IF (zcor .GT. 1.) shigh(jn,jk) = damhih
 125:               END DO
                 ENDIF
              END DO
           ELSE
            ! same loop but WITH statistics
 130:         DO i=1,nns
                 jn=nindex(i)
                 IF (jn >= 3) THEN
                    znn = jn-1
                    DO jk = 1,nlev
 135:                  zcor = (1.+(znn*vmax(jk)-vcrit)) 
                       IF (zcor .GT. 1.) THEN
                          shigh(jn,jk) = 1000.
                          WRITE(ndiahdf,*) nstep, jk, jn-1, vmax(jk)   
                       ENDIF
 140:               END DO
                 ENDIF
              END DO
           END IF
        END IF
 145: 
        ! Vertical loop
      
        DO jk = 1, nlev
      
 150:      jndif = ncdif(jk)
           zztemp = (nkp1*nk*zaa)**(1-iq(jk))*ztwodt
           zdifd = difd*zztemp
           zdifvo = difvo*zztemp
           zdift = dift*zztemp
 155:      
           DO i=1,nns
              jn=nindex(i)
              zndifd(jn,jk) = 1.
              zndifvo(jn,jk) = 1.
 160:         zndift(jn,jk) = 1.
           END DO
           
           DO i=1,nns
              jn=nindex(i)
 165:         IF (jn >= jndif+1) THEN
                 zzzn = jn - 1 - jndif
                 zztemp = (zaa*zzzn*(zzzn+1.))**iq(jk)*zzus(jn,jk)
                 IF (lmidatm) THEN
                    zndifd(jn,jk)  = 1. + shigh(jn,jk)*zdifd*zztemp
 170:               zndifvo(jn,jk) = 1. + shigh(jn,jk)*zdifvo*zztemp
                    zndift(jn,jk)  = 1. + shigh(jn,jk)*zdift*zztemp
                 ELSE
                    zndifd(jn,jk)  = 1. + zdifd*zztemp
                    zndifvo(jn,jk) = 1. + zdifvo*zztemp
 175:               zndift(jn,jk)  = 1. + zdift*zztemp
                 ENDIF
              ENDIF
           END DO
           
 180:      IF ( .NOT. lmidatm) THEN
              IF ( .NOT. ldiahdf) THEN
                 DO i=1,nns
                    jn=nindex(i)
                    IF (jn >= 3) THEN
 185:                  znn = jn - 1
                       zcor = (1.+zcons*(znn*vmax(jk)-vcrit))
                       IF (zcor>1.) THEN
                          zndifd(jn,jk) = zcor*zndifd(jn,jk)
                          zndifvo(jn,jk) = zcor*zndifvo(jn,jk)
 190:                  END IF
                    ENDIF
                 END DO
              ELSE
                 
 195:            ! Same loop but with statistics
                 
                 DO i=1,nns
                    jn=nindex(i)
                    IF (jn >= 3) THEN
 200:                  znn = jn - 1
                       zcor = (1.+zcons*(znn*vmax(jk)-vcrit))
                       IF (zcor>1.) THEN
                          WRITE (ndiahdf,*) nstep, jk, jn - 1, zcor, zndifd(jn,jk)
                          zndifd(jn,jk) = zcor*zndifd(jn,jk)
 205:                     zndifvo(jn,jk) = zcor*zndifvo(jn,jk)
                       END IF
                    ENDIF
                 END DO
              END IF
 210:      END IF
           
           DO i=1,nns
              jn=nindex(i)
              IF (jn >= 3) THEN
 215:            znd(jn,jk) = 1./zndifd(jn,jk)
                 znd(jn,jk+nlev) = znd(jn,jk)
                 znvo(jn,jk) = 1./zndifvo(jn,jk)
                 znvo(jn,jk+nlev) = znvo(jn,jk)
              ELSE
 220:            znd(jn,jk) = 1.
                 znd(jn,jk+nlev) = 1.
                 znvo(jn,jk) = 1.
                 znvo(jn,jk+nlev) = 1.
              ENDIF
 225:      END DO
           
           ! Temperature
           
           DO i=1,nns
 230:         jn=nindex(i)
              IF (jn >= 2) THEN
                 zzcor = 1.
                 IF (jn-1>ntrn(jk)) zzcor = 0.
                 znn = jn - 1.
 235:            IF (lmidatm) THEN
                    zcor = 1.
                 ELSE
                    zcor = (1.+zcons*(znn*vmax(jk)-vcrit))
                    IF (zcor>1.) zndift(jn,jk) = zcor*zndift(jn,jk)
 240:               IF (zcor<=1.) zcor = 1.
                 END IF
                 zcor2(jn,jk) = zzcor*diftcor(jk)
                 zcor1(jn,jk) = zcor2(jn,jk)/zcor
              ELSE
 245:            zcor1(jn,jk) = 0.
                 zcor2(jn,jk) = 0.
              ENDIF
           END DO
           
 250:      DO i=1,nns
              jn=nindex(i)
              znt(jn,jk) = 1./zndift(jn,jk)
           END DO
        END DO
 255:   
      !-- 2. Modify fields
      
        DO is = 1, snsp
      
 260:      ic = np1(is)
      
            IF (ltdiag) THEN
               ! remove vorticity and divergence without diffusion part
               pddiv(:,:,is,8) = pddiv(:,:,is,8) - sd (:,:,is)
 265:          pdvor(:,:,is,8) = pdvor(:,:,is,8) - svo(:,:,is)
            ENDIF
      
            sd (:,1,is) = sd (:,1,is)*znd (ic,:nlev)
            sd (:,2,is) = sd (:,2,is)*znd (ic, nlev+1:)
 270:       svo(:,1,is) = svo(:,1,is)*znvo(ic,:nlev)
            svo(:,2,is) = svo(:,2,is)*znvo(ic, nlev+1:)
      
            IF (ltdiag) THEN
               ! store diffusion part of vorticity and divergence
 275:          pddiv(:,:,is,8) = pddiv(:,:,is,8) + sd (:,:,is)
               pdvor(:,:,is,8) = pdvor(:,:,is,8) + svo(:,:,is)
            ENDIF
      
            DO jr = 1, 2
 280:          ! remove temperature without diffusion part
               IF (ltdiag) pdtem(:,jr,is,11) = pdtem(:,jr,is,11) - stp(:,jr,is)
      !DIR$ IVDEP
      !CDIR NODEP
      !OCL NOVREC,NOALIAS
 285:         DO jlev = 1, nlev
                stp(jlev,jr,is) = stp(nlevp1,jr,is)*zcor1(ic,jlev) + &
      &              (stp(jlev,jr,is)-zcor2(ic,jlev)*stp(nlevp1,jr,is))* &
      &              znt(ic,jlev)
              END DO
 290: 
              ! add temperature with diffusion part
              IF (ltdiag) pdtem(:,jr,is,11) = pdtem(:,jr,is,11) + stp(:,jr,is)
      
            END DO
 295: 
         END DO
      
        RETURN
      END SUBROUTINE hdiff


Info Section
         uses:  mo_constants, mo_control, mo_decomposition, mo_diag_tendency, mo_diff 
                mo_hdiff, mo_memory_sp, mo_midatm, mo_semi_impl, mo_start_dataset 
                mo_truncation 

HTML derived from FORTRAN source by f2html.pl v0.3 (C) 1997,98 Beroud Jean-Marc.