MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
refine.t
Go to the documentation of this file.
1 !> refine one block to its children blocks
2 subroutine refine_grids(child_igrid,child_ipe,igrid,ipe,active)
4 
5  integer, dimension(2^D&), intent(in) :: child_igrid, child_ipe
6  integer, intent(in) :: igrid, ipe
7  logical, intent(in) :: active
8 
9  integer :: ic^D
10 
11  ! allocate solution space for new children
12  {do ic^db=1,2\}
13  call alloc_node(child_igrid(ic^d))
14  {end do\}
15 
16  if ((time_advance .and. active).or.convert.or.firstprocess) then
17  ! prolong igrid to new children
18  call prolong_grid(child_igrid,child_ipe,igrid,ipe)
19  else
20  ! Fill new created children with initial condition
21  {do ic^db=1,2\}
22  call initial_condition(child_igrid(ic^d))
23  {end do\}
24  end if
25 
26  ! remove solution space of igrid
27  !call dealloc_node(igrid)
28 end subroutine refine_grids
29 
30 !> prolong one block
31 subroutine prolong_grid(child_igrid,child_ipe,igrid,ipe)
34  use mod_amr_fct, only: old_neighbors
35 
36  integer, dimension(2^D&), intent(in) :: child_igrid, child_ipe
37  integer, intent(in) :: igrid, ipe
38 
39  integer :: ix^L, ichild, ixCo^L, ic^D
40  double precision :: dxCo^D, xComin^D, dxFi^D, xFimin^D
41 
42  if (prolongation_method=="linear") then
43  ^d&dxlevel(^d)=rnode(rpdx^d_,igrid);
44  {#IFDEF EVOLVINGBOUNDARY
45  if (phyboundblock(igrid)) then
46  ix^l=ixg^ll;
47  else
48  ix^l=ixm^ll^ladd1;
49  end if
50  }{#IFNDEF EVOLVINGBOUNDARY
51  ix^l=ixm^ll^ladd1;
52  }
53 
54  if(prolongprimitive) call phys_to_primitive(ixg^ll,ix^l,ps(igrid)%w,ps(igrid)%x)
55 
56  xcomin^d=rnode(rpxmin^d_,igrid)\
57  dxco^d=rnode(rpdx^d_,igrid)\
58  end if
59 
60  if(stagger_grid) call old_neighbors(child_igrid,child_ipe,igrid,ipe)
61 
62  {do ic^db=1,2\}
63  ichild=child_igrid(ic^d)
64 
65  ixcomin^d=ixmlo^d+(ic^d-1)*block_nx^d/2\
66  ixcomax^d=ixmhi^d+(ic^d-2)*block_nx^d/2\
67 
68  if (prolongation_method=="linear") then
69  xfimin^d=rnode(rpxmin^d_,ichild)\
70  dxfi^d=rnode(rpdx^d_,ichild)\
71  call prolong_2nd(ps(igrid),ixco^l,ps(ichild), &
72  dxco^d,xcomin^d,dxfi^d,xfimin^d,igrid,ichild)
73  else
74  call prolong_1st(ps(igrid)%w,ixco^l,ps(ichild)%w,ps(ichild)%x)
75  end if
76  {end do\}
77 
78  if (prolongation_method=="linear" .and. prolongprimitive) then
79  call phys_to_conserved(ixg^ll,ix^l,ps(igrid)%w,ps(igrid)%x)
80  end if
81 
82 end subroutine prolong_grid
83 
84 !> do 2nd order prolongation
85 subroutine prolong_2nd(sCo,ixCo^L,sFi,dxCo^D,xComin^D,dxFi^D,xFimin^D,igridCo,igridFi)
89 
90  integer, intent(in) :: ixCo^L, igridFi, igridCo
91  double precision, intent(in) :: dxCo^D, xComin^D, dxFi^D, xFimin^D
92  type(state), intent(in) :: sCo
93  type(state), intent(inout) :: sFi
94 
95  integer :: ixCo^D, jxCo^D, hxCo^D, ixFi^D, ix^D, idim, iw, ixCg^L, el
96  double precision :: slopeL, slopeR, slopeC, signC, signR
97  double precision :: slope(nw,ndim)
98  double precision :: eta^D
99  logical :: fine_^L
100 
101  associate(wco=>sco%w, wfi=>sfi%w)
102  ixcg^l=ixco^l;
103  {#IFDEF EVOLVINGBOUNDARY
104  if (phyboundblock(ichild)) then
105  el=ceiling(real(nghostcells)/2.d0)
106  ixcgmin^d=ixcomin^d-el\
107  ixcgmax^d=ixcomax^d+el\
108  end if
109  }
110  {do ixco^db = ixcg^lim^db
111  ! lower left grid index in finer child block
112  ixfi^db=2*(ixco^db-ixcomin^db)+ixmlo^db\}
113 
114  do idim=1,ndim
115  hxco^d=ixco^d-kr(^d,idim)\
116  jxco^d=ixco^d+kr(^d,idim)\
117 
118  do iw=1,nw
119  slopel=wco(ixco^d,iw)-wco(hxco^d,iw)
120  sloper=wco(jxco^d,iw)-wco(ixco^d,iw)
121  slopec=half*(sloper+slopel)
122 
123  ! get limited slope
124  signr=sign(one,sloper)
125  signc=sign(one,slopec)
126  select case(typeprolonglimit)
127  case('unlimit')
128  slope(iw,idim)=slopec
129  case('minmod')
130  slope(iw,idim)=signr*max(zero,min(dabs(sloper), &
131  signr*slopel))
132  case('woodward')
133  slope(iw,idim)=two*signr*max(zero,min(dabs(sloper), &
134  signr*slopel,signr*half*slopec))
135  case('koren')
136  slope(iw,idim)=signr*max(zero,min(two*signr*slopel, &
137  (dabs(sloper)+two*slopel*signr)*third,two*dabs(sloper)))
138  case default
139  slope(iw,idim)=signc*max(zero,min(dabs(slopec), &
140  signc*slopel,signc*sloper))
141  end select
142  end do
143  end do
144  ! cell-centered coordinates of coarse grid point
145  !^D&xCo^D=xCo({ixCo^DD},^D)\
146  {do ix^db=ixfi^db,ixfi^db+1 \}
147  ! cell-centered coordinates of fine grid point
148  !^D&xFi^D=xFi({ix^DD},^D)\
149  if(slab_uniform) then
150  ! normalized distance between fine/coarse cell center
151  ! in coarse cell: ranges from -0.5 to 0.5 in each direction
152  ! (origin is coarse cell center)
153  ! hence this is +1/4 or -1/4 on cartesian mesh
154  !eta^D=(xFi^D-xCo^D)*invdxCo^D;
155  eta^d=0.5d0*(dble(ix^d-ixfi^d)-0.5d0);
156  else
157  {! forefactor is -0.5d0 when ix=ixFi and +0.5d0 for ixFi+1
158  eta^d=(dble(ix^d-ixfi^d)-0.5d0)*(one-sfi%dvolume(ix^dd) &
159  /sum(sfi%dvolume(ixfi^d:ixfi^d+1^d%ix^dd))) \}
160  end if
161  wfi(ix^d,1:nw) = wco(ixco^d,1:nw) &
162  + {(slope(1:nw,^d)*eta^d)+}
163  {end do\}
164  {end do\}
165  if(stagger_grid) then
166  call already_fine(sfi,igridfi,fine_^l)
167  call prolong_2nd_stg(sco,sfi,ixco^l,ixm^ll,dxco^d,xcomin^d,dxfi^d,xfimin^d,.false.,fine_^l)
168  end if
169 
170  if(prolongprimitive) call phys_to_conserved(ixg^ll,ixm^ll,wfi,sfi%x)
171  end associate
172 
173 end subroutine prolong_2nd
174 
175 !> do 1st order prolongation
176 subroutine prolong_1st(wCo,ixCo^L,wFi,xFi)
178 
179  integer, intent(in) :: ixCo^L
180  double precision, intent(in) :: wCo(ixg^t,nw), xFi(ixg^t,1:ndim)
181  double precision, intent(out) :: wFi(ixg^t,nw)
182 
183  integer :: ixCo^D, ixFi^D, iw
184  integer :: ixFi^L
185 
186  {do ixco^db = ixco^lim^db
187  ixfi^db=2*(ixco^db-ixcomin^db)+ixmlo^db\}
188  forall(iw=1:nw) wfi(ixfi^d:ixfi^d+1,iw)=wco(ixco^d,iw)
189  {end do\}
190 
191 end subroutine prolong_1st
This module contains definitions of global parameters and variables and some generic functions/subrou...
character(len=std_len) prolongation_method
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
character(len=std_len) typeprolonglimit
Limiter used for prolongation to refined grids and ghost cells.
logical convert
If true and restart_from_file is given, convert snapshots to other file formats.
integer nghostcells
Number of ghost cells surrounding a grid.
logical stagger_grid
True for using stagger grid.
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:60
integer ixm
the mesh range (within a block with ghost cells)
subroutine, public prolong_2nd_stg(sCo, sFi, ixCoLin, ixFiLin, dxCoD, xCominD, dxFiD, xFiminD, ghost, fine_Lin)
This subroutine performs a 2nd order prolongation for a staggered field F, preserving the divergence ...
Definition: mod_amr_fct.t:41
subroutine alloc_node(igrid)
allocate arrays on igrid node
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer block_nx
number of cells for each dimension in grid block excluding ghostcells
subroutine refine_grids(child_igrid, child_ipe, igrid, ipe, active)
refine one block to its children blocks
Definition: refine.t:3
logical firstprocess
If true, call initonegrid_usr upon restarting.
subroutine prolong_1st(wCo, ixCoL, wFi, xFi)
do 1st order prolongation
Definition: refine.t:177
double precision, dimension(ndim) dxlevel
subroutine prolong_2nd(sCo, ixCoL, sFi, dxCoD, xCominD, dxFiD, xFiminD, igridCo, igridFi)
do 2nd order prolongation
Definition: refine.t:86
subroutine initial_condition(igrid)
fill in initial condition
Definition: amrini.t:39
integer, parameter ndim
Number of spatial dimensions for grid variables.
subroutine prolong_grid(child_igrid, child_ipe, igrid, ipe)
prolong one block
Definition: refine.t:32
subroutine, public old_neighbors(child_igrid, child_ipe, igrid, ipe)
Definition: mod_amr_fct.t:668
This module defines the procedures of a physics module. It contains function pointers for the various...
Definition: mod_physics.t:4
procedure(sub_convert), pointer phys_to_primitive
Definition: mod_physics.t:61
subroutine, public already_fine(sFi, ichild, fine_L)
This routine fills the fine faces before prolonging. It is the face equivalent of fix_conserve...
Definition: mod_amr_fct.t:691
logical, dimension(:), allocatable phyboundblock
True if a block has any physical boundary.
double precision, dimension(:,:), allocatable rnode
Corner coordinates.