tke_eps_advect_3d.F90 4.29 KB
Newer Older
hb's avatar
hb committed
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142
#include "cppdefs.h"
!-----------------------------------------------------------------------
!BOP
!
! !ROUTINE: tke_eps_advect_3d - 3D turbulence advection 
!
! !INTERFACE:
   subroutine tke_eps_advect_3d(hor_adv,ver_adv,adv_split)
!
! !DESCRIPTION:
!
! This routine carries out advection of the prognostic turbulence quantities
! {\tt tke} (turbuent kinetic energy, $k$) and {\tt eps} (lenght scale related
! turbulence quantity, e.g.\ dissipation rate of $k$, $\varepsilon$, or 
! turbulent frequency, $\omega=\varepsilon/k$. Here, the TVD advection 
! schemes are used which are also used for the momentum advection. 
!
! !
! !USES:
   use domain, only: imin,imax,jmin,jmax,kmax,az,au,av
#if defined(SPHERICAL) || defined(CURVILINEAR)
   use domain, only: dxu,dxv,dyu,dyv,arcd1
#else
   use domain, only: dx,dy,ard1
#endif
   use variables_3d, only: dt,kumin,kvmin,uu,vv,ww,hun,hvn,ho,hn,uuEx,vvEx
#ifdef UV_TVD
   use variables_3d, only: uadv,vadv,wadv,huadv,hvadv,hoadv,hnadv
#endif
   use variables_3d, only: tke,eps
   use advection_3d, only: do_advection_3d
   use halo_zones, only: update_3d_halo,wait_halo,U_TAG,V_TAG
   IMPLICIT NONE
!
! !INPUT PARAMETERS:
   integer, intent(in)  :: hor_adv,ver_adv,adv_split
!
! !INPUT/OUTPUT PARAMETERS:
!
! !OUTPUT PARAMETERS:
!
! !REVISION HISTORY:
!  Original author(s): Hans Burchard & Karsten Bolding
!
! !LOCAL VARIABLES:
   integer                   :: i,j,k
#ifdef UV_TVD
   REALTYPE                  :: dxuadv(I2DFIELD)
   REALTYPE                  :: dxvadv(I2DFIELD)
   REALTYPE                  :: dyuadv(I2DFIELD)
   REALTYPE                  :: dyvadv(I2DFIELD)
   REALTYPE                  :: area_inv(I2DFIELD)
   REALTYPE                  :: AH=_ZERO_
   REALTYPE                  :: dxdyi
#endif

!EOP
!-----------------------------------------------------------------------
!BOC
#ifdef DEBUG
   integer, save :: Ncall = 0
   Ncall = Ncall+1
   write(debug,*) 'tke_eps_advect_3d() # ',Ncall
#endif

#ifndef UV_TVD
   STDERR 'Do not use tke_eps_advect_3d() without the compiler option UV_TVD'
   stop 'tke_eps_advect_3d()'
#else


#if defined(SPHERICAL) || defined(CURVILINEAR) 
#else
   dxdyi = _ONE_/(dx*dy)
#endif

! Here begins dimensional split advection for tke and eps
   do k=1,kmax-1
     do j=jmin-HALO,jmax+HALO
         do i=imin-HALO,imax+HALO
            uadv(i,j,k)=_HALF_*(uu(i,j,k)+uu(i,j,k+1))
            vadv(i,j,k)=_HALF_*(vv(i,j,k)+vv(i,j,k+1))
            wadv(i,j,k)=_HALF_*(ww(i,j,k)+ww(i,j,k+1))
            huadv(i,j,k)=_HALF_*(hun(i,j,k)+hun(i,j,k+1))
            hvadv(i,j,k)=_HALF_*(hvn(i,j,k)+hvn(i,j,k+1))
            hoadv(i,j,k)=_HALF_*(ho(i,j,k)+ho(i,j,k+1))
            hnadv(i,j,k)=_HALF_*(hn(i,j,k)+hn(i,j,k+1))
         end do
      end do
   end do
   k=kmax  ! Only needed to allow for advection of k=kmax
   do j=jmin-HALO,jmax+HALO
      do i=imin-HALO,imax+HALO
         uadv(i,j,k)=_HALF_*(uu(i,j,k))
         vadv(i,j,k)=_HALF_*(vv(i,j,k))
         wadv(i,j,k)=_HALF_*(ww(i,j,k))
         huadv(i,j,k)=_HALF_*(hun(i,j,k))
         hvadv(i,j,k)=_HALF_*(hvn(i,j,k))
         hoadv(i,j,k)=_HALF_*(ho(i,j,k))
         hnadv(i,j,k)=_HALF_*(hn(i,j,k))
      end do
   end do

   do j=jmin-HALO,jmax+HALO
      do i=imin-HALO,imax+HALO
#if defined(SPHERICAL) || defined(CURVILINEAR) 
         dxuadv(i,j)=dxu(i,j)
         dxvadv(i,j)=dxv(i,j)
         dyuadv(i,j)=dyu(i,j)
         dyvadv(i,j)=dyv(i,j)
         area_inv(i,j)=arcd1(i,j)
#else
         dxuadv(i,j)=dx
         dxvadv(i,j)=dx
         dyuadv(i,j)=dy
         dyvadv(i,j)=dy
         area_inv(i,j)=dxdyi
#endif
      end do
   end do

   call do_advection_3d(dt,tke,uadv,vadv,wadv,huadv,hvadv,hoadv,hnadv, &
                        dxuadv,dxvadv,dyuadv,dyvadv,area_inv,          &
                        az,au,av,hor_adv,ver_adv,adv_split,AH)

   call do_advection_3d(dt,eps,uadv,vadv,wadv,huadv,hvadv,hoadv,hnadv, &
                        dxuadv,dxvadv,dyuadv,dyvadv,area_inv,          &
                        az,au,av,hor_adv,ver_adv,adv_split,AH)

#endif

#ifdef DEBUG
   write(debug,*) 'Leaving tke_eps_advect_3d()'
   write(debug,*)
#endif
   return
   end subroutine tke_eps_advect_3d
!EOC

!-----------------------------------------------------------------------
! Copyright (C) 2010 - Hans Burchard and Karsten Bolding               !
!-----------------------------------------------------------------------