MPI-AMRVAC  3.0
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
7 
8  ! create and initialize grids on all levels > 1. On entry, all
9  ! level=1 grids have been formed and initialized. This subroutine
10  ! creates and initializes new level grids
11 
12  integer :: igrid, iigrid, levnew
13 
14  ! when only one level allowed, there is nothing to do anymore
15  if (refine_max_level == 1) return
16 
17  do levnew=2,refine_max_level
18  ! advanced needed for refinement based on comparing w_n, w_n-1
19  if(refine_criterion==1) then
20  call setdt
21  call advance(0)
22  end if
23 
24  call errest
25 
26  if(refine_criterion==1) then
27  do iigrid=1,igridstail; igrid=igrids(iigrid);
28  ps1(igrid)%w = pso(igrid)%w
29  pso(igrid)%w = ps(igrid)%w
30  ps(igrid)%w = ps1(igrid)%w
31  end do
32  end if
33 
35 
36  if (.not.reset_grid) then
37  ! if no finer level grids created: exit
38  if (levmax/=levnew) exit
39  end if
40  end do
41 
42  ! set up boundary flux conservation arrays
43  if(levmax>levmin) call allocatebflux
44 
45 end subroutine settree
46 
47 !> reset AMR and (de)allocate boundary flux storage at level changes
48 subroutine resettree
51  use mod_amr_fct
53 
56 
57  call errest
58 
60 
61  ! set up boundary flux conservation arrays
62  if(levmax>levmin) call allocatebflux
63 
64 end subroutine resettree
65 
66 !> Force the tree to desired level(s) from level_io(_min/_max)
67 !> used for conversion to vtk output
72 
73  integer :: igrid,iigrid, my_levmin, my_levmax
74 
75  if(level_io > 0) then
76  my_levmin = level_io
77  my_levmax = level_io
78  else
79  my_levmin = max(1,level_io_min)
80  my_levmax = min(refine_max_level,level_io_max)
81  end if
82 
83  do while(levmin<my_levmin.or.levmax>my_levmax)
84  call getbc(global_time,0.d0,ps,iwstart,nwgc)
85  do iigrid=1,igridstail; igrid=igrids(iigrid);
86  call forcedrefine_grid_io(igrid,ps(igrid)%w)
87  end do
88 
90  end do
91 
92 end subroutine resettree_convert
93 
97 
98  if(phys_trac) then
99  if(phys_trac_type .eq. 3) then
100  if(mype .eq. 0) write(*,*) 'Using TRACL(ine) global method'
101  if(mype .eq. 0) write(*,*) 'By default, magnetic field lines are traced every 4 grid cells'
102  call init_trac_line(.false.)
103  end if
104  if(phys_trac_type .eq. 4) then
105  if(mype .eq. 0) write(*,*) 'Using TRACB(lock) global method'
106  if(mype .eq. 0) write(*,*) 'Currently, only valid in Cartesian uniform settings'
107  if(mype .eq. 0) write(*,*) 'By default, magnetic field lines are traced every 4 grid cells'
108  call init_trac_block(.false.)
109  end if
110  if(phys_trac_type .eq. 5) then
111  if(mype .eq. 0) write(*,*) 'Using TRACL(ine) method with a mask'
112  if(mype .eq. 0) write(*,*) 'By default, magnetic field lines are traced every 4 grid cells'
113  call init_trac_line(.true.)
114  end if
115  if(phys_trac_type .eq. 6) then
116  if(mype .eq. 0) write(*,*) 'Using TRACB(lock) method with a mask'
117  if(mype .eq. 0) write(*,*) 'Currently, only valid in Cartesian uniform settings'
118  if(mype .eq. 0) write(*,*) 'By default, magnetic field lines are traced every 4 grid cells'
119  call init_trac_block(.true.)
120  end if
121  end if
122 
123 end subroutine initialize_after_settree
subroutine resettree_convert
Force the tree to desired level(s) from level_io(_min/_max) used for conversion to vtk output.
Definition: amrgrid.t:69
subroutine initialize_after_settree
Definition: amrgrid.t:95
subroutine settree
Build up AMR.
Definition: amrgrid.t:3
subroutine resettree
reset AMR and (de)allocate boundary flux storage at level changes
Definition: amrgrid.t:49
subroutine errest
Do all local error estimation which determines (de)refinement.
Definition: errest.t:3
subroutine forcedrefine_grid_io(igrid, w)
Definition: errest.t:431
Module containing all the time stepping schemes.
Definition: mod_advance.t:2
subroutine, public advance(iit)
Advance all the grids over one time step, including all sources.
Definition: mod_advance.t:18
subroutine, public deallocatebfaces
Definition: mod_amr_fct.t:654
Module to coarsen and refine grids for AMR.
subroutine, public amr_coarsen_refine
coarsen and refine blocks to update AMR grid
Module for flux conservation near refinement boundaries.
subroutine, public deallocatebflux
subroutine, public allocatebflux
update ghost cells of all blocks including physical boundaries
subroutine getbc(time, qdt, psb, nwstart, nwbc, req_diag)
do update ghost cells of all blocks including physical boundaries
This module contains definitions of global parameters and variables and some generic functions/subrou...
double precision global_time
The global simulation time.
logical stagger_grid
True for using stagger grid.
integer mype
The rank of the current MPI task.
logical phys_trac
Use TRAC (Johnston 2019 ApJL, 873, L22) for MHD or 1D HD.
logical reset_grid
If true, rebuild the AMR grid upon restarting.
integer refine_max_level
Maximal number of AMR levels.
subroutine init_trac_block(mask)
Definition: mod_trac.t:78
subroutine init_trac_line(mask)
Definition: mod_trac.t:21
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