2D AMRCLAW
filval.f90
Go to the documentation of this file.
1 
3 subroutine filval(val, mitot, mjtot, dx, dy, level, time, mic, &
4  mjc, xleft, xright, ybot, ytop, nvar, mptr, ilo, ihi, &
5  jlo, jhi, aux, naux)
6 
10  use amr_module, only: newstl, iregsz, jregsz
11 
12  implicit none
13 
14  ! Input
15  integer, intent(in) :: mitot, mjtot, level, mic, mjc, nvar, mptr, ilo, ihi
16  integer, intent(in) :: jlo, jhi, naux
17  real(kind=8), intent(in) :: dx, dy, time, xleft, xright, ybot, ytop
18 
19  ! Output
20  real(kind=8), intent(in out) :: val(nvar,mitot,mjtot), aux(naux,mitot,mjtot)
21 
22  ! Local storage
23  integer :: refinement_ratio_x, refinement_ratio_y, iclo, jclo, ichi, jchi, ng
24  integer :: ivar, i, j, ico, jco, ifine, jfine, nx, ny
25  real(kind=8) :: valc(nvar,mic,mjc), auxc(naux,mic,mjc)
26  real(kind=8) :: dx_coarse, dy_coarse, xl, xr, yb, yt, area
27  real(kind=8) :: s1m, s1p, slopex, slopey, xoff, yoff
28  real(kind=8) :: fliparray((mitot+mjtot)*nghost*(nvar+naux))
29  real(kind=8) :: setflags(mitot,mjtot),maxauxdif,aux2(naux,mitot,mjtot)
30  integer :: mjb
31  integer :: jm, im, nm
32  logical :: sticksoutxfine, sticksoutyfine,sticksoutxcrse,sticksoutycrse
33 
34  !for setaux timing
35  integer :: clock_start, clock_finish, clock_rate
36  real(kind=8) :: cpu_start, cpu_finish
37 
38 
39  ! External function definitions
40  real(kind=8) :: get_max_speed
41 
42 
43  refinement_ratio_x = intratx(level-1)
44  refinement_ratio_y = intraty(level-1)
45  dx_coarse = dx * refinement_ratio_x
46  dy_coarse = dy * refinement_ratio_y
47  xl = xleft - dx_coarse
48  xr = xright + dx_coarse
49  yb = ybot - dy_coarse
50  yt = ytop + dy_coarse
51 
52  ! if topo not yet final then aux is set outside filval (in gfixup)
53  ! and so aux has real data already, (ie dont overwrite here)
54 
55  ! set integer indices for coarser patch enlarged by 1 cell
56  ! (can stick out of domain). proper nesting will insure this one
57  ! call is sufficient.
58  iclo = ilo / refinement_ratio_x - 1
59  jclo = jlo / refinement_ratio_y - 1
60  ichi = (ihi + 1) / refinement_ratio_x - 1 + 1
61  jchi = (jhi + 1) / refinement_ratio_y - 1 + 1
62  ng = 0
63 
64  sticksoutxfine = ( (ilo .lt. 0) .or. (ihi .ge. iregsz(level)))
65  sticksoutyfine = ( (jlo .lt. 0) .or. (jhi .ge. jregsz(level)))
66  sticksoutxcrse = ((iclo .lt. 0) .or. (ichi .ge. iregsz(level-1)))
67  sticksoutycrse = ((jclo .lt. 0) .or. (jchi .ge. jregsz(level-1)))
68 
69  if (naux == 0) then
70  if ((xperdom .and. sticksoutxcrse).or. (yperdom.and.sticksoutycrse) .or. spheredom) then
71  call preintcopy(valc,mic,mjc,nvar,iclo,ichi,jclo,jchi,level-1,fliparray)
72  else
73  call intcopy(valc,mic,mjc,nvar,iclo,ichi,jclo,jchi,level-1,1,1)
74  endif
75  else
76  ! intersect grids and copy all (soln and aux)
77  auxc(1,:,:) = needs_to_be_set
78  if ((xperdom .and.sticksoutxcrse) .or. (yperdom.and.sticksoutycrse) .or. spheredom) then
79  call preicall(valc,auxc,mic,mjc,nvar,naux,iclo,ichi,jclo,jchi, &
80  level-1,fliparray)
81  else
82  call icall(valc,auxc,mic,mjc,nvar,naux,iclo,ichi,jclo,jchi,level-1,1,1)
83  endif
84 !!$ do i = 1, mic
85 !!$ do j = 1, mjc
86 !!$ if (auxc(1,:,:) == NEEDS_TO_BE_SET) then
87 !!$ write(*,*)" *** coarsenened new fine grid not completely set from previously" &
88 !!$ "existing coarse grids ***"
89 !!$ stop
90 !!$ endif
91 !!$ end do
92 !!$ end do
93  ! no ghost cells on coarse enlarged patch. set any remaining
94  ! vals. should only be bcs that stick out of domain.
95  call setaux(ng,mic,mjc,xl,yb,dx_coarse,dy_coarse,naux,auxc)
96  endif
97 
98  call bc2amr(valc,auxc,mic,mjc,nvar,naux,dx_coarse,dy_coarse,level-1,time,xl,xr,yb, &
99  yt,xlower,ylower,xupper,yupper,xperdom,yperdom,spheredom)
100 
101 
102 
103 ! NOTE change in order of code. Since the interp from coarse to fine needs the aux
104 ! arrays set already, the fine copy is done first, to set up the aux arrays.
105 ! we can do this since we have the flag array to test where to overwrite.
106 
107 ! SO this is no longer overwriting but setting for the first time.
108 ! overwrite interpolated values with fine grid values, if available.
109  nx = mitot - 2*nghost
110  ny = mjtot - 2*nghost
111 
112  if (naux .gt. 0) then
113 ! ## NEEDS_TO_BE_SET is signal that aux array not set.
114 ! ## after calling icall to copy aux from other grids
115 ! ## any remaining NEEDS_TO_BE_SET signals will be set in setaux.
116 ! ## it also signals where soln was copied, so it wont be
117 ! ## overwritten with coarse grid interpolation
118  aux(1,:,:) = needs_to_be_set ! indicates fine cells not yet set.
119 
120  if ((xperdom.and.sticksoutxfine) .or. (yperdom.and.sticksoutyfine) .or. spheredom) then
121  call preicall(val,aux,mitot,mjtot,nvar,naux,ilo-nghost,ihi+nghost, &
122  jlo-nghost,jhi+nghost,level,fliparray)
123  else
124  call icall(val,aux,mitot,mjtot,nvar,naux,ilo-nghost,ihi+nghost, &
125  jlo-nghost,jhi+nghost,level,1,1)
126  endif
127  setflags = aux(1,:,:) ! save since will overwrite in setaux when setting all aux vals
128  ! need this so we know where to use coarse grid to set fine solution w/o overwriting
129  ! set remaining aux vals not set by copying from prev existing grids
130  call setaux(nghost,nx,ny,xleft,ybot,dx,dy,naux,aux)
131  else ! either no aux exists, or cant reuse yet
132  ! so only call intcopy (which copies soln) and not icall.
133  ! in this case flag q(1,:) to NEEDS_TO_BE_SET flag so wont be overwritten
134  ! by coarse grid interp. this is needed due to reversing order of
135  ! work - first copy from fine grids, then interpolate from coarse grids
136  val(1,:,:) = needs_to_be_set
137  if ((xperdom.and.sticksoutxfine) .or. (yperdom.and.sticksoutyfine) .or. spheredom) then
138  call preintcopy(val,mitot,mjtot,nvar,ilo-nghost,ihi+nghost, &
139  jlo-nghost,jhi+nghost,level,fliparray)
140  else
141  call intcopy(val,mitot,mjtot,nvar,ilo-nghost,ihi+nghost, &
142  jlo-nghost,jhi+nghost,level,1,1)
143  endif
144  setflags = val(1,:,:) ! remaining flags signals need to set
145  endif
146 
147 
148  ! Prepare slopes - use min-mod limiters
149 
150  do j=2, mjc-1
151  do i=2, mic-1
152  do ivar = 1, nvar
153 
154  s1p = valc(ivar,i+1,j) - valc(ivar,i,j)
155  s1m = valc(ivar,i,j) - valc(ivar,i-1,j)
156  slopex = min(abs(s1p), abs(s1m)) &
157  * sign(1.d0,valc(ivar,i+1,j) - valc(ivar,i-1,j))
158  ! if there's a sign change, set slope to 0.
159  if ( s1m*s1p <= 0.d0) slopex = 0.d0
160 
161  s1p = valc(ivar,i,j+1) - valc(ivar,i,j)
162  s1m = valc(ivar,i,j) - valc(ivar,i,j-1)
163  slopey = min(abs(s1p), abs(s1m)) &
164  * sign(1.0d0, valc(ivar,i,j+1) - valc(ivar,i,j-1))
165  if ( s1m*s1p <= 0.d0) slopey = 0.d0
166 
167  ! Interpolate from coarse cells to fine grid to find depth
168  do jco = 1,refinement_ratio_y
169  yoff = (real(jco,kind=8) - 0.5d0) / refinement_ratio_y - 0.5d0
170  jfine = (j-2) * refinement_ratio_y + nghost + jco
171 
172  do ico = 1,refinement_ratio_x
173  xoff = (real(ico,kind=8) - 0.5d0) / refinement_ratio_x - 0.5d0
174  ifine = (i-2) * refinement_ratio_x + nghost + ico
175 
176  if (setflags(ifine,jfine) .eq. needs_to_be_set) then
177  val(ivar,ifine,jfine) = valc(ivar,i,j) + xoff*slopex + yoff*slopey
178  endif
179 
180  end do
181  end do
182 
183  enddo !end of ivar loop
184  enddo !end of coarse i loop
185  enddo !end of coarse j loop
186 
187  ! adjust to conserve kappa*q, but only where coarse grid was interpolated
188  ! so now need to pass setflags to this subr.
189  if (mcapa .ne. 0) then
190  call fixcapaq(val,aux,mitot,mjtot,valc,auxc,mic,mjc,nvar,naux,level-1,setflags)
191  endif
192 
193 !!$! CHECK BY CALLING SETAUX AND SETTING ALL, THEN DIFFING
194 !!$ mjb = 0
195 !!$ if (naux .gt. 0 .and. mjb .eq. 1) then
196 !!$ aux2(1,:,:) = NEEDS_TO_BE_SET ! indicates fine cells not yet set
197 !!$ call setaux(nghost,nx,ny,xleft,ybot,dx,dy,naux,aux2)
198 !!$ maxauxdif = 1.d-13
199 !!$ do i = 1, mitot
200 !!$ do j = 1, mjtot
201 !!$ if (abs(aux(1,i,j)-aux2(1,i,j)) .gt. maxauxdif) then
202 !!$ maxauxdif = abs(aux(1,i,j)-aux2(1,i,j))
203 !!$ write(*,444)i,j,aux(1,i,j),aux2(1,i,j),maxauxdif
204 !!$444 format("i,j = ",2i4," auxs ",2e15.7," maxauxdif ",e12.5)
205 !!$ endif
206 !!$ end do
207 !!$ end do
208 !!$ if (maxauxdif .gt. 2.d-13) then
209 !!$ write(*,*)" maxauxdif = ",maxauxdif," with mitot,mjtot ",mitot,mjtot, &
210 !!$ " on grid ",mptr," level ",level
211 !!$ endif
212 !!$ endif
213 
214 end subroutine filval
215 
216 ! ::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
217 
218 subroutine dumpaux(aux,naux,mitot,mjtot)
219  implicit none
220  real(kind=8) :: aux(naux,mitot,mjtot)
221  integer :: naux,mitot,mjtot,i,j,iaux
222 
223  do j = 1, mjtot
224  do i = 1, mitot
225  write(*,444) i,j,(aux(iaux,i,j),iaux=1,naux)
226  444 format(2i4,5e12.5)
227  end do
228  end do
229 
230 end subroutine dumpaux
subroutine bc2amr(val, aux, nrow, ncol, meqn, naux,
Take a grid patch with mesh widths hx,hy, of dimensions nrow by ncol, and set the values of any piece...
Definition: bc2amr.f:88
subroutine filval(val, mitot, mjtot, dx, dy, level, time, mic, mjc, xleft, xright, ybot, ytop, nvar, mptr, ilo, ihi, jlo, jhi, aux, naux)
Fill grid mptr on level level by copying values from OLD level level grids if available, otherwise by interpolating values from coarser grids.
Definition: filval.f90:6
integer, dimension(maxlv) iregsz
Definition: amr_module.f90:198
integer, dimension(maxlv) newstl
Definition: amr_module.f90:198
integer, dimension(maxlv) jregsz
Definition: amr_module.f90:198
real(kind=8) xupper
Definition: amr_module.f90:231
real(kind=8) xlower
Definition: amr_module.f90:231
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
subroutine intcopy(val, mitot, mjtot, nvar, ilo, ihi, jlo, jhi, level, iputst, jputst)
For a rectangle that is on level level, described by ilo, ihi, jlo, jhi and made up by mitot mjtot c...
Definition: intcopy.f:14
real(kind=8) yupper
Definition: amr_module.f90:231
real(kind=8) ylower
Definition: amr_module.f90:231
real(kind=8), parameter needs_to_be_set
Definition: amr_module.f90:168
subroutine fixcapaq(val, aux, mitot, mjtot, valc, auxc, mic, mjc, nvar, naux, levc, setflags)
Definition: fixcapaq.f:6
logical spheredom
Definition: amr_module.f90:230
subroutine preintcopy(val, mitot, mjtot, nvar, ilo, ihi, jlo, jhi, level, fliparray)
Definition: preintcopy.f:6
integer mcapa
Definition: amr_module.f90:253
integer, dimension(maxlv) intraty
Definition: amr_module.f90:198
subroutine dumpaux(aux, naux, mitot, mjtot)
Definition: filval.f90:219
subroutine preicall(val, aux, nrow, ncol, nvar, naux, ilo, ihi, jlo, jhi, level, fliparray)
Definition: preicall.f:6
integer, parameter outunit
Definition: amr_module.f90:290
logical xperdom
Definition: amr_module.f90:230
integer, dimension(maxlv) intratx
Definition: amr_module.f90:198
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
real(kind=8), dimension(:), allocatable alloc
Definition: amr_module.f90:218