diffu_corr Subroutine

public subroutine diffu_corr()

Uses

  • proc~~diffu_corr~~UsesGraph proc~diffu_corr diffu_corr decomp_2d decomp_2d proc~diffu_corr->decomp_2d module~modfields modfields proc~diffu_corr->module~modfields module~modglobal modglobal proc~diffu_corr->module~modglobal module~modsubgriddata modsubgriddata proc~diffu_corr->module~modsubgriddata module~modfields->decomp_2d

Arguments

None

Calls

proc~~diffu_corr~~CallsGraph proc~diffu_corr diffu_corr zstart zstart proc~diffu_corr->zstart

Called by

proc~~diffu_corr~~CalledByGraph proc~diffu_corr diffu_corr proc~ibmwallfun ibmwallfun proc~ibmwallfun->proc~diffu_corr program~dalesurban DALESURBAN program~dalesurban->proc~ibmwallfun

Source Code

   subroutine diffu_corr
     ! Negate subgrid rhs contributions from solid points (added by diffu in modsubgrid)
     use modglobal,      only : eps1, ib, ie, ih, jb, je, jh, kb, ke, kh, &
                                dx2i, dxi5, dy2i, dyi5, dzf, dzh2i, dzfi, dzhi, dzfi5, dzhiq
     use modfields,      only : u0, up
     use modsubgriddata, only : ekm
     use decomp_2d,      only : zstart

     real :: empo, emmo, emop, emom
     integer :: i, j, k, n, m

     do n = 1,bound_info_u%nbndptsrank
      !n = bound_info_u%bndptsrank(m)
         i = bound_info_u%bndpts_loc(n,1) - zstart(1) + 1
         j = bound_info_u%bndpts_loc(n,2) - zstart(2) + 1
         k = bound_info_u%bndpts_loc(n,3) - zstart(3) + 1

         if (abs(mask_u(i,j+1,k)) < eps1) then
           empo = 0.25 * ((ekm(i,j,k) + ekm(i,j+1,k)) + (ekm(i-1,j,k) + ekm(i-1,j+1,k)))
           up(i,j,k) = up(i,j,k) - empo * (u0(i,j+1,k) - u0(i,j,k))*dy2i
         end if

         if (abs(mask_u(i,j-1,k)) < eps1) then
           emmo = 0.25 * ((ekm(i,j,k) + ekm(i,j-1,k)) + (ekm(i-1,j-1,k) + ekm(i-1,j,k)))
           up(i,j,k) = up(i,j,k) + emmo * (u0(i,j,k) - u0(i,j-1,k))*dy2i
         end if

         if (abs(mask_u(i,j,k+1)) < eps1) then
           emop = (dzf(k+1) * ( ekm(i,j,k)   + ekm(i-1,j,k  ))  + &
                   dzf(k)   * ( ekm(i,j,k+1) + ekm(i-1,j,k+1))) * dzhiq(k+1)
           up(i,j,k) = up(i,j,k) - emop * (u0(i,j,k+1) - u0(i,j,k))*dzhi(k+1)*dzfi(k)
         end if

         if (abs(mask_u(i,j,k-1)) < eps1) then
           emom = (dzf(k-1) * (ekm(i,j,k  ) + ekm(i-1,j,k  ))  + &
                   dzf(k)   * (ekm(i,j,k-1) + ekm(i-1,j,k-1))) * dzhiq(k)
           up(i,j,k) = up(i,j,k) + emom * (u0(i,j,k) - u0(i,j,k-1))*dzhi(k)*dzfi(k)
         end if

     end do


   end subroutine diffu_corr