c-----------------------------------------------------------------------
c Three-dimensional restriction operator for AMROC.
c A coarse cell value is overwritten by the mean value
c of all refined cells within this particular coarse cell.
c
c Interface:
c mfx,mfy,mfz := shape of fine grid
c mcx,mcy,mcz := shape of coarse grid
c
c uf(,,) := fine grid
c uc(,,) := coarse grid
c
c lbc(3) := lower bound for coarse grid
c ubc(3) := upper bound for coarse grid
c lbf(3) := lower bound for fine grid
c ubf(3) := upper bound for fine grid
c lbr(3) := lower bound for region restriction desired
c ufr(3) := upper bound for region restriction desired
c shaper(3) := shape of region restriction desired
c
c Copyright (C) 2002 Ralf Deiterding
c Brandenburgische Universitaet Cottbus
c
c Copyright (C) 2003-2007 California Institute of Technology
c Ralf Deiterding, ralf@cacr.caltech.edu
c-----------------------------------------------------------------------
subroutine restrict3(uf,mfx,mfy,mfz,lbf,ubf,
& uc,mcx,mcy,mcz,lbc,ubc,
& lbr,ubr,shaper,meqn,mbc)
implicit none
integer meqn,mbc,mfx,mfy,mfz,mcx,mcy,mcz
integer shaper(3)
double precision uf(meqn,mfx,mfy,mfz), uc(meqn,mcx,mcy,mcz)
integer lbf(3), ubf(3),
& lbc(3), ubc(3),
& lbr(3), ubr(3)
c Local variables
integer i, j, k, m, imin, imax, jmin, jmax, kmin, kmax,
& ii, jj, kk, ifine, icoarse, jfine, jcoarse, kfine, kcoarse,
& refine, stridec, stridef, getindx, mbcf
c- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
c See definition of member-function extents() in BBox.h
c for calculation of stride
c- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
stridec = (ubc(1) - lbc(1))/(mcx-1)
stridef = (ubf(1) - lbf(1))/(mfx-1)
refine = stridec/stridef
mbcf = mbc * stridef
c- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
c Find coarse domain over which to refine
c Take three regions and select out intersection
c- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
imin = max(lbf(1)+mbcf, lbc(1), lbr(1))
imax = min(ubf(1)-mbcf, ubc(1), ubr(1))
jmin = max(lbf(2)+mbcf, lbc(2), lbr(2))
jmax = min(ubf(2)-mbcf, ubc(2), ubr(2))
kmin = max(lbf(3)+mbcf, lbc(3), lbr(3))
kmax = min(ubf(3)-mbcf, ubc(3), ubr(3))
if (mod(imin-lbc(1),stridec) .ne. 0) then
imin = imin + stridec - mod(imin-lbc(1),stridec)
endif
if (mod(jmin-lbc(2),stridec) .ne. 0) then
jmin = jmin + stridec - mod(jmin-lbc(2),stridec)
endif
if (mod(kmin-lbc(3),stridec) .ne. 0) then
kmin = kmin + stridec - mod(kmin-lbc(3),stridec)
endif
c- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
c Inject points to coarse grid from fine grid
c Loop from lower bound to upper bound with stride of refine.
c Convert the integer coordinates to fine and coarse grid absolute
c coordinates...
c- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
do 10 k=kmin, kmax, stridec
kfine = getindx(k, lbf(3), stridef)
kcoarse = getindx(k, lbc(3), stridec)
do 10 j=jmin, jmax, stridec
jfine = getindx(j, lbf(2), stridef)
jcoarse = getindx(j, lbc(2), stridec)
do 10 i=imin, imax, stridec
ifine = getindx(i, lbf(1), stridef)
icoarse = getindx(i, lbc(1), stridec)
! if (icoarse .gt. mcx .or.
! & jcoarse .gt. mcy .or.
! & kcoarse .gt. mcz) then
! write(0,*)'ERROR in restriction: ',
! & icoarse,jcoarse,kcoarse
! end if
do 10 m=1, meqn
uc(m,icoarse,jcoarse,kcoarse) = 0
do 20 kk=0, refine-1
do 20 jj=0, refine-1
do 20 ii=0, refine-1
uc(m,icoarse,jcoarse,kcoarse) =
& uc(m,icoarse,jcoarse,kcoarse) +
& uf(m,ifine+ii,jfine+jj,kfine+kk)
20 continue
uc(m,icoarse,jcoarse,kcoarse) =
& uc(m,icoarse,jcoarse,kcoarse) / refine**3
10 continue
return
end