diffv_corr Subroutine

public subroutine diffv_corr()

Uses

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

Arguments

None

Calls

proc~~diffv_corr~~CallsGraph proc~diffv_corr diffv_corr zstart zstart proc~diffv_corr->zstart

Called by

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

Source Code

   subroutine diffv_corr
     ! Negate subgrid rhs contributions from solid points (added by diffv 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 : v0, vp
     use modsubgriddata, only : ekm
     use decomp_2d,      only : zstart

     real :: epmo, emmo, eomp, eomm
     integer :: i, j, k, n, m

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

         if (abs(mask_v(i+1,j,k)) < eps1) then
           epmo = 0.25 * (ekm(i,j,k) + ekm(i,j-1,k) + ekm(i+1,j-1,k) + ekm(i+1,j,k))
           vp(i,j,k) = vp(i,j,k) - epmo * (v0(i+1,j,k) - v0(i,j,k))*dx2i
         end if

         if (abs(mask_v(i-1,j,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))
           vp(i,j,k) = vp(i,j,k) + emmo * (v0(i,j,k) - v0(i-1,j,k))*dx2i
         end if

         if (abs(mask_v(i,j,k+1)) < eps1) then
           eomp = ( dzf(k+1) * ( ekm(i,j,k)   + ekm(i,j-1,k)  )  + &
                    dzf(k  ) * ( ekm(i,j,k+1) + ekm(i,j-1,k+1))) * dzhiq(k+1)
           vp(i,j,k) = vp(i,j,k) - eomp * (v0(i,j,k+1) - v0(i,j,k))*dzhi(k+1)*dzfi(k)
         end if

         if (abs(mask_v(i,j,k-1)) < eps1) then
           eomm = ( dzf(k-1) * ( ekm(i,j,k  )  + ekm(i,j-1,k)   ) + &
                    dzf(k)   * ( ekm(i,j,k-1)  + ekm(i,j-1,k-1))) * dzhiq(k)
           vp(i,j,k) = vp(i,j,k) + eomm * (v0(i,j,k) - v0(i,j,k-1))*dzhi(k)*dzfi(k)
         end if

     end do

   end subroutine diffv_corr