subroutine diffv (putout)
use modglobal, only : ib,ie,ih,jb,je,jh,kb,ke,kh,dxi,dzf,dzfi,dyi,&
dy2i,dzhi,dzhiq,jmax,numol,lles
use modfields, only : u0,v0,w0
use modsurfdata,only : ustar
use modmpi, only : myid
implicit none
real, intent(inout) :: putout(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh)
real :: emmo, eomm,eomp,epmo
real :: fv, vcv,vpcv
integer :: i,j,k,jm,jp,km,kp
if (lles) then
do k=kb,ke
kp=k+1
km=k-1
do j=jb,je
jp=j+1
jm=j-1
do i=ib,ie
eomm = ( dzf(km) * ( ekm(i,j,k) + ekm(i,jm,k) ) + &
dzf(k) * ( ekm(i,j,km) + ekm(i,jm,km) ) ) * dzhiq(k)
eomp = ( dzf(kp) * ( ekm(i,j,k) + ekm(i,jm,k) ) + &
dzf(k) * ( ekm(i,j,kp) + ekm(i,jm,kp) ) ) * dzhiq(kp)
emmo = 0.25 * ( ekm(i,j,k)+ekm(i,jm,k)+ekm(i-1,jm,k)+ekm(i-1,j,k) )
epmo = 0.25 * ( ekm(i,j,k)+ekm(i,jm,k)+ekm(i+1,jm,k)+ekm(i+1,j,k) )
! discretized diffusion term
putout(i,j,k) = putout(i,j,k) &
+ ( &
epmo * ( &
(v0(i+1,j,k)-v0(i,j,k)) *dxi &
+ (u0(i+1,j,k)-u0(i+1,jm,k))*dyi &
) &
-emmo * ( &
(v0(i,j,k)-v0(i-1,j,k)) *dxi &
+ (u0(i,j,k)-u0(i,jm,k)) *dyi &
) &
) * dxi & ! = d/dx( Km*(dv/dx + du/dy) )
+ ( &
ekm(i,j,k) * (v0(i,jp,k)-v0(i,j,k)) &
- ekm(i,jm,k)* (v0(i,j,k)-v0(i,jm,k)) &
) * 2. * dy2i & ! = d/dy( 2*Km*(dv/dy) )
+ ( &
eomp * ( &
(v0(i,j,kp)-v0(i,j,k)) *dzhi(kp) &
+ (w0(i,j,kp)-w0(i,jm,kp)) *dyi &
) &
- eomm * ( &
(v0(i,j,k)-v0(i,j,km)) *dzhi(k) &
+ (w0(i,j,k)-w0(i,jm,k)) *dyi &
) &
) * dzfi(k) ! = d/dz( Km*(dv/dz + dw/dy) )
end do
end do
end do
else ! DNS
do k=kb,ke
kp=k+1
km=k-1
do j=jb,je
jp=j+1
jm=j-1
do i=ib,ie
putout(i,j,k) = putout(i,j,k) &
+ ( &
numol * ( &
(v0(i+1,j,k)-v0(i,j,k)) *dxi &
+ (u0(i+1,j,k)-u0(i+1,jm,k))*dyi &
) &
- numol * ( &
(v0(i,j,k)-v0(i-1,j,k)) *dxi &
+ (u0(i,j,k)-u0(i,jm,k)) *dyi &
) &
) * dxi & ! = d/dx( Km*(dv/dx + du/dy) )
+ ( &
numol * (v0(i,jp,k)-v0(i,j,k)) &
- numol * (v0(i,j,k)-v0(i,jm,k)) &
) * 2. * dy2i & ! = d/dy( 2*Km*(dv/dy) )
+ ( &
numol * ( &
(v0(i,j,kp)-v0(i,j,k)) *dzhi(kp) &
+ (w0(i,j,kp)-w0(i,jm,kp)) *dyi) &
- numol * ( &
(v0(i,j,k)-v0(i,j,km)) *dzhi(k) &
+ (w0(i,j,k)-w0(i,jm,k)) *dyi &
) &
) * dzfi(k) ! = d/dz( Km*(dv/dz + dw/dy) )
end do
end do
end do
end if
end subroutine diffv