2D AMRCLAW
Functions/Subroutines
setgrd.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine setgrd (nvar, cut, naux, dtinit, start_time)
 

Function/Subroutine Documentation

◆ setgrd()

subroutine setgrd (   nvar,
  cut,
  naux,
  dtinit,
  start_time 
)

set up the entire tree/grid structure. only at this time t = 0 can we take advantage of initialization routines. remember that regridding/error estimation needs to have two time steps of soln. values.

Definition at line 5 of file setgrd.f.

References advanc(), arrangegrids(), amr_module::avenumgrids, amr_module::evol, amr_module::flag_richardson, ginit(), grdfit(), amr_module::iregridcount, amr_module::kratio, amr_module::levelptr, amr_module::lfine, amr_module::lstart, makebndrylist(), makegridlist(), amr_module::method, amr_module::mxnest, amr_module::ndihi, amr_module::ndilo, amr_module::ndjhi, amr_module::ndjlo, amr_module::newstl, amr_module::node, amr_module::numcells, amr_module::numgrids, amr_module::possk, prepc(), prepf(), amr_module::rnode, amr_module::rvol, amr_module::rvoll, amr_module::store1, amr_module::store2, amr_module::timemult, update(), and amr_module::verbosity_regrid.

Referenced by amr2().

5 c
6  use amr_module
7  implicit double precision (a-h,o-z)
8 
9  logical vtime
10  data vtime/.false./
11 
12 c # may as well not bother to calculate time step for error est.
13 c
14 c :::::::::::::::::::::::::::: SETGRD :::::::::::::::::::::::::::::::;
19 c 6/21/05: added dtinit arg. to allow for better choice of initial timestep
20 c as discovered by advance/setgrd in first step.
21 c ::::::::::::::::::::::::::::::::::::::;::::::::::::::::::::::::::::;
22 c
23  dtinit = possk(1)
24  if (mxnest .eq. 1) go to 99
25 c
26  levnew = 2
27  time = start_time
29 c
30  10 if (levnew .gt. mxnest) go to 30
31  levold = levnew - 1
32  if (lstart(levold) .eq. 0) go to 30
33  lbase = levold
34  lfnew = lbase
35 c
36 c set up level to be flagged. need a solution t=0,and if
37 c using richardson error estimation then one at dt too.
38 c then richardson makes one at 2*dt. for just flag2refine
39 c only need to set boundary conditions, done from grdfit->flagger.
40 c
41  if (flag_richardson) then
42  call advanc(levold,nvar,dtlev,vtime,naux)
43  evol = evol + rvol !time stepping stats for error estimation
44  rvol = 0.d0 ! reset since not 'real' time stepping stats
45  kfac = 1
46  do i = 1, levold-1
47  kfac = kfac * kratio(i)
48  end do
49  dtinit = min(dtinit, dtlev*kfac)
50 
51 c dont count it in real integration stats
52  do 20 level=1,mxnest
53  20 rvoll(level) = 0.d0
54  endif
55 c
56 c flag, cluster, and make new grids. grdfit set bcs, controls flagging,
57 c colating and making grids. But advanc called above since if using
58 c richardson, it assumes two levels of solution already exist
59 c
60  call grdfit(lbase,levold,nvar,naux,cut,time,start_time)
61  if (newstl(levnew) .ne. 0) lfnew = levnew
62 c
63 c init new level. after each iteration. fix the data structure
64 c also reinitalize coarser grids so fine grids can be advanced
65 c and interpolate correctly for their bndry vals from coarser grids.
66 c
67  call ginit(newstl(levnew),.true., nvar, naux,start_time)
68  lstart(levnew) = newstl(levnew)
69  lfine = lfnew
70  call ginit(lstart(levold),.false., nvar, naux,start_time)
71 c
72 c count number of grids on newly created levels (needed for openmp
73 c parallelization). this is also done in regridding.
74 c set up numgrids now for new level, since advanc will use it for parallel execution
75 c
76  mptr = lstart(levnew)
77  ngrids = 0
78  ncells = 0
79  do while (mptr .gt. 0)
80  ngrids = ngrids + 1
81  ncells = ncells + (node(ndihi,mptr)-node(ndilo,mptr)+1)
82  . * (node(ndjhi,mptr)-node(ndjlo,mptr)+1)
83  mptr = node(levelptr, mptr)
84  end do
85  numgrids(levnew) = ngrids
86  numcells(levnew) = ncells
87  avenumgrids(levnew) = avenumgrids(levnew) + ngrids
88  iregridcount(levnew) = iregridcount(levnew) + 1
89  if (ngrids .gt. 1) call arrangegrids(levnew,ngrids)
90  if (verbosity_regrid .ge. levnew) then
91  write(*,100) ngrids,ncells,levnew
92  100 format("there are ",i6," grids with ",i10,
93  & " cells at level ", i3)
94  endif
95 c
96 c need to make gridList here before calling again to make finer grids.
97 c This is because ths list if used in advanc. level 1 is called from domain
98  call makegridlist(lbase)
99  call makebndrylist(levnew)
100 c
101  levnew = levnew + 1
102  go to 10
103  30 continue
104 c
105 c switch location of old and new storage for soln. vals,
106 c and reset time to 0.0 (or initial time start_time)
107 c
108 c if (mxnest .eq. 1) go to 99 shouldnt be nec. tested above.
109 c
110  lev = 1
111  40 if ((lev .eq. mxnest) .or. (lev .gt. lfine)) go to 60
112  mptr = lstart(lev)
113  50 itemp = node(store1,mptr)
114  node(store1,mptr) = node(store2,mptr)
115  node(store2,mptr) = itemp
116  rnode(timemult,mptr) = start_time
117  mptr = node(levelptr,mptr)
118  if (mptr .ne. 0) go to 50
119  lev = lev + 1
120  go to 40
121  60 continue
122 c
123 c initial updating so can do conservation check. can do before
124 c bndry flux arrays set, since dont need them for this
125 c
126  do 65 level = 1, lfine-1
127  call update(lfine-level,nvar,naux)
128  65 continue
129 c
130 c set up boundary flux conservation arrays
131 c
132  do 70 level = 1, lfine-1
133  call prepf(level+1,nvar,naux)
134  call prepc(level,nvar)
135  70 continue
136 c
137 c
138 c
139 c
140  99 continue
141  return
subroutine makebndrylist(level)
Preprocess each grid on level level to have a linked list of other grids at the same level that suppl...
Definition: nodget.f:183
integer, dimension(7) method
Definition: amr_module.f90:253
integer, parameter timemult
current simulation time on this grid
Definition: amr_module.f90:151
integer, dimension(maxlv) kratio
Definition: amr_module.f90:198
integer, dimension(maxlv) numcells
Definition: amr_module.f90:198
subroutine prepc(level, nvar)
This routine is called because regridding just changed the fine grids.
Definition: prepc.f:12
real(kind=8) evol
Definition: amr_module.f90:237
subroutine makegridlist(lbase)
Make (modify) array of grid numbers, listOfGrids, (after sorting them in the linked list so they are ...
Definition: nodget.f:109
subroutine ginit(msave, first, nvar, naux, start_time)
Initializes solution on all grids at level of grid msave, by calling qinit().
Definition: ginit.f:18
integer, dimension(maxlv) newstl
Definition: amr_module.f90:198
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, dimension(maxlv) numgrids
Definition: amr_module.f90:198
real(kind=8), dimension(rsize, maxgr) rnode
Definition: amr_module.f90:193
integer, dimension(maxlv) iregridcount
Definition: amr_module.f90:238
logical flag_richardson
Definition: amr_module.f90:258
real(kind=8), dimension(maxlv) avenumgrids
Definition: amr_module.f90:237
subroutine update(level, nvar, naux)
Synchronize between all grids on level level and grids on level level+1.
Definition: update.f:17
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
integer, parameter store1
pointer to the address of memory storing the first copy of solution data on this grid, usually for storing new solution
Definition: amr_module.f90:101
subroutine advanc(level, nvar, dtlevnew, vtime, naux)
Integrate all grids at the input level by one step of its delta(t)
Definition: advanc.f:11
integer, dimension(maxlv) lstart
Definition: amr_module.f90:198
real(kind=8) rvol
Definition: amr_module.f90:237
subroutine grdfit(lbase, lcheck, nvar, naux, cut, time, start_time)
This subroutine is called by setgrd and regrid to actually fit the new grids on each level...
Definition: grdfit2.f:22
integer, parameter store2
pointer to the address of memory storing the second copy of solution data on this grid...
Definition: amr_module.f90:105
real(kind=8), dimension(maxlv) rvoll
Definition: amr_module.f90:237
real(kind=8), dimension(maxlv) possk
Definition: amr_module.f90:193
integer, parameter ndjhi
global j index of upper border of this grid
Definition: amr_module.f90:117
subroutine arrangegrids(level, numg)
Sort all grids at level level.
Definition: regrid.f:133
integer, parameter levelptr
node number (index) of next grid on the same level
Definition: amr_module.f90:35
integer mxnest
Definition: amr_module.f90:198
integer verbosity_regrid
Definition: amr_module.f90:259
The module contains the definition of a "node descriptor" as well as other global variables used duri...
Definition: amr_module.f90:21
integer lfine
Definition: amr_module.f90:198
subroutine prepf(level, nvar, naux)
For new fine grids (grids on level level), allocate space for saving fluxes at boundary (even if it's...
Definition: prepf.f:19
Here is the call graph for this function:
Here is the caller graph for this function: