BOUNDARIES: Rewrite block_corner_prolong().
Signed-off-by: Grzegorz Kowal <grzegorz@amuncode.org>
This commit is contained in:
parent
115659486b
commit
9b030f0b13
@ -946,7 +946,7 @@ module boundaries
|
||||
|
||||
! prolong corner boundaries from blocks at lower levels
|
||||
!
|
||||
call boundaries_corner_prolong()
|
||||
call boundaries_corner_prolong(status)
|
||||
|
||||
end if ! minlev /= maxlev
|
||||
|
||||
@ -4388,7 +4388,7 @@ module boundaries
|
||||
!
|
||||
!===============================================================================
|
||||
!
|
||||
subroutine boundaries_corner_prolong()
|
||||
subroutine boundaries_corner_prolong(status)
|
||||
|
||||
! import external procedures and variables
|
||||
!
|
||||
@ -4410,6 +4410,8 @@ module boundaries
|
||||
!
|
||||
implicit none
|
||||
|
||||
integer, intent(out) :: status
|
||||
|
||||
! local pointers
|
||||
!
|
||||
type(block_meta), pointer :: pmeta, pneigh
|
||||
@ -4529,12 +4531,12 @@ module boundaries
|
||||
#if NDIMS == 2
|
||||
call block_corner_prolong((/ i, j, k /) &
|
||||
, pneigh%data%q(1:nv, : , : , : ) &
|
||||
, pmeta%data%q(1:nv,il:iu,jl:ju, : ))
|
||||
, pmeta%data%q(1:nv,il:iu,jl:ju, : ), status)
|
||||
#endif /* NDIMS == 2 */
|
||||
#if NDIMS == 3
|
||||
call block_corner_prolong((/ i, j, k /) &
|
||||
, pneigh%data%q(1:nv, : , : , : ) &
|
||||
, pmeta%data%q(1:nv,il:iu,jl:ju,kl:ku))
|
||||
, pmeta%data%q(1:nv,il:iu,jl:ju,kl:ku), status)
|
||||
#endif /* NDIMS == 3 */
|
||||
|
||||
#ifdef MPI
|
||||
@ -4644,12 +4646,12 @@ module boundaries
|
||||
#if NDIMS == 2
|
||||
call block_corner_prolong((/ i, j, k /) &
|
||||
, pneigh%data%q(1:nv, : , : , : ) &
|
||||
, sbuf(l,1:nv,1:ng,1:ng, : ))
|
||||
, sbuf(l,1:nv,1:ng,1:ng, : ), status)
|
||||
#endif /* NDIMS == 2 */
|
||||
#if NDIMS == 3
|
||||
call block_corner_prolong((/ i, j, k /) &
|
||||
, pneigh%data%q(1:nv, : , : , : ) &
|
||||
, sbuf(l,1:nv,1:ng,1:ng,1:ng))
|
||||
, sbuf(l,1:nv,1:ng,1:ng,1:ng), status)
|
||||
#endif /* NDIMS == 3 */
|
||||
|
||||
! associate the pointer with the next block
|
||||
@ -5996,38 +5998,45 @@ module boundaries
|
||||
! subroutine BLOCK_CORNER_PROLONG:
|
||||
! -------------------------------
|
||||
!
|
||||
! Subroutine returns the corner boundary region by prolongating
|
||||
! the corresponding region from the provided input variable array.
|
||||
! Subroutine prolongs the region of qn specified by the corner position
|
||||
! and inserts it in the corresponding corner boundary region of qb.
|
||||
!
|
||||
! Remarks:
|
||||
!
|
||||
! This subroutine requires the ghost zone regions of the qn already
|
||||
! updated in order to perform a consistent interpolation.
|
||||
!
|
||||
! Arguments:
|
||||
!
|
||||
! pos - the corner position;
|
||||
! qn - the input neighbor variable array;
|
||||
! qb - the output corner boundary array;
|
||||
! pos - the corner position;
|
||||
! qn - the array to prolong from;
|
||||
! qb - the array of the corner boundary region;
|
||||
! status - the call status;
|
||||
!
|
||||
!===============================================================================
|
||||
!
|
||||
subroutine block_corner_prolong(pos, qn, qb)
|
||||
subroutine block_corner_prolong(pos, qn, qb, status)
|
||||
|
||||
use coordinates , only : corners_dp
|
||||
use coordinates , only : nb, ne, nh => nghosts_half
|
||||
use equations , only : nv, positive
|
||||
use helpers , only : print_message
|
||||
use interpolations, only : limiter_prol
|
||||
|
||||
implicit none
|
||||
|
||||
integer , dimension(3) , intent(in) :: pos
|
||||
real(kind=8), dimension(:,:,:,:), intent(in) :: qn
|
||||
real(kind=8), dimension(:,:,:,:), intent(out) :: qb
|
||||
integer , dimension(3) , intent(in) :: pos
|
||||
real(kind=8), dimension(:,:,:,:), intent(in) :: qn
|
||||
real(kind=8), dimension(:,:,:,:), intent(inout) :: qb
|
||||
integer , intent(out) :: status
|
||||
|
||||
integer :: p
|
||||
integer :: i, il, iu, is, it, im1, ip1
|
||||
integer :: j, jl, ju, js, jt, jm1, jp1
|
||||
integer :: k, kt
|
||||
integer :: i, il, iu, im, ip, is, it
|
||||
integer :: j, jl, ju, jm, jp, js, jt
|
||||
#if NDIMS == 3
|
||||
integer :: kl, ku, ks, km1, kp1
|
||||
integer :: k, kl, ku, km, kp, ks, kt
|
||||
#else /* NDIMS == 3 */
|
||||
integer :: k
|
||||
#endif /* NDIMS == 3 */
|
||||
real(kind=8) :: dql, dqr
|
||||
real(kind=8) :: dq1, dq2
|
||||
#if NDIMS == 3
|
||||
real(kind=8) :: dq3, dq4
|
||||
@ -6035,75 +6044,74 @@ module boundaries
|
||||
|
||||
real(kind=8), dimension(NDIMS) :: dq
|
||||
|
||||
character(len=80) :: msg
|
||||
logical , save :: first = .true.
|
||||
integer, dimension(2,2), save :: nlims
|
||||
|
||||
character(len=*), parameter :: loc = 'BOUNDARIES::block_corner_prolong()'
|
||||
!
|
||||
|
||||
!-------------------------------------------------------------------------------
|
||||
!
|
||||
#if NDIMS == 2
|
||||
il = corners_dp(pos(1),pos(2))%l(1)
|
||||
jl = corners_dp(pos(1),pos(2))%l(2)
|
||||
iu = corners_dp(pos(1),pos(2))%u(1)
|
||||
ju = corners_dp(pos(1),pos(2))%u(2)
|
||||
#endif /* NDIMS == 2 */
|
||||
#if NDIMS == 3
|
||||
il = corners_dp(pos(1),pos(2),pos(3))%l(1)
|
||||
jl = corners_dp(pos(1),pos(2),pos(3))%l(2)
|
||||
kl = corners_dp(pos(1),pos(2),pos(3))%l(3)
|
||||
iu = corners_dp(pos(1),pos(2),pos(3))%u(1)
|
||||
ju = corners_dp(pos(1),pos(2),pos(3))%u(2)
|
||||
ku = corners_dp(pos(1),pos(2),pos(3))%u(3)
|
||||
#endif /* NDIMS == 3 */
|
||||
if (first) then
|
||||
nlims(1,1) = ne - nh + 1
|
||||
nlims(1,2) = ne
|
||||
nlims(2,1) = nb
|
||||
nlims(2,2) = nb + nh - 1
|
||||
|
||||
#if NDIMS == 2
|
||||
k = 1
|
||||
kt = 1
|
||||
#endif /* NDIMS == 2 */
|
||||
first = .false.
|
||||
end if
|
||||
|
||||
il = nlims(pos(1),1)
|
||||
iu = nlims(pos(1),2)
|
||||
jl = nlims(pos(2),1)
|
||||
ju = nlims(pos(2),2)
|
||||
#if NDIMS == 3
|
||||
kl = nlims(pos(3),1)
|
||||
ku = nlims(pos(3),2)
|
||||
|
||||
do k = kl, ku
|
||||
km1 = k - 1
|
||||
kp1 = k + 1
|
||||
ks = 2 * (k - kl) + 1
|
||||
kt = ks + 1
|
||||
km = k - 1
|
||||
kp = k + 1
|
||||
ks = 2 * (k - kl) + 1
|
||||
kt = ks + 1
|
||||
#else /* NDIMS == 3 */
|
||||
k = 1
|
||||
#endif /* NDIMS == 3 */
|
||||
do j = jl, ju
|
||||
jm1 = j - 1
|
||||
jp1 = j + 1
|
||||
js = 2 * (j - jl) + 1
|
||||
jt = js + 1
|
||||
jm = j - 1
|
||||
jp = j + 1
|
||||
js = 2 * (j - jl) + 1
|
||||
jt = js + 1
|
||||
do i = il, iu
|
||||
im1 = i - 1
|
||||
ip1 = i + 1
|
||||
is = 2 * (i - il) + 1
|
||||
it = is + 1
|
||||
im = i - 1
|
||||
ip = i + 1
|
||||
is = 2 * (i - il) + 1
|
||||
it = is + 1
|
||||
|
||||
do p = 1, nv
|
||||
|
||||
dql = qn(p,i ,j,k) - qn(p,im1,j,k)
|
||||
dqr = qn(p,ip1,j,k) - qn(p,i ,j,k)
|
||||
dq(1) = limiter_prol(2.5d-01, dql, dqr)
|
||||
dq1 = qn(p,i ,j,k) - qn(p,im,j,k)
|
||||
dq2 = qn(p,ip,j,k) - qn(p,i ,j,k)
|
||||
dq(1) = limiter_prol(2.5d-01, dq1, dq2)
|
||||
|
||||
dql = qn(p,i,j ,k) - qn(p,i,jm1,k)
|
||||
dqr = qn(p,i,jp1,k) - qn(p,i,j ,k)
|
||||
dq(2) = limiter_prol(2.5d-01, dql, dqr)
|
||||
dq1 = qn(p,i,j ,k) - qn(p,i,jm,k)
|
||||
dq2 = qn(p,i,jp,k) - qn(p,i,j ,k)
|
||||
dq(2) = limiter_prol(2.5d-01, dq1, dq2)
|
||||
|
||||
#if NDIMS == 3
|
||||
dql = qn(p,i,j,k ) - qn(p,i,j,km1)
|
||||
dqr = qn(p,i,j,kp1) - qn(p,i,j,k )
|
||||
dq(3) = limiter_prol(2.5d-01, dql, dqr)
|
||||
dq1 = qn(p,i,j,k ) - qn(p,i,j,km)
|
||||
dq2 = qn(p,i,j,kp) - qn(p,i,j,k )
|
||||
dq(3) = limiter_prol(2.5d-01, dq1, dq2)
|
||||
#endif /* NDIMS == 3 */
|
||||
|
||||
if (positive(p) .and. qn(p,i,j,k) < sum(abs(dq(1:NDIMS)))) then
|
||||
if (positive(p)) then
|
||||
if (qn(p,i,j,k) > 0.0d+00) then
|
||||
do while (qn(p,i,j,k) <= sum(abs(dq(1:NDIMS))))
|
||||
dq(:) = 0.5d+00 * dq(:)
|
||||
end do
|
||||
dq1 = sum(abs(dq(1:NDIMS)))
|
||||
if (qn(p,i,j,k) <= dq1) &
|
||||
dq(:) = (5.0d-01 * qn(p,i,j,k) / dq1) * dq(:)
|
||||
else
|
||||
write(msg,"(a,3i4,a)") &
|
||||
"Positive variable is not positive at (", i, j, k, " )"
|
||||
call print_message(loc, msg)
|
||||
dq(:) = 0.0d+00
|
||||
call print_message(loc, "Positive variable is not positive!")
|
||||
status = 1
|
||||
return
|
||||
end if
|
||||
end if
|
||||
|
||||
@ -6132,12 +6140,12 @@ module boundaries
|
||||
qb(p,it,jt,kt) = qn(p,i,j,k) + dq1
|
||||
#endif /* NDIMS == 3 */
|
||||
|
||||
end do ! p
|
||||
end do
|
||||
|
||||
end do ! i
|
||||
end do ! j
|
||||
end do
|
||||
end do
|
||||
#if NDIMS == 3
|
||||
end do ! k
|
||||
end do
|
||||
#endif /* NDIMS == 3 */
|
||||
|
||||
!-------------------------------------------------------------------------------
|
||||
|
||||
Loading…
x
Reference in New Issue
Block a user