2D AMRCLAW
Functions/Subroutines
preicall.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine preicall (val, aux, nrow, ncol, nvar, naux, ilo, ihi, jlo, jhi, level, fliparray)
 

Function/Subroutine Documentation

◆ preicall()

subroutine preicall ( dimension(nvar,nrow,ncol)  val,
dimension(naux,nrow,ncol)  aux,
  nrow,
  ncol,
  nvar,
  naux,
  ilo,
  ihi,
  jlo,
  jhi,
  level,
dimension((nrow+ncol)*nghost*(nvar+naux))  fliparray 
)

Definition at line 6 of file preicall.f.

References amr_module::hxposs, amr_module::hyposs, iadd(), icall(), amr_module::iregsz, amr_module::jregsz, amr_module::needs_to_be_set, amr_module::nghost, setaux(), amr_module::spheredom, amr_module::xlower, amr_module::xperdom, amr_module::ylower, and amr_module::yperdom.

Referenced by filval(), and saveqc().

6 c
7  use amr_module
8  implicit double precision (a-h,o-z)
9 
10  dimension fliparray((nrow+ncol)*nghost*(nvar+naux))
11  dimension val(nvar,nrow,ncol)
12  dimension aux(naux,nrow,ncol)
13 
14  dimension ist(3), iend(3), jst(3), jend(3), ishift(3), jshift(3)
15  logical xint, yint
16 
17 c
18 c NEW INDEXING - ORDER SWITCHED
19  iadd(ivar,i,j) = locflip + ivar-1 + nvar*((j-1)*nc+i-1)
20  iaddaux(iaux,i,j) = locflipaux + iaux-1 + naux*((j-1)*nc+i-1)
21 
22 c
23 c :::::::::::::: PREICALL :::::::::::::::::::::::::::::::::::::::::::
24 c For periodic boundary conditions more work needed to initialize a
25 c new grid that sticks out. This routine was
26 c called because the patch sticks out of the domain,
27 c and has periodic bc.s preprocess the patch before calling
28 c icall to shift the patch periodically back into the domain.
29 c
30 c Inputs to this routine:
31 c ilo,ihi,jlo,jhi = the location in index space of this patch.
32 c
33 c Outputs from this routine:
34 c The values of the grid are inserted
35 c directly into the enlarged val array for this piece.
36 c
37 c :::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
38 c
39 ! ## fliparray now passed in; index into it below
40  locflip = 1
41  locflipaux = 1 + nvar*(ncol+nrow)
42 c
43 c ## will divide patch into 9 possibilities (some empty):
44 c ## x sticks out left, x interior, x sticks out right
45 c ## same for y. for example, the max. would be
46 c ## i from (ilo,-1), (0,iregsz(level)-1), (iregsz(level),ihi)
47 
48  if (xperdom) then
49  ist(1) = ilo
50  ist(2) = 0
51  ist(3) = iregsz(level)
52  iend(1) = -1
53  iend(2) = iregsz(level)-1
54  iend(3) = ihi
55  ishift(1) = iregsz(level)
56  ishift(2) = 0
57  ishift(3) = -iregsz(level)
58  else ! if not periodic, set vals to only have nonnull intersection for interior regoin
59  ist(1) = iregsz(level)
60  ist(2) = ilo
61  ist(3) = iregsz(level)
62  iend(1) = -nghost
63  iend(2) = ihi
64  iend(3) = -nghost
65  ishift(1) = 0
66  ishift(2) = 0
67  ishift(3) = 0
68  endif
69 
70 
71  if (yperdom .or. spheredom) then
72  jst(1) = jlo
73  jst(2) = 0
74  jst(3) = jregsz(level)
75  jend(1) = -1
76  jend(2) = jregsz(level)-1
77  jend(3) = jhi
78  jshift(1) = jregsz(level)
79  jshift(2) = 0
80  jshift(3) = -jregsz(level)
81  else
82  jst(1) = jregsz(level)
83  jst(2) = jlo
84  jst(3) = jregsz(level)
85  jend(1) = -nghost
86  jend(2) = jhi
87  jend(3) = -nghost
88  jshift(1) = 0
89  jshift(2) = 0
90  jshift(3) = 0
91  endif
92 
93 c ## loop over the 9 regions (in 2D) of the patch - the interior is i=j=2 plus
94 c ## the ghost cell regions. If any parts stick out of domain and are periodic
95 c ## map indices periodically, but stick the values in the correct place
96 c ## in the orig grid (indicated by iputst,jputst.
97 c ## if a region sticks out of domain but is not periodic, not handled in (pre)-icall
98 c ## but in setaux/bcamr (not called from here).
99  do 20 i = 1, 3
100  i1 = max(ilo, ist(i))
101  i2 = min(ihi, iend(i))
102  if (i1 .gt. i2) go to 20 ! non-empty intersection not possible
103  do 10 j = 1, 3
104  j1 = max(jlo, jst(j))
105  j2 = min(jhi, jend(j))
106 
107  if (j1 .le. j2) then ! part of patch in this region
108 c
109 c check if special mapping needed for spherical bc.
110 c (j=2 is interior,nothing special needed)
111  if (.not. spheredom .or. j .eq. 2) then
112  iputst = i1 - ilo + 1
113  jputst = j1 - jlo + 1
114  call icall(val,aux,nrow,ncol,nvar,naux,
115  1 i1+ishift(i),i2+ishift(i),
116  2 j1+jshift(j),j2+jshift(j),level,
117  3 iputst,jputst)
118  else
119  nr = i2 - i1 + 1
120  nc = j2 - j1 + 1
121  ng = 0 ! no ghost cells in this little patch. fill everything.
122 
123  jbump = 0
124  if (j1 < 0) jbump = abs(j1) ! bump up so new bottom is 0
125  if (j2 >= jregsz(level)) jbump = -(j2+1-jregsz(level)) ! bump down
126 
127 c next 2 lines would take care of periodicity in x
128  iwrap1 = i1 + ishift(i)
129  iwrap2 = i2 + ishift(i)
130 c next 2 lines take care of mapped sphere bcs
131  iwrap1 = iregsz(level) - iwrap1 -1
132  iwrap2 = iregsz(level) - iwrap2 -1
133 c swap so that smaller one is left index, etc since mapping reflects
134  tmp = iwrap1
135  iwrap1 = iwrap2
136  iwrap2 = tmp
137 
138  jwrap1 = j1 + jbump
139  xlwrap = xlower + iwrap1*hxposs(level)
140  ybwrap = ylower + jwrap1*hyposs(level)
141  jwrap2 = j2 + jbump
142 
143  if (naux>0) then
144 ! fliparray(locflipaux:locflipaux+naux*(ncol+nrow)-1) =
145  iflipchunksize = naux*nc*nr - 1 + nvar*(ncol+nrow)
146  idimen = (nrow+ncol)*nghost*(nvar+naux)
147  if (iflipchunksize .gt. idimen) then
148  write(*,*) "Error in fliparray size: asking for ",
149  . iflipchunksize," but dimension is",idimen
150  stop
151  endif
152  fliparray(locflipaux:locflipaux+naux*nc*nr - 1) =
154  call setaux(ng,nr,nc,xlwrap,ybwrap,
155  1 hxposs(level),hyposs(level),naux,
156  2 fliparray(locflipaux))
157  endif
158 
159 c write(dbugunit,101) i1,i2,j1,j2
160 c write(dbugunit6,102) iwrap1,iwrap2,j1+jbump,j2+jbump
161  101 format(" actual patch from i:",2i5," j :",2i5)
162  102 format(" icall called w i:",2i5," j :",2i5)
163  call icall(fliparray(locflip),fliparray(locflipaux),
164  1 nr,nc,nvar, naux,iwrap1,iwrap2,jwrap1,jwrap2,
165  2 level,1,1)
166 
167 c copy back using weird mapping for spherical folding
168  nrowst = 1 ! start filling up val at (1,1) - no additional offset
169  ncolst = 1
170  do 15 ii = i1, i2
171  do 15 jj = j1, j2
172 c write(dbugunit6,100)nrowst+ii-ilo,ncolst+jj-jlo,nr-(ii-i1),
173 c 1 nc-jj+j1
174  100 format(" filling loc ",2i5," with ",2i5)
175 
176  do 17 ivar = 1, nvar
177  val(ivar,nrowst+(ii-ilo),ncolst+(jj-jlo)) =
178  1 fliparray(iadd(ivar,nr-(ii-i1),nc-(jj-j1)))
179  17 continue
180 
181  do 16 iaux = 1, naux
182  aux(iaux,nrowst+(ii-ilo),ncolst+(jj-jlo)) =
183  1 fliparray(iaddaux(iaux,nr-(ii-i1),nc-(jj-j1)))
184  16 continue
185  15 continue
186 
187  endif
188 
189  endif
190 
191 
192  10 continue
193  20 continue
194 
195 
196  return
integer, dimension(maxlv) iregsz
Definition: amr_module.f90:198
real(kind=8), dimension(maxlv) hyposs
Definition: amr_module.f90:193
real(kind=8), dimension(maxlv) hxposs
Definition: amr_module.f90:193
integer, dimension(maxlv) jregsz
Definition: amr_module.f90:198
real(kind=8) xlower
Definition: amr_module.f90:231
integer pure function iadd(ivar, i, j)
Definition: intfil.f90:294
logical yperdom
Definition: amr_module.f90:230
subroutine icall(val, aux, nrow, ncol, nvar, naux, ilo, ihi, jlo, jhi, level, iputst, jputst)
For a rectangle defined on level level and bound by ilo, ihi, jlo, jhi, find intersecting grids at th...
Definition: icall.f:10
real(kind=8) ylower
Definition: amr_module.f90:231
real(kind=8), parameter needs_to_be_set
Definition: amr_module.f90:168
logical spheredom
Definition: amr_module.f90:230
logical xperdom
Definition: amr_module.f90:230
integer nghost
Definition: amr_module.f90:232
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
subroutine setaux(mbc, mx, my, xlower, ylower, dx, dy, maux, aux)
Definition: setaux.f90:2
Here is the call graph for this function:
Here is the caller graph for this function: