2D AMRCLAW
Functions/Subroutines
bound_orig.f90 File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine bound (time, nvar, ng, valbig, mitot, mjtot, mptr, aux, naux)
 

Function/Subroutine Documentation

◆ bound()

subroutine bound ( real(kind=8), intent(in)  time,
integer, intent(in)  nvar,
integer, intent(in)  ng,
real(kind=8), dimension(nvar,mitot,mjtot), intent(inout)  valbig,
integer, intent(in)  mitot,
integer, intent(in)  mjtot,
integer, intent(in)  mptr,
real(kind=8), dimension(naux,mitot,mjtot), intent(inout)  aux,
integer, intent(in)  naux 
)

Definition at line 20 of file bound_orig.f90.

References amr_module::cornxhi, amr_module::cornxlo, amr_module::cornyhi, amr_module::cornylo, filrecur(), amr_module::hxposs, amr_module::hyposs, amr_module::ndihi, amr_module::ndilo, amr_module::ndjhi, amr_module::ndjlo, amr_module::nestlevel, amr_module::node, prefilrecur(), amr_module::rnode, amr_module::spheredom, amr_module::xlower, amr_module::xperdom, amr_module::xupper, amr_module::ylower, amr_module::yperdom, and amr_module::yupper.

20 
24  use amr_module, only: xperdom, yperdom, spheredom
25 
26  implicit none
27 
28  ! Input
29  integer, intent(in) :: nvar, ng, mitot, mjtot, mptr, naux
30  real(kind=8), intent(in) :: time
31  real(kind=8), intent(in out) :: valbig(nvar,mitot,mjtot)
32  real(kind=8), intent(in out) :: aux(naux,mitot,mjtot)
33 
34  ! Locals
35  integer :: ilo, ihi, jlo, jhi, level,i,j ! rect(4)
36  real(kind=8) :: xleft, xright, ybot, ytop, hx, hy, xl, xr, yb, yt
37 
38  xleft = rnode(cornxlo, mptr)
39  xright = rnode(cornxhi, mptr)
40  ybot = rnode(cornylo, mptr)
41  ytop = rnode(cornyhi, mptr)
42  ilo = node(ndilo, mptr)
43  ihi = node(ndihi, mptr)
44  jlo = node(ndjlo, mptr)
45  jhi = node(ndjhi, mptr)
46  level = node(nestlevel, mptr)
47  hx = hxposs(level)
48  hy = hyposs(level)
49 
50 
51  ! left boundary
52  xl = xleft - ng*hx
53  xr = xleft
54  yb = ybot
55  yt = ytop
56 
57  !rect = [ilo-ng,ilo-1,jlo,jhi]
58  if ((xl < xlower) .and. xperdom) then
59  call prefilrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot,1,ng+1, &
60  ilo-ng,ilo-1,jlo,jhi)
61  else
62  call filrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot,1,ng+1,ilo-ng,ilo-1,jlo,jhi)
63  endif
64 
65  ! right boundary
66  xl = xright
67  xr = xright + ng*hx
68  yb = ybot
69  yt = ytop
70 
71  !rect = [ihi+1,ihi+ng,jlo,jhi]
72  if ((xr .gt. xupper) .and. xperdom) then
73  call prefilrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot, &
74  mitot-ng+1,ng+1,ihi+1,ihi+ng,jlo,jhi)
75  else
76  call filrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot, &
77  mitot-ng+1,ng+1,ihi+1,ihi+ng,jlo,jhi)
78  endif
79 
80 
81  ! bottom boundary
82  xl = xleft - ng*hx
83  xr = xright + ng*hx
84  yb = ybot - ng*hy
85  yt = ybot
86 
87  !rect = [ilo-ng,ihi+ng,jlo-ng,jlo-1]
88  if ( ((yb < ylower) .and. (yperdom .or. spheredom)) .or. &
89  ( ((xl < xlower) .or. (xr > xupper)) .and. xperdom) ) then
90  call prefilrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot,1,1,ilo-ng,ihi+ng,jlo-ng,jlo-1)
91  else
92  call filrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot,1,1,ilo-ng,ihi+ng,jlo-ng,jlo-1)
93  endif
94 
95  ! top boundary
96  xl = xleft - ng*hx
97  xr = xright + ng*hx
98  yb = ytop
99  yt = ytop + ng*hy
100 
101  !rect = [ilo-ng,ihi+ng,jhi+1,jhi+ng]
102  if ( ((yt .gt. yupper) .and. (yperdom .or. spheredom)) .or. &
103  (((xl .lt. xlower) .or. (xr .gt. xupper)) .and. xperdom) ) then
104  call prefilrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot, &
105  1,mjtot-ng+1,ilo-ng,ihi+ng,jhi+1,jhi+ng)
106  else
107  call filrecur(level,nvar,valbig,aux,naux,time,mitot,mjtot, &
108  1,mjtot-ng+1,ilo-ng,ihi+ng,jhi+1,jhi+ng)
109  endif
110 
111 ! external boundary conditions THIS IS NOW DONE IN THE FILPATCHES
112 ! (in the recursive filrecur.f, since filpatches had to call bc2amr, have them
113 ! all do it).
114 
integer, parameter cornxlo
x-coordinate of the left border of this grid
Definition: amr_module.f90:143
recursive subroutine prefilrecur(level, nvar, valbig, auxbig, naux, time, mitot, mjtot, nrowst, ncolst, ilo, ihi, jlo, jhi, iglo, ighi, jglo, jghi, patchOnly)
For periodic boundary conditions more work needed to fill the piece of the boundary.
Definition: prefilp.f90:20
real(kind=8), dimension(maxlv) hyposs
Definition: amr_module.f90:193
real(kind=8), dimension(maxlv) hxposs
Definition: amr_module.f90:193
integer, parameter ndihi
global i index of right border of this grid
Definition: amr_module.f90:111
integer, dimension(nsize, maxgr) node
Definition: amr_module.f90:198
integer, parameter nestlevel
AMR level of the grid.
Definition: amr_module.f90:44
real(kind=8) xupper
Definition: amr_module.f90:231
real(kind=8), dimension(rsize, maxgr) rnode
Definition: amr_module.f90:193
real(kind=8) xlower
Definition: amr_module.f90:231
logical yperdom
Definition: amr_module.f90:230
real(kind=8) yupper
Definition: amr_module.f90:231
real(kind=8) ylower
Definition: amr_module.f90:231
integer, parameter ndilo
global i index of left border of this grid
Definition: amr_module.f90:108
integer, parameter ndjlo
global j index of lower border of this grid
Definition: amr_module.f90:114
logical spheredom
Definition: amr_module.f90:230
recursive subroutine filrecur(level, nvar, valbig, aux, naux, t, mitot, mjtot, nrowst, ncolst, ilo, ihi, jlo, jhi, patchOnly, msrc)
Fill a region (patch) described by:
Definition: filpatch.f90:73
integer, parameter ndjhi
global j index of upper border of this grid
Definition: amr_module.f90:117
logical xperdom
Definition: amr_module.f90:230
integer, parameter cornylo
y-coordinate of the lower border of this grid
Definition: amr_module.f90:145
integer, parameter cornxhi
x-coordinate of the right border of this grid
Definition: amr_module.f90:147
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
integer, parameter cornyhi
y-coordinate of the upper border of this grid
Definition: amr_module.f90:149
Here is the call graph for this function: