MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
amrgrid.t
Go to the documentation of this file.
1 !> Build up AMR
2 subroutine settree
5  use mod_advance, only: advance
6 
7  ! create and initialize grids on all levels > 1. On entry, all
8  ! level=1 grids have been formed and initialized. This subroutine
9  ! creates and initializes new level grids
10 
11  integer :: igrid, iigrid, levnew
12 
13  ! when only one level allowed, there is nothing to do anymore
14  if (refine_max_level == 1) return
15 
16  do levnew=2,refine_max_level
17  ! advanced needed for refinement based on comparing w_n, w_n-1
18  if(refine_criterion==1) then
19  call setdt
20  call advance(0)
21  end if
22 
23  call errest
24 
25  if(refine_criterion==1) then
26  do iigrid=1,igridstail; igrid=igrids(iigrid);
27  ps1(igrid)%w = pso(igrid)%w
28  pso(igrid)%w = ps(igrid)%w
29  ps(igrid)%w = ps1(igrid)%w
30  end do
31  end if
32 
34 
35  if (.not.reset_grid) then
36  ! if no finer level grids created: exit
37  if (levmax/=levnew) exit
38  end if
39  end do
40 
41  ! set up boundary flux conservation arrays
42  if(levmax>levmin) call allocatebflux
43 
44 end subroutine settree
45 
46 !> reset AMR and (de)allocate boundary flux storage at level changes
47 subroutine resettree
50  use mod_amr_fct
51 
54 
55  call errest
56 
58 
59  ! set up boundary flux conservation arrays
60  if(levmax>levmin) call allocatebflux
61 
62 end subroutine resettree
63 
64 !> Force the tree to desired level(s) from level_io(_min/_max)
65 !> used for conversion to vtk output
66 subroutine resettree_convert
69  integer :: igrid,iigrid, my_levmin, my_levmax
70 
71  if(level_io > 0) then
72  my_levmin = level_io
73  my_levmax = level_io
74  else
75  my_levmin = max(1,level_io_min)
76  my_levmax = min(refine_max_level,level_io_max)
77  end if
78 
79  do while(levmin<my_levmin.or.levmax>my_levmax)
80  call getbc(global_time,0.d0,ps,1,nwflux+nwaux)
81  do iigrid=1,igridstail; igrid=igrids(iigrid);
82  call forcedrefine_grid_io(igrid,ps(igrid)%w)
83  end do
84 
86  end do
87 
88 end subroutine resettree_convert
subroutine, public deallocatebflux
This module contains definitions of global parameters and variables and some generic functions/subrou...
update ghost cells of all blocks including physical boundaries
subroutine deallocatebfaces
Definition: mod_amr_fct.t:635
Module containing all the time stepping schemes.
Definition: mod_advance.t:2
subroutine forcedrefine_grid_io(igrid, w)
Definition: errest.t:426
subroutine errest
Do all local error estimation which determines (de)refinement.
Definition: errest.t:3
Module for flux conservation near refinement boundaries.
subroutine resettree
reset AMR and (de)allocate boundary flux storage at level changes
Definition: amrgrid.t:48
subroutine, public advance(iit)
Advance all the grids over one time step, including all sources.
Definition: mod_advance.t:20
subroutine amr_coarsen_refine
coarsen and refine blocks to update AMR grid
logical stagger_grid
True for using stagger grid.
subroutine settree
Build up AMR.
Definition: amrgrid.t:3
subroutine getbc(time, qdt, psb, nwstart, nwbc, req_diag)
do update ghost cells of all blocks including physical boundaries
double precision global_time
The global simulation time.
subroutine resettree_convert
Force the tree to desired level(s) from level_io(_min/_max) used for conversion to vtk output...
Definition: amrgrid.t:67
integer refine_max_level
Maximal number of AMR levels.
logical reset_grid
If true, rebuild the AMR grid upon restarting.
subroutine setdt()
setdt - set dt for all levels between levmin and levmax. dtpar>0 –> use fixed dtpar for all level dt...
Definition: setdt.t:5
subroutine, public allocatebflux