diffw_corr Subroutine

public subroutine diffw_corr()

Uses

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

Arguments

None

Calls

proc~~diffw_corr~~CallsGraph proc~diffw_corr diffw_corr zstart zstart proc~diffw_corr->zstart

Called by

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

Source Code

   subroutine diffw_corr
     ! Negate subgrid rhs contributions from solid points (added by diffw 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 : w0, wp
     use modsubgriddata, only : ekm
     use decomp_2d,      only : zstart

     real :: epom, emom, eopm, eomm
     integer :: i, j, k, n, m

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

         ! Account for solid w points
         if (abs(mask_w(i+1,j,k)) < eps1) then
           epom = ( 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)
           wp(i,j,k) = wp(i,j,k) - epom * (w0(i+1,j,k) - w0(i,j,k))*dx2i
         end if

         if (abs(mask_w(i-1,j,k)) < 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)
           wp(i,j,k) = wp(i,j,k) + emom * (w0(i,j,k) - w0(i-1,j,k))*dx2i
         end if

         if (abs(mask_w(i,j+1,k)) < eps1) then
           eopm = ( 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)
           wp(i,j,k) = wp(i,j,k) - eopm * (w0(i,j+1,k) - w0(i,j,k))*dy2i
         end if

         if (abs(mask_w(i,j-1,k)) < 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)
           wp(i,j,k) = wp(i,j,k) + eomm * (w0(i,j,k) - w0(i,j-1,k))*dy2i
         end if

     end do

   end subroutine diffw_corr