MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
forest.t
Go to the documentation of this file.
1 !> build root forest
2 subroutine init_forest_root
6 
7  integer :: ig^D, level, igrid, ipe
8  integer :: iside, i^D, Morton_no, isfc
9 
10  integer, external :: getnode
11 
12  level=1
13  morton_no=0
14  ipe=0
15  nparents=0
16  nleafs={ng^d(1)*}
18  nleafs_level(1)={ng^d(1)*}
20  call get_morton_range
21  ! Generate Morton-order space-filling curve to connect level 1 blocks
23  do isfc=1,nglev1
24  ig^d=sfc_iglevel1(^d,isfc)\
25  morton_no=morton_no+1
26  if (morton_no>morton_stop(ipe)) ipe=ipe+1
27  igrid=getnode(ipe)
28  if (ipe==mype) then
29  sfc_to_igrid(morton_no)=igrid
30  igrid_to_sfc(igrid)=morton_no
31  end if
32  call init_tree_leaf(tree_root(ig^d),ig^d,level,igrid,ipe,.true.)
33  end do
34 
35  ! update root neighbor
36  {do ig^db=1,ng^db(1)\}
37  {do iside=1,2
38  i^dd=kr(^dd,^d)*(2*iside-3);
39  call find_root_neighbor(tree_root(ig^dd)%node%neighbor(iside,^d), &
40  tree_root(ig^dd),i^dd)
41  end do\}
42  {end do\}
43 
44  ! This call is here to ensure the sfc array is initialized
45  call amr_morton_order()
46 
47 end subroutine init_forest_root
48 
49 subroutine init_tree_leaf(tree,ig^D,level,igrid,ipe,active)
51  implicit none
52 
53  type(tree_node_ptr) :: tree
54  integer, intent(in) :: ig^D, level, igrid, ipe
55  logical, intent(in) :: active
56  integer :: ic^D
57 
58  allocate(tree%node)
59 
60  tree%node%ig^d=ig^d;
61  tree%node%level=level
62  tree%node%igrid=igrid
63  tree%node%ipe=ipe
64 
65  tree%node%leaf=.true.
66  tree%node%active=active
67 
68  nullify(tree%node%parent%node)
69  {do ic^db=1,2\}
70  nullify(tree%node%child(ic^d)%node)
71  {end do\}
72 
73  call add_to_linked_list(level,tree)
74 
75  ! initialize neighbor pointers
76  nullify({tree%node%neighbor(1,^d)%node},{tree%node%neighbor(2,^d)%node})
77 
78  igrid_to_node(igrid,ipe)%node => tree%node
79 
80 end subroutine init_tree_leaf
81 
82 subroutine coarsen_tree_leaf(igrid,ipe,child_igrid,child_ipe,active)
84  implicit none
85 
86  integer, intent(in) :: igrid, ipe
87  integer, dimension(2^D&), intent(in) :: child_igrid, child_ipe
88  logical, intent(out) :: active
89 
90 
91  integer :: level, ic^D, child_level, iside, iotherside, vote
92  type(tree_node_ptr) :: tree, child, child_neighbor
93 
94  tree%node => igrid_to_node(child_igrid(1^d&),child_ipe(1^d&))%node%parent%node
95  level=tree%node%level
96 
97  call add_to_linked_list(level,tree)
98 
99  child_level=level+1
100  vote=0
101 
102  {do ic^db=1,2\}
103  child%node => tree%node%child(ic^d)%node
104 
105  ! vote for active:
106  if(child%node%active) vote=vote+1
107 
108  call delete_from_linked_list(child_level,child)
109 
110  ! update neighbor pointers
111  {iside=ic^d
112  child_neighbor%node => child%node%neighbor(iside,^d)%node
113  if (associated(child_neighbor%node)) then
114  if (child%node%ig^d==child_neighbor%node%ig^d) then ! pole
115  nullify(child_neighbor%node%neighbor(iside,^d)%node)
116  else
117  iotherside=3-iside
118  nullify(child_neighbor%node%neighbor(iotherside,^d)%node)
119  end if
120  end if\}
121 
122  nullify(tree%node%child(ic^d)%node)
123  deallocate(igrid_to_node(child_igrid(ic^d),child_ipe(ic^d))%node)
124  {end do\}
125 
126  tree%node%leaf=.true.
127  tree%node%igrid=igrid
128  tree%node%ipe=ipe
129  igrid_to_node(igrid,ipe)%node => tree%node
130 
131  ! Count the vote and set active/passive state:
132 
133  if (vote /= 2**^nd) then
134  !if (vote == 0) then
135  tree%node%active = .false.
136  nleafs_active = nleafs_active - vote
137  else
138  tree%node%active = .true.
139  nleafs_active = nleafs_active - vote + 1
140  end if
141  active = tree%node%active
142 
143  nleafs=nleafs-2**^nd+1
144  nparents=nparents-1
145  nleafs_level(child_level)=nleafs_level(child_level)-2**^nd
146  nleafs_level(level)=nleafs_level(level)+1
147 
148 end subroutine coarsen_tree_leaf
149 
150 subroutine refine_tree_leaf(child_igrid,child_ipe,igrid,ipe,active)
153 
154  integer, dimension(2^D&), intent(in) :: child_igrid, child_ipe
155  integer, intent(in) :: igrid, ipe
156  logical, intent(out):: active
157 
158  integer :: ig^D, level, i^D, ic^D, child_ig^D, child_level, iside
159  integer :: my_neighbor_type
160  logical, dimension(ndim) :: pole
161  type(tree_node_ptr) :: tree, child, my_neighbor
162 
163  tree%node => igrid_to_node(igrid,ipe)%node
164  ig^d=tree%node%ig^d;
165  level=tree%node%level
166  active=tree%node%active
167 
168  tree%node%ipe=-1
169  tree%node%igrid=0
170  tree%node%leaf=.false.
171  tree%node%active=.true.
172 
173  call delete_from_linked_list(level,tree)
174 
175  child_level=level+1
176 
177  {do ic^db=1,2\}
178  child_ig^d=2*(ig^d-1)+ic^d;
179  call init_tree_leaf(child,child_ig^d,child_level, &
180  child_igrid(ic^d),child_ipe(ic^d),active)
181 
182  igrid_to_node(child_igrid(ic^d),child_ipe(ic^d))%node => child%node
183 
184  tree%node%child(ic^d)%node => child%node
185  child%node%parent%node => tree%node
186  {end do\}
187 
188  ! update neighbor pointers
189  {do ic^db=1,2\}
190  child%node => tree%node%child(ic^d)%node
191  {iside=ic^d
192  i^dd=kr(^dd,^d)*(2*iside-3);
193  call find_neighbor(my_neighbor,my_neighbor_type,child,i^dd,pole)
194  select case (my_neighbor_type)
195  case (neighbor_sibling, neighbor_fine)
196  child%node%neighbor(iside,^d)%node => my_neighbor%node
197  if (pole(^d)) then
198  my_neighbor%node%neighbor(iside,^d)%node => child%node
199  else
200  my_neighbor%node%neighbor(3-iside,^d)%node => child%node
201  end if
202  case default
203  nullify(child%node%neighbor(iside,^d)%node)
204  end select
205  child%node%neighbor(3-ic^d,^d)%node=>tree%node%child(3-ic^d^d%ic^dd)%node\}
206  {end do\}
207 
208  nleafs=nleafs+2**^nd-1
210  nleafs_level(child_level)=nleafs_level(child_level)+2**^nd
211  nleafs_level(level)=nleafs_level(level)-1
212 
213  if (active) nleafs_active = nleafs_active + 2**^nd-1
214 
215 end subroutine refine_tree_leaf
216 
217 subroutine change_ipe_tree_leaf(recv_igrid,recv_ipe,send_igrid,send_ipe)
219  implicit none
220 
221  integer, intent(in) :: recv_igrid, recv_ipe, send_igrid, send_ipe
222 
223  type(tree_node_ptr) :: tree
224 
225  tree%node => igrid_to_node(send_igrid,send_ipe)%node
226 
227  tree%node%igrid=recv_igrid
228  tree%node%ipe=recv_ipe
229 
230  nullify(igrid_to_node(send_igrid,send_ipe)%node)
231  igrid_to_node(recv_igrid,recv_ipe)%node => tree%node
232 
233 end subroutine change_ipe_tree_leaf
234 
235 subroutine add_to_linked_list(level,tree)
237  implicit none
238 
239  integer, intent(in) :: level
240  type(tree_node_ptr) :: tree
241 
242  nullify(tree%node%next%node)
243  if (associated(level_head(level)%node)) then
244  tree%node%prev%node => level_tail(level)%node
245  level_tail(level)%node%next%node => tree%node
246  level_tail(level)%node => tree%node
247  else
248  level_head(level)%node => tree%node
249  level_tail(level)%node => tree%node
250  nullify(tree%node%prev%node)
251  end if
252 
253 end subroutine add_to_linked_list
254 
255 subroutine delete_from_linked_list(level,tree)
257  implicit none
258 
259  integer, intent(in) :: level
260  type(tree_node_ptr) :: tree
261 
262  type(tree_node_ptr) :: next, prev
263 
264  prev%node => tree%node%prev%node
265  next%node => tree%node%next%node
266  if (associated(next%node).and.associated(prev%node)) then
267  prev%node%next%node => next%node
268  next%node%prev%node => prev%node
269  else if (associated(prev%node)) then
270  level_tail(level)%node => prev%node
271  nullify(prev%node%next%node)
272  else if (associated(next%node)) then
273  level_head(level)%node => next%node
274  nullify(next%node%prev%node)
275  else
276  nullify(level_head(level)%node)
277  nullify(level_tail(level)%node)
278  end if
279 
280 end subroutine delete_from_linked_list
281 
282 subroutine write_forest(file_handle)
285 
286  integer, intent(in) :: file_handle
287 
288  integer, dimension(MPI_STATUS_SIZE) :: status
289  integer :: ig^D,isfc
290 
291  do isfc=1,nglev1
292  ig^d=sfc_iglevel1(^d,isfc)\
293  call write_node(tree_root(ig^d))
294  end do
295 
296  contains
297 
298  recursive subroutine write_node(tree)
299  implicit none
300 
301  type(tree_node_ptr) :: tree
302 
303  integer :: ic^D
304 
305  call mpi_file_write(file_handle,tree%node%leaf,1,mpi_logical,status,ierrmpi)
306 
307  if (.not.tree%node%leaf) then
308  {do ic^db=1,2\}
309  call write_node(tree%node%child(ic^d))
310  {end do\}
311  end if
312 
313  end subroutine write_node
314 
315 end subroutine write_forest
316 
317 subroutine read_forest(file_handle)
321 
322  integer, intent(in) :: file_handle
323 
324  integer, dimension(MPI_STATUS_SIZE) :: status
325  !integer :: ig^D, level, size_logical, Morton_no, igrid, ipe
326  integer :: ig^D, level, Morton_no, igrid, ipe, isfc
327  integer, external :: getnode
328 
329  morton_no=0
330  ipe=0
331  level=1
332  nleafs_level(1:nlevelshi) = 0
333  nparents = 0
334 
335  call get_morton_range
337  do isfc=1,nglev1
338  ig^d=sfc_iglevel1(^d,isfc)\
339  allocate(tree_root(ig^d)%node)
340  nullify(tree_root(ig^d)%node%parent%node)
341  call read_node(tree_root(ig^d),ig^d,level)
342  end do
343 
344  call get_level_range
345 
346  ! Rebuild tree connectivity
347  call getigrids
348  call build_connectivity
349 
350  ! This call is here to ensure the sfc array is initialized
351  call amr_morton_order()
352 
353  contains
354 
355  recursive subroutine read_node(tree,ig^D,level)
356  implicit none
357 
358  type(tree_node_ptr) :: tree
359  integer, intent(in) :: ig^D, level
360 
361  logical :: leaf
362  integer :: ic^D, child_ig^D, child_level
363 
364  if (mype==0) then
365  call mpi_file_read(file_handle,leaf,1,mpi_logical, &
366  status,ierrmpi)
367  end if
368  if (npe>1) call mpi_bcast(leaf,1,mpi_logical,0,icomm,ierrmpi)
369 
370  tree%node%leaf=leaf
371  tree%node%ig^d=ig^d;
372  tree%node%level=level
373  tree%node%active=.true. .and. leaf
374 
375  {do ic^db=1,2\}
376  nullify(tree%node%child(ic^d)%node)
377  {end do\}
378  nullify({tree%node%neighbor(1,^d)%node},{tree%node%neighbor(2,^d)%node})
379  nullify(tree%node%next%node,tree%node%prev%node)
380 
381  call asign_tree_neighbor(tree)
382 
383  if (leaf) then
384  call add_to_linked_list(level,tree)
385  nleafs_level(level) = nleafs_level(level) + 1
386 
387  morton_no=morton_no+1
388  if (morton_no>morton_stop(ipe)) ipe=ipe+1
389  igrid=getnode(ipe)
390  tree%node%igrid=igrid
391  tree%node%ipe=ipe
392  igrid_to_node(igrid,ipe)%node => tree%node
393  if (ipe==mype) sfc_to_igrid(morton_no)=igrid
394  else
395  nparents = nparents + 1
396  tree%node%igrid=0
397  tree%node%ipe=-1
398  child_level=level+1
399  {do ic^db=1,2\}
400  child_ig^d=2*(ig^d-1)+ic^d;
401  allocate(tree%node%child(ic^d)%node)
402  tree%node%child(ic^d)%node%parent%node => tree%node
403  call read_node(tree%node%child(ic^d),child_ig^d,child_level)
404  {end do\}
405  end if
406 
407  end subroutine read_node
408 
409 end subroutine read_forest
This module contains definitions of global parameters and variables and some generic functions/subrou...
subroutine get_morton_range
Set the Morton range for each processor.
type(tree_node_ptr), dimension(:,:), allocatable, save igrid_to_node
Array to go from an [igrid, ipe] index to a node pointer.
Definition: mod_forest.t:32
subroutine init_tree_leaf(tree, igD, level, igrid, ipe, active)
Definition: forest.t:50
integer nglev1
Definition: mod_forest.t:78
Module with basic grid data structures.
Definition: mod_forest.t:2
subroutine build_connectivity
Definition: connectivity.t:45
subroutine delete_from_linked_list(level, tree)
Definition: forest.t:256
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
integer npe
The number of MPI tasks.
subroutine init_forest_root
build root forest
Definition: forest.t:3
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
subroutine write_forest(file_handle)
Definition: forest.t:283
integer, dimension(:), allocatable, save sfc_to_igrid
Go from a Morton number to an igrid index (for a single processor)
Definition: mod_forest.t:53
integer, dimension(:,:), allocatable, save sfc_iglevel1
Space filling curve for level 1 grid. sfc_iglevel1(^D, MN) gives ig^D (the spatial index of the grid)...
Definition: mod_forest.t:47
subroutine coarsen_tree_leaf(igrid, ipe, child_igrid, child_ipe, active)
Definition: forest.t:83
type(tree_node_ptr), dimension(:^d &), allocatable, save tree_root
Pointers to the coarse grid.
Definition: mod_forest.t:29
integer ierrmpi
A global MPI error return code.
type(tree_node_ptr), dimension(:), allocatable, save level_head
The head pointer of the linked list per refinement level.
Definition: mod_forest.t:35
recursive subroutine write_node(tree)
Definition: forest.t:299
integer, dimension(3, 3) kr
Kronecker delta tensor.
integer nleafs_active
Definition: mod_forest.t:78
recursive subroutine read_node(tree, igD, level)
Definition: forest.t:356
integer, dimension(:), allocatable, save morton_stop
Last Morton number per processor.
Definition: mod_forest.t:65
subroutine getigrids
Definition: connectivity.t:26
integer mype
The rank of the current MPI task.
subroutine read_forest(file_handle)
Definition: forest.t:318
integer, dimension(:), allocatable, save nleafs_level
How many leaves are present per refinement level.
Definition: mod_forest.t:81
subroutine level1_morton_order
build Morton space filling curve for level 1 grid
subroutine change_ipe_tree_leaf(recv_igrid, recv_ipe, send_igrid, send_ipe)
Definition: forest.t:218
subroutine refine_tree_leaf(child_igrid, child_ipe, igrid, ipe, active)
Definition: forest.t:151
type(tree_node_ptr), dimension(:), allocatable, save level_tail
The tail pointer of the linked list per refinement level.
Definition: mod_forest.t:38
subroutine get_level_range
Definition: connectivity.t:2
integer icomm
The MPI communicator.
subroutine amr_morton_order
Construct Morton-order as a global recursive lexicographic ordering.
integer, parameter nlevelshi
The maximum number of levels in the grid refinement.
integer, save nparents
Number of parent blocks.
Definition: mod_forest.t:73
subroutine add_to_linked_list(level, tree)
Definition: forest.t:236
integer, dimension(:), allocatable, save igrid_to_sfc
Go from a grid index to Morton number (for a single processor)
Definition: mod_forest.t:56
integer, save nleafs
Number of leaf block.
Definition: mod_forest.t:76
subroutine find_root_neighbor(tree_neighbor, tree, iD)
find neighors of level-one root blocks
Definition: amr_neighbors.t:3