MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
amr_neighbors.t
Go to the documentation of this file.
1 !> find neighors of level-one root blocks
2 subroutine find_root_neighbor(tree_neighbor,tree,i^D)
5  use mod_geometry
6 
7  type(tree_node_ptr) :: tree_neighbor, tree
8  integer, intent(in) :: i^D
9 
10  integer :: jg^D
11 
12  jg^d=tree%node%ig^d+i^d;
13 
14  ! find the periodic grid indices, modulo(-1,10)=9
15  {if (periodb(^d)) jg^d=1+modulo(jg^d-1,ng^d(1))\}
16 
17  ! pi-periodicity at pole
18  select case (coordinate)
19  case (spherical) {^ifthreed
20  if (poleb(1,2).and.jg2==0) then ! northpole (theta=0)
21  jg2=1; jg3=1+modulo(jg3+ng3(1)/2-1,ng3(1))
22  end if
23  if (poleb(2,2).and.jg2==ng2(1)+1) then ! southpole (theta=pi)
24  jg2=ng2(1); jg3=1+modulo(jg3+ng3(1)/2-1,ng3(1))
25  end if}
26  case (cylindrical)
27  if (poleb(1,1).and.jg1==0) then ! cylindrical axis
28  jg1=1
29  {if (^d==phi_) jg^d=1+modulo(jg^d+ng^d(1)/2-1,ng^d(1))\}
30  end if
31  end select
32 
33  if (jg^d>=1.and.jg^d<=ng^d(1)|.and.) then
34  tree_neighbor%node => tree_root(jg^d)%node
35  else
36  nullify(tree_neighbor%node)
37  end if
38 
39 end subroutine find_root_neighbor
40 
41 !> find neighors of all blocks
42 subroutine find_neighbor(my_neighbor,my_neighbor_type,tree,i^D,pole)
45 
46  type(tree_node_ptr) :: tree, my_neighbor
47  integer, intent(in) :: i^D
48  integer, intent(out) :: my_neighbor_type
49  logical, dimension(ndim), intent(out) :: pole
50 
51  integer :: level, ig^D, ic^D, n_ic^D, inp^D
52 
53  pole=.false.
54  level=tree%node%level
55  if (level==1) then
56  call find_root_neighbor(my_neighbor,tree,i^d)
57  if (associated(my_neighbor%node)) then
58  if (phi_ > 0) then
59  ig^d=tree%node%ig^d;
60  {if ((poleb(2,^d).and.ig^d==ng^d(1).and.i^d==1) .or. &
61  (poleb(1,^d).and.ig^d==1.and.i^d==-1)) pole(^d)=.true.\}
62  end if
63  if (my_neighbor%node%leaf) then
64  my_neighbor_type=3
65  else
66  my_neighbor_type=4
67  end if
68  else
69  my_neighbor_type=1
70  return
71  end if
72  else
73  ig^d=tree%node%ig^d;
74 
75  if (phi_ > 0) then
76  {if ((poleb(2,^d).and.ig^d==ng^d(level).and.i^d==1) .or. &
77  (poleb(1,^d).and.ig^d==1.and.i^d==-1)) pole(^d)=.true.\}
78  end if
79 
80  ! ic^D is 1 when ig^D is odd, is 2 when ig^D is even
81  ic^d=1+modulo(ig^d-1,2);
82  inp^d=int((ic^d+i^d+1)/2)-1;
83  my_neighbor%node => tree%node%parent%node
84  {if (inp^d/=0) then
85  my_neighbor%node => my_neighbor%node%neighbor(ic^d,^d)%node
86  if (.not.associated(my_neighbor%node)) then
87  my_neighbor_type=1
88  return
89  end if
90  end if\}
91  if (my_neighbor%node%leaf) then
92  my_neighbor_type=2
93  else
94  {if (i^d==0 .or. pole(^d)) then
95  n_ic^d=ic^d
96  else
97  n_ic^d=3-ic^d ! switch 1 <--> 2
98  end if\}
99  my_neighbor%node => my_neighbor%node%child(n_ic^d)%node
100  if (associated(my_neighbor%node)) then
101  if (my_neighbor%node%leaf) then
102  my_neighbor_type=3
103  else
104  my_neighbor_type=4
105  end if
106  else
107  my_neighbor_type=0
108  end if
109  end if
110  end if
111 
112 end subroutine find_neighbor
113 
114 !> asign tree node neighor
115 subroutine asign_tree_neighbor(tree)
118 
119  type(tree_node_ptr) :: tree
120 
121  logical, dimension(ndim) :: pole
122  integer :: my_neighbor_type, i^D, iside
123  type(tree_node_ptr) :: my_neighbor
124 
125  {do iside=1,2
126  i^dd=kr(^dd,^d)*(2*iside-3);
127  call find_neighbor(my_neighbor,my_neighbor_type,tree,i^dd,pole)
128  select case (my_neighbor_type)
129  case (neighbor_sibling, neighbor_fine)
130  tree%node%neighbor(iside,^d)%node => my_neighbor%node
131  if (associated(my_neighbor%node)) then
132  if (pole(^d)) then
133  my_neighbor%node%neighbor(iside,^d)%node => tree%node
134  else
135  my_neighbor%node%neighbor(3-iside,^d)%node => tree%node
136  end if
137  end if
138  case default
139  nullify(tree%node%neighbor(iside,^d)%node)
140  end select
141  end do\}
142 
143 end subroutine asign_tree_neighbor
This module contains definitions of global parameters and variables and some generic functions/subrou...
Module with geometry-related routines (e.g., divergence, curl)
Definition: mod_geometry.t:2
integer coordinate
Definition: mod_geometry.t:6
Module with basic grid data structures.
Definition: mod_forest.t:2
subroutine find_neighbor(my_neighbor, my_neighbor_type, tree, iD, pole)
find neighors of all blocks
Definition: amr_neighbors.t:43
subroutine asign_tree_neighbor(tree)
asign tree node neighor
Pointer to a tree_node.
Definition: mod_forest.t:6
integer, dimension(:), allocatable ng
number of grid blocks in domain per dimension, in array over levels
logical, dimension(2, ndim) poleb
Indicates whether there is a pole at a boundary.
integer, dimension(3, 3) kr
Kronecker delta tensor.
logical, dimension(ndim) periodb
True for dimensions with periodic boundaries.
integer, parameter spherical
Definition: mod_geometry.t:10
subroutine find_root_neighbor(tree_neighbor, tree, iD)
find neighors of level-one root blocks
Definition: amr_neighbors.t:3