2D AMRCLAW
saveqc.f
Go to the documentation of this file.
1 c
2 c ================================================================
3  subroutine saveqc(level,nvar,naux)
4 c ================================================================
5 c
6  use amr_module
7  implicit double precision (a-h,o-z)
8 
9  logical sticksout, found
10 ! make fliparray largest possible grid size
11  dimension fliparray(2*max1d*nghost*(nvar+naux))
12 c
13 c ::::::::::::::::::::::::: SAVEQC :::::::::::::::::::::::::::::::::
14 c prepare new fine grids to save fluxes after each integration step
15 c for future conservative fix-up on coarse grids.
16 c save all boundary fluxes of fine grid (even if on a phys. bndry.) -
17 c but only save space for every intrat of them.
18 c:::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::::
19 c
20  levc = level - 1
21  hxc = hxposs(levc)
22  hyc = hyposs(levc)
23  ng = 0 ! no ghost cells on coarsened enlarged patch
24 
25  mkid = lstart(level)
26  10 if (mkid .eq. 0) go to 99
27  nx = node(ndihi,mkid)-node(ndilo,mkid) + 1
28  ny = node(ndjhi,mkid)-node(ndjlo,mkid) + 1
29  ikeep = nx/intratx(level-1)
30  jkeep = ny/intraty(level-1)
31  lenbc = 2*(ikeep+jkeep)
32  ist = node(ffluxptr,mkid)
33  time = rnode(timemult,mkid)
34 
35 c make coarsened enlarged patch for conservative fixup
36  ilo = node(ndilo,mkid)
37  jlo = node(ndjlo,mkid)
38  ihi = node(ndihi,mkid)
39  jhi = node(ndjhi,mkid)
40  iclo = ilo/intratx(level-1) - 1
41  jclo = jlo/intraty(level-1) - 1
42  ichi = (ihi+1)/intratx(level-1)
43  jchi = (jhi+1)/intraty(level-1)
44  nrow = ichi-iclo+1
45  ncol = jchi-jclo+1
46  xl = rnode(cornxlo,mkid) - hxc
47  yb = rnode(cornylo,mkid) - hyc
48  xr = rnode(cornxhi,mkid) + hxc
49  yt = rnode(cornyhi,mkid) + hyc
50  loctmp = igetsp(nrow*ncol*(nvar+naux))
51  loctx = loctmp + nrow*ncol*nvar
52  do i = 1, nrow*ncol*naux
53  alloc(loctx+i-1) = needs_to_be_set
54  end do
55  locaux = node(storeaux,mkid)
56 
57  if (iclo .lt. 0 .or. ichi .eq. iregsz(levc) .or.
58  & jclo .lt. 0 .or. jchi .eq. jregsz(levc)) then
59  sticksout = .true.
60  else
61  sticksout = .false.
62  endif
63 
64  if (sticksout .and. (xperdom.or.yperdom.or.spheredom)) then
65  !iperim = nrow+ncol
66  !locflip = igetsp(iperim*nghost*(nvar+naux))
67  call preicall(alloc(loctmp),alloc(loctx),nrow,ncol,nvar,
68  . naux,iclo,ichi,jclo,jchi,level-1,
69  . fliparray)
70 ! . alloc(locflip))
71 ! call reclam(locflip,iperim*nghost*(nvar+naux))
72  else
73  call icall(alloc(loctmp),alloc(loctx),nrow,ncol,nvar,naux,
74  . iclo,ichi,jclo,jchi,level-1,1,1)
75  endif
76 ! in case any part sticks out of domain still need to set remaining aux
77 ! cells
78  if (naux .gt. 0 .and. sticksout) then
79  call setaux(ng,nrow,ncol,xl,yb,hxc,hyc,naux,alloc(loctx))
80  endif
81 !-- found = .false.
82 !-- do i = 1, naux*nrow*ncol, naux
83 !-- if (alloc(loctx+i-1) .eq. NEEDS_TO_BE_SET) then
84 !-- found = .true.
85 !-- endif
86 !-- end do
87 !-- if (found) write(*,*) "still have unset aux vals in qad"
88  call bc2amr(alloc(loctmp),alloc(loctx),nrow,ncol,nvar,naux,
89  . hxc,hyc,level,time,
90  . xl,xr,yb,yt)
91  call cstore(alloc(loctmp),nrow,ncol,nvar,
92  . alloc(ist+nvar*lenbc),lenbc,naux,alloc(loctx),
93  . alloc(ist+2*nvar*lenbc))
94  call reclam(loctmp,nrow*ncol*(nvar+naux))
95 
96  mkid = node(levelptr,mkid)
97  go to 10
98  99 return
99  end
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:87
subroutine saveqc(level, nvar, naux)
Definition: saveqc.f:3
function igetsp(nwords)
Definition: igetsp.f:4
logical pure function sticksout(iplo, iphi, jplo, jphi)
Definition: filpatch.f90:348
subroutine reclam(index, nwords)
Definition: reclam.f:4
subroutine cstore(qc, nrow, ncol, nvar, qc1d, lenbc, naux, auxc, auxc1d)
Definition: cstore.f:4
subroutine icall(val, aux, nrow, ncol, nvar, naux, ilo, ihi, jlo, jhi, level, iputst, jputst)
Definition: icall.f:4
subroutine preicall(val, aux, nrow, ncol, nvar, naux, ilo, ihi, jlo, jhi, level, fliparray)
Definition: preicall.f:4
subroutine setaux(mbc, mx, my, xlower, ylower, dx, dy, maux, aux)
Definition: setaux.f90:1