MPI-AMRVAC  3.1
The MPI - Adaptive Mesh Refinement - Versatile Advection Code (development version)
mod_refine.t
Go to the documentation of this file.
1 module mod_refine
2 
3  implicit none
4  private
5 
6  public :: refine_grids
7 
8 contains
9 
10 
11  !> refine one block to its children blocks
12  subroutine refine_grids(child_igrid,child_ipe,igrid,ipe,active)
16 
17  integer, dimension(2^D&), intent(in) :: child_igrid, child_ipe
18  integer, intent(in) :: igrid, ipe
19  logical, intent(in) :: active
20 
21  integer :: ic^d
22 
23  ! allocate solution space for new children
24  {do ic^db=1,2\}
25  call alloc_node(child_igrid(ic^d))
26  {end do\}
27 
28  if ((time_advance .and. active).or.convert.or.reset_grid) then
29  ! prolong igrid to new children
30  call prolong_grid(child_igrid,child_ipe,igrid,ipe)
31  else
32  ! Fill new created children with initial condition
33  {do ic^db=1,2\}
34  call initial_condition(child_igrid(ic^d))
35  {end do\}
36  end if
37 
38  ! remove solution space of igrid to save memory when converting data
39  if(convert) call dealloc_node(igrid)
40  end subroutine refine_grids
41 
42  !> prolong one block
43  subroutine prolong_grid(child_igrid,child_ipe,igrid,ipe)
46  use mod_amr_fct, only: old_neighbors
47 
48  integer, dimension(2^D&), intent(in) :: child_igrid, child_ipe
49  integer, intent(in) :: igrid, ipe
50 
51  integer :: ix^L, ichild, ixCo^L, ic^D
52  double precision :: dxCo^D, xComin^D, dxFi^D, xFimin^D
53 
54  ix^l=ixm^ll^ladd1;
55 
56  if(prolongprimitive) call phys_to_primitive(ixg^ll,ix^l,ps(igrid)%w,ps(igrid)%x)
57 
58  xcomin^d=rnode(rpxmin^d_,igrid)\
59  dxco^d=rnode(rpdx^d_,igrid)\
60 
61  if(stagger_grid) call old_neighbors(child_igrid,child_ipe,igrid,ipe)
62 
63  {do ic^db=1,2\}
64  ichild=child_igrid(ic^d)
65 
66  ixcomin^d=ixmlo^d+(ic^d-1)*block_nx^d/2\
67  ixcomax^d=ixmhi^d+(ic^d-2)*block_nx^d/2\
68 
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  !call prolong_1st(ps(igrid)%w,ixCo^L,ps(ichild)%w,ps(ichild)%x)
74  {end do\}
75 
76  if (prolongprimitive) call phys_to_conserved(ixg^ll,ix^l,ps(igrid)%w,ps(igrid)%x)
77 
78  end subroutine prolong_grid
79 
80  !> do 2nd order prolongation
81  subroutine prolong_2nd(sCo,ixCo^L,sFi,dxCo^D,xComin^D,dxFi^D,xFimin^D,igridCo,igridFi)
85 
86  integer, intent(in) :: ixCo^L, igridFi, igridCo
87  double precision, intent(in) :: dxCo^D, xComin^D, dxFi^D, xFimin^D
88  type(state), intent(in) :: sCo
89  type(state), intent(inout) :: sFi
90 
91  integer :: ixCo^D, jxCo^D, hxCo^D, ixFi^D, ix^D, idim, iw, ixCg^L, el
92  double precision :: slopeL, slopeR, slopeC, signC, signR
93  double precision :: slope(nw,ndim)
94  double precision :: eta^D
95  logical :: fine_^L
96 
97  associate(wco=>sco%w, wfi=>sfi%w)
98  ixcg^l=ixco^l;
99  {do ixco^db = ixcg^lim^db
100  ! lower left grid index in finer child block
101  ixfi^db=2*(ixco^db-ixcomin^db)+ixmlo^db\}
102 
103  do idim=1,ndim
104  hxco^d=ixco^d-kr(^d,idim)\
105  jxco^d=ixco^d+kr(^d,idim)\
106 
107  do iw=1,nw
108  slopel=wco(ixco^d,iw)-wco(hxco^d,iw)
109  sloper=wco(jxco^d,iw)-wco(ixco^d,iw)
110  slopec=half*(sloper+slopel)
111 
112  ! get limited slope
113  signr=sign(one,sloper)
114  signc=sign(one,slopec)
115  !select case(prolong_limiter)
116  !case(1)
117  ! ! unlimited
118  ! slope(iw,idim)=slopeC
119  !case(2)
120  ! ! minmod
121  ! slope(iw,idim)=signR*max(zero,min(dabs(slopeR), &
122  ! signR*slopeL))
123  !case(3)
124  ! ! woodward
125  ! slope(iw,idim)=two*signR*max(zero,min(dabs(slopeR), &
126  ! signR*slopeL,signR*half*slopeC))
127  !case(4)
128  ! ! koren
129  ! slope(iw,idim)=signR*max(zero,min(two*signR*slopeL, &
130  ! (dabs(slopeR)+two*slopeL*signR)*third,two*dabs(slopeR)))
131  !case default
132  slope(iw,idim)=signc*max(zero,min(dabs(slopec), &
133  signc*slopel,signc*sloper))
134  !end select
135  end do
136  end do
137  ! cell-centered coordinates of coarse grid point
138  !^D&xCo^D=xCo({ixCo^DD},^D)
139  {do ix^db=ixfi^db,ixfi^db+1 \}
140  ! cell-centered coordinates of fine grid point
141  !^D&xFi^D=xFi({ix^DD},^D)
142  if(slab_uniform) then
143  ! normalized distance between fine/coarse cell center
144  ! in coarse cell: ranges from -0.5 to 0.5 in each direction
145  ! (origin is coarse cell center)
146  ! hence this is +1/4 or -1/4 on cartesian mesh
147  !eta^D=(xFi^D-xCo^D)*invdxCo^D;
148  eta^d=0.5d0*(dble(ix^d-ixfi^d)-0.5d0);
149  else
150  {! forefactor is -0.5d0 when ix=ixFi and +0.5d0 for ixFi+1
151  eta^d=(dble(ix^d-ixfi^d)-0.5d0)*(one-sfi%dvolume(ix^dd) &
152  /sum(sfi%dvolume(ixfi^d:ixfi^d+1^d%ix^dd))) \}
153  end if
154  wfi(ix^d,1:nw) = wco(ixco^d,1:nw) &
155  + {(slope(1:nw,^d)*eta^d)+}
156  {end do\}
157  {end do\}
158  if(stagger_grid) then
159  call already_fine(sfi,igridfi,fine_^l)
160  call prolong_2nd_stg(sco,sfi,ixco^l,ixm^ll,dxco^d,xcomin^d,dxfi^d,xfimin^d,.false.,fine_^l)
161  end if
162 
163  if(fix_small_values) call phys_handle_small_values(prolongprimitive,wfi,sfi%x,ixg^ll,ixm^ll,'prolong_2nd')
164  if(prolongprimitive) call phys_to_conserved(ixg^ll,ixm^ll,wfi,sfi%x)
165  end associate
166 
167  end subroutine prolong_2nd
168 
169  !> do 1st order prolongation
170  subroutine prolong_1st(wCo,ixCo^L,wFi,xFi)
172 
173  integer, intent(in) :: ixCo^L
174  double precision, intent(in) :: wCo(ixG^T,nw), xFi(ixG^T,1:ndim)
175  double precision, intent(out) :: wFi(ixG^T,nw)
176 
177  integer :: ixCo^D, ixFi^D, iw
178  integer :: ixFi^L
179 
180  {do ixco^db = ixco^lim^db
181  ixfi^db=2*(ixco^db-ixcomin^db)+ixmlo^db\}
182  forall(iw=1:nw) wfi(ixfi^d:ixfi^d+1,iw)=wco(ixco^d,iw)
183  {end do\}
184 
185  end subroutine prolong_1st
186 
187 end module mod_refine
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
subroutine, public dealloc_node(igrid)
subroutine, public alloc_node(igrid)
allocate arrays on igrid node
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.
integer block_nx
number of cells for each dimension in grid block excluding ghostcells
integer, dimension(:), allocatable, parameter d
integer ixm
the mesh range of a physical block without ghost cells
double precision, dimension(:,:), allocatable rnode
Corner coordinates.
logical slab_uniform
uniform Cartesian geometry or not (stretched Cartesian)
subroutine, public initial_condition(igrid)
fill in initial condition
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:82
procedure(sub_convert), pointer phys_to_conserved
Definition: mod_physics.t:58
subroutine prolong_1st(wCo, ixCoL, wFi, xFi)
do 1st order prolongation
Definition: mod_refine.t:171
subroutine prolong_grid(child_igrid, child_ipe, igrid, ipe)
prolong one block
Definition: mod_refine.t:44
subroutine prolong_2nd(sCo, ixCoL, sFi, dxCoD, xCominD, dxFiD, xFiminD, igridCo, igridFi)
do 2nd order prolongation
Definition: mod_refine.t:82
subroutine, public refine_grids(child_igrid, child_ipe, igrid, ipe, active)
refine one block to its children blocks
Definition: mod_refine.t:13