2D AMRCLAW
Functions/Subroutines
domain.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine domain (nvar, vtime, nx, ny, naux, start_time)
 

Function/Subroutine Documentation

◆ domain()

subroutine domain (   nvar,
logical  vtime,
  nx,
  ny,
  naux,
  start_time 
)

allocate initial coarse grid domain. set node info & initialize grid initial space and time step set here too

Definition at line 5 of file domain.f.

References amr_module::alloc, arrangegrids(), amr_module::avenumgrids, birect(), amr_module::cfl, amr_module::cornxhi, amr_module::cornxlo, amr_module::cornyhi, amr_module::cornylo, estdt(), amr_module::flag_richardson, ginit(), amr_module::hxposs, amr_module::hyposs, amr_module::intratx, amr_module::intraty, amr_module::iregend, amr_module::iregridcount, amr_module::iregst, amr_module::iregsz, amr_module::jregend, amr_module::jregst, amr_module::jregsz, amr_module::kratio, amr_module::levelptr, amr_module::lfine, amr_module::liststart, amr_module::lstart, makebndrylist(), makegridlist(), amr_module::mstart, amr_module::mxnest, amr_module::ndihi, amr_module::ndilo, amr_module::ndjhi, amr_module::ndjlo, amr_module::nestlevel, amr_module::nghost, amr_module::node, nodget(), amr_module::numcells, amr_module::numgrids, amr_module::outunit, amr_module::possk, amr_module::rnode, amr_module::store1, amr_module::storeaux, amr_module::xlower, amr_module::xupper, amr_module::ylower, and amr_module::yupper.

Referenced by amr2(), drivesort(), and flagcheck().

5 c
6  use amr_module
7  implicit double precision (a-h,o-z)
8  logical vtime
9 
10 c
13 c
14  mstart = nodget()
15 c
16 c code assumes in many places that lower left corner at (0,0)
17 c this initial code sets the domain - assumed rectangular
18 c if it is too large, birect will chop it up into several rectangular
19 c pieces
20 c
25  node(nestlevel,mstart) = 1
26  node(levelptr,mstart) = 0
27  lstart(1) = mstart
28 
29  if (flag_richardson) then
30  if (((nx/2)*2 .ne. nx) .or. (ny/2)*2 .ne. ny) then
31  write(outunit,*)" must have even number of cells"
32  write(*,*) " must have even number of cells"
33  stop
34  endif
35  endif
36 
37  node(ndilo,mstart) = 0
38  node(ndjlo,mstart) = 0
39  node(ndihi,mstart) = nx-1
40  node(ndjhi,mstart) = ny-1
41 
42  lfine = 1
43  call birect(mstart)
44  call ginit (mstart, .true., nvar, naux, start_time)
45 c
46 c compute number of grids at level 1 (may have been bi-rected above)
47 c needs to be done here since this is used when calling advanc for
48 c parallelization
49  ngrids = 0
50  ncells = 0
51  mptr = lstart(1)
52  do while (mptr .gt. 0)
53  ngrids = ngrids + 1
54  ncells = ncells + (node(ndihi,mptr)-node(ndilo,mptr)+1)
55  & * (node(ndjhi,mptr)-node(ndjlo,mptr)+1)
56  mptr = node(levelptr, mptr)
57  end do
58  numgrids(1) = ngrids
59  numcells(1) = ncells
60  avenumgrids(1) = avenumgrids(1) + ngrids
61  iregridcount(1) = 1
62  if (ngrids .gt. 1) call arrangegrids(1,ngrids)
63 
64  write(*,100) ngrids,ncells
65  100 format("there are ",i4," grids with ",i8," cells at level 1")
66 
67 c set lbase to 1 here, to put domain 1 grids in lsit
68 c once and for all. Only here, this once, (and if restarting)
69 c does listStart have to be set outside of makeGridList
70 c but call it with lbase 0 to make grid 1
71  liststart(1) = 1
72  call makegridlist(0)
73  call makebndrylist(1) ! 1 means level 1
74 c
75 c set stable initial time step using coarse grid data
76 c
77  if (vtime) then
78  mptr = lstart(1)
79  dx = hxposs(1)
80  dy = hyposs(1)
81  dt = possk(1)
82  dtgrid = dt
83  60 mitot = node(ndihi,mptr)-node(ndilo,mptr) + 1 + 2*nghost
84  mjtot = node(ndjhi,mptr)-node(ndjlo,mptr) + 1 + 2*nghost
85  locaux = node(storeaux,mptr)
86 c # added cfl to call to estdt so call.i isnt needed in estdt:
87  call estdt(alloc(node(store1,mptr)),mitot,mjtot,nvar,
88  1 dx,dy,dtgrid,nghost,alloc(locaux),naux,cfl)
89  dt = dmin1(dt,dtgrid)
90  mptr = node(levelptr,mptr)
91  if (mptr .ne. 0) go to 60
92  possk(1) = dt
93  endif
94 c
95 c set rest of possk array for refined timesteps
96 c
97  iregsz(1) = nx
98  jregsz(1) = ny
99  iregst(1) = 0
100  jregst(1) = 0
101  iregend(1) = nx-1
102  jregend(1) = ny-1
103  do 70 level = 2, mxnest
104  iregsz(level) = iregsz(level-1) * intratx(level-1)
105  jregsz(level) = jregsz(level-1) * intraty(level-1)
106  possk(level) = possk(level-1)/dble(kratio(level-1))
107  70 continue
108 c
109  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(maxlv) kratio
Definition: amr_module.f90:198
integer, parameter cornxlo
x-coordinate of the left border of this grid
Definition: amr_module.f90:143
integer, dimension(maxlv) numcells
Definition: amr_module.f90:198
real(kind=8) cfl
Definition: amr_module.f90:255
subroutine makegridlist(lbase)
Make (modify) array of grid numbers, listOfGrids, (after sorting them in the linked list so they are ...
Definition: nodget.f:109
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
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
subroutine estdt(val, mitot, mjtot, nvar, dx, dy, dt, nghost, aux, naux, cfl)
Definition: estdt.f:6
integer, parameter ndihi
global i index of right border of this grid
Definition: amr_module.f90:111
integer, dimension(maxlv) jregsz
Definition: amr_module.f90:198
integer, dimension(nsize, maxgr) node
Definition: amr_module.f90:198
integer, dimension(maxlv) iregst
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
integer, dimension(0:maxlv+1) liststart
Definition: amr_module.f90:189
integer, dimension(maxlv) numgrids
Definition: amr_module.f90:198
real(kind=8), dimension(rsize, maxgr) rnode
Definition: amr_module.f90:193
real(kind=8) xlower
Definition: amr_module.f90:231
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
real(kind=8) yupper
Definition: amr_module.f90:231
real(kind=8) ylower
Definition: amr_module.f90:231
integer, dimension(maxlv) jregst
Definition: amr_module.f90:198
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
integer, dimension(maxlv) lstart
Definition: amr_module.f90:198
integer, dimension(maxlv) intraty
Definition: amr_module.f90:198
integer function nodget()
Get first free node of the linked list kept in node array.
Definition: nodget.f:11
integer, parameter outunit
Definition: amr_module.f90:290
real(kind=8), dimension(maxlv) possk
Definition: amr_module.f90:193
integer, dimension(maxlv) jregend
Definition: amr_module.f90:198
integer, parameter ndjhi
global j index of upper border of this grid
Definition: amr_module.f90:117
subroutine birect(mptr1)
Check each grid, starting with mptr1 (either newstl or lstart) to see that it has no more than max1d ...
Definition: birect.f:11
subroutine arrangegrids(level, numg)
Sort all grids at level level.
Definition: regrid.f:133
integer, parameter cornylo
y-coordinate of the lower border of this grid
Definition: amr_module.f90:145
integer, dimension(maxlv) intratx
Definition: amr_module.f90:198
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, parameter cornxhi
x-coordinate of the right border of this grid
Definition: amr_module.f90:147
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
integer, dimension(maxlv) iregend
Definition: amr_module.f90:198
integer mstart
Definition: amr_module.f90:198
integer lfine
Definition: amr_module.f90:198
integer, parameter cornyhi
y-coordinate of the upper border of this grid
Definition: amr_module.f90:149
integer, parameter storeaux
pointer to the address of memory storing auxiliary data on this grid
Definition: amr_module.f90:120
real(kind=8), dimension(:), allocatable alloc
Definition: amr_module.f90:218
Here is the call graph for this function:
Here is the caller graph for this function: