!-----------------------------------------------------------------------
! Physical boundary conditions
! Interface:
! mx,my := shape of grid function
!
! u(,) := grid function
!
! lb(1) := lower bound for grid
! ub(1) := upper bound for grid
! lbbnd(1) := lower bound for boundary region
! ubnd(1) := upper bound for boundary region
! shapebnd(1) := shape of boundary region
! xc(1) := lower left corner of grid
! dx(1) := grid spacing
! dir := at which side of the grid is the boundary?
! bnd(,2,1) := lower left and upper right corner of global grid and
! of mb-1 internal boundary regions
!-----------------------------------------------------------------------
subroutine physbd(u,mx,lb,ub,lbbnd,ubbnd,shapebnd,
& xc,dx,dir,bnd,mb,time,meqn)
implicit double precision (a-h,o-z)
c
include "cuser.i"
c
integer meqn, mx, mb, dir
integer lb(1), ub(1), lbbnd(1), ubbnd(1), shapebnd(1)
double precision u(meqn,mx), xc(1), dx(1), bnd(mb,2,1), time
! Local variables
integer i, imin, imax, m, isym
integer stride, getindx
c
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! See definition of member-function extents() in BBox.h
! for calculation of stride
stride = (ub(1) - lb(1))/(mx-1)
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
! Find working domain
imin = getindx(max(lbbnd(1), lb(1)), lb(1), stride)
imax = getindx(min(ubbnd(1), ub(1)), lb(1), stride)
if(imax .gt. mx .or. imin .lt. 1) then
write(0,*)'INDEX ERROR in physbd'
end if
go to (100,200) dir+1
! Left Side --- Fixed wall
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
100 continue
do 110 i = imax, imin, -1
isym = 2*imax+1-i
do m = 1, meqn
u(m,i) = u(m,isym)
enddo
u(3,i) = -u(3,isym)
110 continue
return
! Right Side --- Symmetry or Fixed Walls
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
200 continue
if (Nbnd.eq.1) then
do 210 i = imin, imax
isym = 2*imin-1-i
do m = 1, meqn
u(m,i) = u(m,isym)
enddo
u(3,i) = -u(3,isym)
210 continue
else
! Right Side --- Outflow
!- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
do 310 i = imin, imax
do m = 1, meqn
u(m,i) = u(m,i-1)
enddo
310 continue
endif
return
c
end