MPI-AMRVAC  3.0
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.reset_grid) 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 to save memory when converting data
27  if(convert) 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  {#IFDEF EVOLVINGBOUNDARY
43  if (phyboundblock(igrid)) then
44  ix^l=ixg^ll;
45  else
46  ix^l=ixm^ll^ladd1;
47  end if
48  }{#IFNDEF EVOLVINGBOUNDARY
49  ix^l=ixm^ll^ladd1;
50  }
51 
52  if(prolongprimitive) call phys_to_primitive(ixg^ll,ix^l,ps(igrid)%w,ps(igrid)%x)
53 
54  xcomin^d=rnode(rpxmin^d_,igrid)\
55  dxco^d=rnode(rpdx^d_,igrid)\
56 
57  if(stagger_grid) call old_neighbors(child_igrid,child_ipe,igrid,ipe)
58 
59  {do ic^db=1,2\}
60  ichild=child_igrid(ic^d)
61 
62  ixcomin^d=ixmlo^d+(ic^d-1)*block_nx^d/2\
63  ixcomax^d=ixmhi^d+(ic^d-2)*block_nx^d/2\
64 
65  xfimin^d=rnode(rpxmin^d_,ichild)\
66  dxfi^d=rnode(rpdx^d_,ichild)\
67  call prolong_2nd(ps(igrid),ixco^l,ps(ichild), &
68  dxco^d,xcomin^d,dxfi^d,xfimin^d,igrid,ichild)
69  !call prolong_1st(ps(igrid)%w,ixCo^L,ps(ichild)%w,ps(ichild)%x)
70  {end do\}
71 
72  if (prolongprimitive) call phys_to_conserved(ixg^ll,ix^l,ps(igrid)%w,ps(igrid)%x)
73 
74 end subroutine prolong_grid
75 
76 !> do 2nd order prolongation
77 subroutine prolong_2nd(sCo,ixCo^L,sFi,dxCo^D,xComin^D,dxFi^D,xFimin^D,igridCo,igridFi)
81 
82  integer, intent(in) :: ixCo^L, igridFi, igridCo
83  double precision, intent(in) :: dxCo^D, xComin^D, dxFi^D, xFimin^D
84  type(state), intent(in) :: sCo
85  type(state), intent(inout) :: sFi
86 
87  integer :: ixCo^D, jxCo^D, hxCo^D, ixFi^D, ix^D, idim, iw, ixCg^L, el
88  double precision :: slopeL, slopeR, slopeC, signC, signR
89  double precision :: slope(nw,ndim)
90  double precision :: eta^D
91  logical :: fine_^L
92 
93  associate(wco=>sco%w, wfi=>sfi%w)
94  ixcg^l=ixco^l;
95  {#IFDEF EVOLVINGBOUNDARY
96  if (phyboundblock(ichild)) then
97  el=ceiling(real(nghostcells)/2.d0)
98  ixcgmin^d=ixcomin^d-el\
99  ixcgmax^d=ixcomax^d+el\
100  end if
101  }
102  {do ixco^db = ixcg^lim^db
103  ! lower left grid index in finer child block
104  ixfi^db=2*(ixco^db-ixcomin^db)+ixmlo^db\}
105 
106  do idim=1,ndim
107  hxco^d=ixco^d-kr(^d,idim)\
108  jxco^d=ixco^d+kr(^d,idim)\
109 
110  do iw=1,nw
111  slopel=wco(ixco^d,iw)-wco(hxco^d,iw)
112  sloper=wco(jxco^d,iw)-wco(ixco^d,iw)
113  slopec=half*(sloper+slopel)
114 
115  ! get limited slope
116  signr=sign(one,sloper)
117  signc=sign(one,slopec)
118  !select case(prolong_limiter)
119  !case(1)
120  ! ! unlimited
121  ! slope(iw,idim)=slopeC
122  !case(2)
123  ! ! minmod
124  ! slope(iw,idim)=signR*max(zero,min(dabs(slopeR), &
125  ! signR*slopeL))
126  !case(3)
127  ! ! woodward
128  ! slope(iw,idim)=two*signR*max(zero,min(dabs(slopeR), &
129  ! signR*slopeL,signR*half*slopeC))
130  !case(4)
131  ! ! koren
132  ! slope(iw,idim)=signR*max(zero,min(two*signR*slopeL, &
133  ! (dabs(slopeR)+two*slopeL*signR)*third,two*dabs(slopeR)))
134  !case default
135  slope(iw,idim)=signc*max(zero,min(dabs(slopec), &
136  signc*slopel,signc*sloper))
137  !end select
138  end do
139  end do
140  ! cell-centered coordinates of coarse grid point
141  !^D&xCo^D=xCo({ixCo^DD},^D)
142  {do ix^db=ixfi^db,ixfi^db+1 \}
143  ! cell-centered coordinates of fine grid point
144  !^D&xFi^D=xFi({ix^DD},^D)
145  if(slab_uniform) then
146  ! normalized distance between fine/coarse cell center
147  ! in coarse cell: ranges from -0.5 to 0.5 in each direction
148  ! (origin is coarse cell center)
149  ! hence this is +1/4 or -1/4 on cartesian mesh
150  !eta^D=(xFi^D-xCo^D)*invdxCo^D;
151  eta^d=0.5d0*(dble(ix^d-ixfi^d)-0.5d0);
152  else
153  {! forefactor is -0.5d0 when ix=ixFi and +0.5d0 for ixFi+1
154  eta^d=(dble(ix^d-ixfi^d)-0.5d0)*(one-sfi%dvolume(ix^dd) &
155  /sum(sfi%dvolume(ixfi^d:ixfi^d+1^d%ix^dd))) \}
156  end if
157  wfi(ix^d,1:nw) = wco(ixco^d,1:nw) &
158  + {(slope(1:nw,^d)*eta^d)+}
159  {end do\}
160  {end do\}
161  if(stagger_grid) then
162  call already_fine(sfi,igridfi,fine_^l)
163  call prolong_2nd_stg(sco,sfi,ixco^l,ixm^ll,dxco^d,xcomin^d,dxfi^d,xfimin^d,.false.,fine_^l)
164  end if
165 
166  if(fix_small_values) call phys_handle_small_values(prolongprimitive,wfi,sfi%x,ixg^ll,ixm^ll,'prolong_2nd')
167  if(prolongprimitive) call phys_to_conserved(ixg^ll,ixm^ll,wfi,sfi%x)
168  end associate
169 
170 end subroutine prolong_2nd
171 
172 !> do 1st order prolongation
173 subroutine prolong_1st(wCo,ixCo^L,wFi,xFi)
175 
176  integer, intent(in) :: ixCo^L
177  double precision, intent(in) :: wCo(ixG^T,nw), xFi(ixG^T,1:ndim)
178  double precision, intent(out) :: wFi(ixG^T,nw)
179 
180  integer :: ixCo^D, ixFi^D, iw
181  integer :: ixFi^L
182 
183  {do ixco^db = ixco^lim^db
184  ixfi^db=2*(ixco^db-ixcomin^db)+ixmlo^db\}
185  forall(iw=1:nw) wfi(ixfi^d:ixfi^d+1,iw)=wco(ixco^d,iw)
186  {end do\}
187 
188 end subroutine prolong_1st
subroutine alloc_node(igrid)
allocate arrays on igrid node
subroutine dealloc_node(igrid)
subroutine initial_condition(igrid)
fill in initial condition
Definition: amrini.t:39
subroutine, public old_neighbors(child_igrid, child_ipe, igrid, ipe)
Definition: mod_amr_fct.t:668
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
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
This module contains definitions of global parameters and variables and some generic functions/subrou...
integer, dimension(3, 3) kr
Kronecker delta tensor.
logical stagger_grid
True for using stagger grid.
logical, dimension(:), allocatable phyboundblock
True if a block has any physical boundary.
integer block_nx
number of cells for each dimension in grid block excluding ghostcells
integer ixm
the mesh range (within a block with ghost cells)
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
integer nghostcells
Number of ghost cells surrounding a grid.
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
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:59
procedure(sub_small_values), pointer phys_handle_small_values
Definition: mod_physics.t:81
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:58
subroutine prolong_2nd(sCo, ixCoL, sFi, dxCoD, xCominD, dxFiD, xFiminD, igridCo, igridFi)
do 2nd order prolongation
Definition: refine.t:78
subroutine prolong_grid(child_igrid, child_ipe, igrid, ipe)
prolong one block
Definition: refine.t:32
subroutine refine_grids(child_igrid, child_ipe, igrid, ipe, active)
refine one block to its children blocks
Definition: refine.t:3
subroutine prolong_1st(wCo, ixCoL, wFi, xFi)
do 1st order prolongation
Definition: refine.t:174