subroutine xwallscalar(hi, hj, hk, putin, putout, bcvaluem, bcvaluep, n)
use modglobal, only:jmax, dxf, dxfi, dxfi5, dxhi, dxh2i, nsv, ib, ie, jb, je, kb, ke, prandtlmoli, numol
use modfields, only:u0
use modmpi, only:myid
use modsubgriddata, only:ekh
use initfac, only:block
integer i, j, k, jl, ju, kl, ku, iww, iee
integer, intent(in) :: hi !<size of halo in i
integer, intent(in) :: hj !<size of halo in j
integer, intent(in) :: hk !<size of halo in k
real, intent(in) :: putin(ib - hi:ie + hi, jb - hj:je + hj, kb - hk:ke + hk)
real, intent(inout) :: putout(ib - hi:ie + hi, jb - hj:je + hj, kb:ke + hk)
real, intent(in) :: bcvaluem, bcvaluep
integer, intent(in) :: n
iww = block(n, 1) - 1
iee = block(n, 2) + 1
jl = MAX(block(n, 3) - myid*jmax, 1) ! starting j-index
ju = MIN(block(n, 4) - myid*jmax, jmax) ! ending j-index
kl = block(n, 5) ! starting k-index
ku = block(n, 6) ! ending k-index
!fixed flux
!remove standard diffusion term, add flux=bcvalue
do k = kl, ku
do j = jl, ju
putout(iee, j, k) = putout(iee, j, k) + ( &
0.5*(ekh(iee, j, k)*dxf(iee - 1) + ekh(iee - 1, j, k)*dxf(iee))* &
(putin(iee, j, k) - putin(iee - 1, j, k))*dxh2i(iee) - &
bcvaluep)*dxfi(iee) !
putout(iww, j, k) = putout(iww, j, k) + ( &
-0.5*(ekh(iww + 1, j, k)*dxf(iww) + ekh(iww, j, k)*dxf(iww + 1))* &
(putin(iww + 1, j, k) - putin(iww, j, k))*dxh2i(iww + 1) - &
bcvaluem)*dxfi(iww) !
end do
end do
end subroutine xwallscalar