subroutine solid(solid_info, var, rhs, val, hi, hj, hk, mask)
use modglobal, only : ib, ie, jb, je, kb, ke, ih, jh, kh, eps1
use decomp_2d, only : zstart
type(solid_info_type), intent(in) :: solid_info
integer, intent(in) :: hi, hj, hk
real, intent(inout) :: var(ib-hi:ie+hi,jb-hj:je+hj,kb-hk:ke+hk)
real, intent(inout) :: rhs(ib-hi:ie+hi,jb-hj:je+hj,kb :ke+hk)
real, intent(in) :: val
real, intent(in), optional :: mask(ib-ih:ie+ih,jb-jh:je+jh,kb-kh:ke+kh)
real :: count
integer :: i, j, k, n, m
if (present(mask) .eqv. .false.) then
do n=1,solid_info%nsolptsrank
!n = solid_info%solptsrank(m)
i = solid_info%solpts_loc(n,1) - zstart(1) + 1
j = solid_info%solpts_loc(n,2) - zstart(2) + 1
k = solid_info%solpts_loc(n,3) - zstart(3) + 1
var(i,j,k) = val
rhs(i,j,k) = 0.
end do
else
do n=1,solid_info%nsolptsrank
!n = solid_info%solptsrank(m)
i = solid_info%solpts_loc(n,1) - zstart(1) + 1
j = solid_info%solpts_loc(n,2) - zstart(2) + 1
k = solid_info%solpts_loc(n,3) - zstart(3) + 1
var(i,j,k) = val
rhs(i,j,k) = 0.
count = 0
! Attempt to set zero flux BC
if (abs(mask(i,j+1,k) - 1.) < eps1) then ! fluid neighbour
count = count + 1
var(i,j,k) = var(i,j,k) + var(i,j+1,k)
rhs(i,j,k) = rhs(i,j,k) + rhs(i,j+1,k)
end if
if (abs(mask(i,j-1,k) - 1.) < eps1) then
count = count + 1
var(i,j,k) = var(i,j,k) + var(i,j-1,k)
rhs(i,j,k) = rhs(i,j,k) + rhs(i,j-1,k)
end if
if (abs(mask(i,j,k+1) - 1.) < eps1) then
count = count + 1
var(i,j,k) = var(i,j,k) + var(i,j,k+1)
rhs(i,j,k) = rhs(i,j,k) + rhs(i,j,k+1)
end if
if (abs(mask(i,j,k-1) - 1.) < eps1) then
count = count + 1
var(i,j,k) = var(i,j,k) + var(i,j,k-1)
rhs(i,j,k) = rhs(i,j,k) + rhs(i,j,k-1)
end if
if (abs(mask(i+1,j,k) - 1.) < eps1) then
count = count + 1
var(i,j,k) = var(i,j,k) + var(i+1,j,k)
rhs(i,j,k) = rhs(i,j,k) + rhs(i+1,j,k)
end if
if (abs(mask(i-1,j,k) - 1.) < eps1) then
count = count + 1
var(i,j,k) = var(i,j,k) + var(i-1,j,k)
rhs(i,j,k) = rhs(i,j,k) + rhs(i-1,j,k)
end if
if (count > 0) then
var(i,j,k) = (var(i,j,k) - val) / count
rhs(i,j,k) = rhs(i,j,k) / count
end if
end do
end if
end subroutine solid