MPI-AMRVAC  2.2
The MPI - Adaptive Mesh Refinement - Versatile Advection Code
mod_multigrid_coupling.t
Go to the documentation of this file.
1 !> Module to couple the octree-mg library to AMRVAC. This file uses the VACPP
2 !> preprocessor, but its use is kept to a minimum.
4  {^ifoned
6  }
7  {^iftwod
9  }
10  {^ifthreed
11  use m_octree_mg_3d
12  }
13 
14 
15  implicit none
16  public
17 
18  !> Data structure containing the multigrid tree.
19  type(mg_t) :: mg
20 
21  !> If defined, this routine is called after a new multigrid tree is
22  !> constructed.
23  procedure(after_new_tree), pointer :: mg_after_new_tree => null()
24 
25  interface
26  subroutine after_new_tree()
27  end subroutine after_new_tree
28  end interface
29 
30 contains
31 
32  !> Setup multigrid for usage
33  subroutine mg_setup_multigrid()
35  use mod_geometry
36 
37  if (ndim /= mg_ndim) &
38  error stop "Multigrid module was compiled for different ndim"
39 
40  select case (coordinate)
41  case (cartesian)
42  continue
43  case (cylindrical)
44  if (ndim == 3) error stop "Multigrid does not support cylindrical 3D"
45  mg%geometry_type = mg_cylindrical
46  case default
47  error stop "Multigrid does not support your geometry"
48  end select
49 
50  if (any([ block_nx^d ] /= block_nx1)) &
51  error stop "Multigrid requires all block_nx to be equal"
52 
53  call mg_comm_init(mg)
54  call mg_set_methods(mg)
56  end subroutine mg_setup_multigrid
57 
58  !> Set multigrid boundary conditions for the solution according to variable iw
59  subroutine mg_copy_boundary_conditions(mg, iw)
61  type(mg_t), intent(inout) :: mg
62  integer, intent(in) :: iw
63  character(len=std_len) :: bnd_name(mg_num_neighbors)
64  integer :: n
65 
66  do n = 1, mg_num_neighbors
67  select case (typeboundary(iw, n))
68  case ('symm')
69  mg%bc(n, mg_iphi)%bc_type = mg_bc_neumann
70  mg%bc(n, mg_iphi)%bc_value = 0.0_dp
71  case ('asymm')
72  mg%bc(n, mg_iphi)%bc_type = mg_bc_dirichlet
73  mg%bc(n, mg_iphi)%bc_value = 0.0_dp
74  case ('cont')
75  mg%bc(n, mg_iphi)%bc_type = mg_bc_continuous
76  mg%bc(n, mg_iphi)%bc_value = 0.0_dp ! Not needed
77  case ('periodic')
78  ! Nothing to do here
79  case default
80  print *, "Not a standard: ", trim(typeboundary(iw, n))
81  error stop "You have to set a user-defined boundary method"
82  end select
83  end do
84  end subroutine mg_copy_boundary_conditions
85 
86  !> If the grid has changed, rebuild the full multigrid tree
87  subroutine mg_update_refinement(n_coarsen, n_refine)
89  integer, intent(in) :: n_coarsen
90  integer, intent(in) :: n_refine
91 
92  ! Don't build multigrid tree while doing initial refinement
93  if (.not. time_advance) return
94 
95  if (.not. mg%is_allocated) then
97  else if (n_coarsen + n_refine > 0) then
100  end if
101  end subroutine mg_update_refinement
102 
103  !> Copy a variable to the multigrid tree, including a layer of ghost cells
104  subroutine mg_copy_to_tree(iw_from, iw_to, restrict, restrict_gc, factor, state_from)
106  use mod_forest
107  integer, intent(in) :: iw_from !< Variable to use as right-hand side
108  integer, intent(in) :: iw_to !< Copy to this variable
109  logical, intent(in), optional :: restrict !< Restrict variable on multigrid tree
110  logical, intent(in), optional :: restrict_gc !< Fill ghost cells after restrict
111  real(dp), intent(in), optional :: factor !< out = factor * in
112  !> If present, use this as the input state
113  type(state), intent(in), optional, target :: state_from(max_blocks)
114  integer :: iigrid, igrid, id
115  integer :: nc, lvl, ilo^D, ihi^D
116  type(tree_node), pointer :: pnode
117  real(dp) :: fac
118  logical :: do_restrict, do_gc
119 
120  if (.not. mg%is_allocated) &
121  error stop "mg_copy_to_tree: tree not allocated yet"
122 
123  fac = 1.0_dp; if (present(factor)) fac = factor
124  do_restrict = .false.; if (present(restrict)) do_restrict = restrict
125  do_gc = .false.; if (present(restrict_gc)) do_gc = restrict_gc
126 
127  ilo^d=ixmlo^d-1;
128  ihi^d=ixmhi^d+1;
129 
130  do iigrid = 1, igridstail
131  igrid = igrids(iigrid);
132  pnode => igrid_to_node(igrid, mype)%node
133  id = pnode%id
134  lvl = mg%boxes(id)%lvl
135  nc = mg%box_size_lvl(lvl)
136 
137  ! Include one layer of ghost cells on grid leaves
138  if (present(state_from)) then
139  mg%boxes(id)%cc({0:nc+1}, iw_to) = &
140  fac * state_from(igrid)%w({ilo^d:ihi^d}, iw_from)
141  else
142  mg%boxes(id)%cc({0:nc+1}, iw_to) = fac * ps(igrid)%w({ilo^d:ihi^d}, iw_from)
143  end if
144  end do
145 
146  if (do_restrict) then
147  call mg_restrict(mg, iw_to)
148  if (do_gc) call mg_fill_ghost_cells(mg, iw_to)
149  end if
150 
151  end subroutine mg_copy_to_tree
152 
153  !> Copy a variable from the multigrid tree
154  subroutine mg_copy_from_tree(iw_from, iw_to, state_to)
156  use mod_forest
157  integer, intent(in) :: iw_from !< Variable to use as right-hand side
158  integer, intent(in) :: iw_to !< Copy to this variable
159  !> If present, use this as the output state
160  type(state), intent(inout), optional, target :: state_to(max_blocks)
161  integer :: iigrid, igrid, id
162  integer :: nc, lvl
163  type(tree_node), pointer :: pnode
164 
165  if (.not. mg%is_allocated) &
166  error stop "mg_copy_from_tree: tree not allocated yet"
167 
168  do iigrid = 1, igridstail
169  igrid = igrids(iigrid);
170  pnode => igrid_to_node(igrid, mype)%node
171  id = pnode%id
172  lvl = mg%boxes(id)%lvl
173  nc = mg%box_size_lvl(lvl)
174 
175  if (present(state_to)) then
176  state_to(igrid)%w(ixm^t, iw_to) = mg%boxes(id)%cc({1:nc}, iw_from)
177  else
178  ps(igrid)%w(ixm^t, iw_to) = mg%boxes(id)%cc({1:nc}, iw_from)
179  end if
180  end do
181  end subroutine mg_copy_from_tree
182 
183  !> Copy from multigrid tree with one layer of ghost cells. Corner ghost cells
184  !> are not used/set.
185  subroutine mg_copy_from_tree_gc(iw_from, iw_to, state_to)
187  use mod_forest
188  integer, intent(in) :: iw_from !< Variable to use as right-hand side
189  integer, intent(in) :: iw_to !< Copy to this variable
190  !> If present, use this as the output state
191  type(state), intent(inout), optional, target :: state_to(max_blocks)
192  integer :: iigrid, igrid, id
193  integer :: nc, lvl, ilo^D, ihi^D
194  type(tree_node), pointer :: pnode
195 
196  if (.not. mg%is_allocated) &
197  error stop "mg_copy_from_tree_gc: tree not allocated yet"
198 
199  ilo^d=ixmlo^d-1;
200  ihi^d=ixmhi^d+1;
201 
202  do iigrid = 1, igridstail
203  igrid = igrids(iigrid);
204  pnode => igrid_to_node(igrid, mype)%node
205  id = pnode%id
206  lvl = mg%boxes(id)%lvl
207  nc = mg%box_size_lvl(lvl)
208 
209  if (present(state_to)) then
210  state_to(igrid)%w({ilo^d:ihi^d}, iw_to) = mg%boxes(id)%cc({0:nc+1}, iw_from)
211  else
212  ps(igrid)%w({ilo^d:ihi^d}, iw_to) = mg%boxes(id)%cc({0:nc+1}, iw_from)
213  end if
214  end do
215  end subroutine mg_copy_from_tree_gc
216 
217  !> Generate a multigrid tree that includes the amrvac tree, but also contains
218  !> coarser grid levels. A number of checks has already been performed in
219  !> mg_setup_multigrid, so we don't repeat these checks here.
220  subroutine mg_tree_from_amrvac(mg)
223  type(mg_t), intent(inout) :: mg
224  integer :: i, n, id, ix(ndim)
225  integer :: n_boxes_total, i_c, c_id, c_ix(ndim)
226  integer :: min_lvl, lvl
227  integer :: nb, nb_ix, nb_dim
228  integer :: n_finer
229  type(tree_node), pointer :: pnode, pnode_ch
230  type(tree_node_ptr), allocatable :: id_to_node(:)
231  real(dp) :: dr_coarse
232 
233  ! Estimate number of finer blocks
234  n_finer = nparents+nleafs
235 
236  call mg_build_rectangle(mg, [ domain_nx^d ], block_nx1, dx(:,1), &
237  [ xprobmin^d ], periodb, n_finer)
238 
239  mg%highest_lvl = levmax
240  n_boxes_total = mg%n_boxes + n_finer
241 
242  ! To link the two trees
243  allocate(id_to_node(n_boxes_total))
244 
245  ! Link base level
246  do i = 1, size(mg%lvls(1)%ids)
247  id = mg%lvls(1)%ids(i)
248  ix = mg%boxes(id)%ix
249 
250  pnode => tree_root({ix(^d)})%node
251  pnode%id = id
252  id_to_node(id)%node => pnode
253  mg%boxes(id)%rank = pnode%ipe
254  end do
255 
256  ! Add refinement
257  do lvl = 1, mg%highest_lvl
258  do i = 1, size(mg%lvls(lvl)%ids)
259  id = mg%lvls(lvl)%ids(i)
260  pnode => id_to_node(id)%node
261 
262  if (.not. pnode%leaf) then
263  call mg_add_children(mg, id)
264 
265  do i_c = 1, mg_num_children
266  c_id = mg%boxes(id)%children(i_c)
267  c_ix = mg_child_dix(:, i_c) + 1
268  pnode_ch => pnode%child({c_ix(^d)})%node
269  id_to_node(c_id)%node => pnode_ch
270  pnode_ch%id = c_id
271  mg%boxes(c_id)%rank = pnode_ch%ipe
272  end do
273  end if
274  end do
275 
276  call mg_set_leaves_parents(mg%boxes, mg%lvls(lvl))
277 
278  if (lvl < mg%highest_lvl) then
279  call mg_set_next_level_ids(mg, lvl)
280  call mg_set_neighbors_lvl(mg, lvl+1)
281  end if
282  end do
283 
284  ! Store boxes with refinement boundaries (from the coarse side)
285  do lvl = 1, mg%highest_lvl
286  call mg_set_refinement_boundaries(mg%boxes, mg%lvls(lvl))
287  end do
288 
289  ! Assign boxes to MPI processes
290  call mg_load_balance_parents(mg)
291 
292  ! Allocate storage for boxes owned by this process
293  call mg_allocate_storage(mg)
294 
295  if (associated(mg_after_new_tree)) then
296  call mg_after_new_tree()
297  end if
298 
299  end subroutine mg_tree_from_amrvac
300 
301 end module mod_multigrid_coupling
302 
integer, parameter cylindrical
Definition: mod_geometry.t:9
subroutine mg_update_refinement(n_coarsen, n_refine)
If the grid has changed, rebuild the full multigrid tree.
This module contains definitions of global parameters and variables and some generic functions/subrou...
subroutine, public mg_deallocate_storage(mg)
Deallocate all allocatable arrays.
integer domain_nx
number of cells for each dimension in level-one mesh
integer max_blocks
The maximum number of grid blocks in a processor.
The data structure that contains information about a tree node/grid block.
Definition: mod_forest.t:11
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 mg_tree_from_amrvac(mg)
Generate a multigrid tree that includes the amrvac tree, but also contains coarser grid levels...
Module with geometry-related routines (e.g., divergence, curl)
Definition: mod_geometry.t:2
integer, dimension(1, 2), parameter, public mg_child_dix
integer coordinate
Definition: mod_geometry.t:6
Module with basic grid data structures.
Definition: mod_forest.t:2
integer, parameter, public mg_num_children
type(mg_t) mg
Data structure containing the multigrid tree.
integer, parameter, public mg_bc_continuous
Value to indicate a continuous boundary condition.
Module to couple the octree-mg library to AMRVAC. This file uses the VACPP preprocessor, but its use is kept to a minimum.
Pointer to a tree_node.
Definition: mod_forest.t:6
subroutine, public mg_restrict(mg, iv)
Restrict all levels.
integer, parameter, public mg_num_neighbors
subroutine mg_copy_from_tree(iw_from, iw_to, state_to)
Copy a variable from the multigrid tree.
subroutine, public mg_set_methods(mg)
subroutine mg_copy_boundary_conditions(mg, iw)
Set multigrid boundary conditions for the solution according to variable iw.
subroutine, public mg_add_children(mg, id)
subroutine, public mg_comm_init(mg, comm)
Prolong from a parent to a child with index offset dix. This method could sometimes work better than ...
type(tree_node_ptr), dimension(:^d &), allocatable, save tree_root
Pointers to the coarse grid.
Definition: mod_forest.t:29
integer, parameter cartesian
Definition: mod_geometry.t:7
integer, parameter, public mg_bc_neumann
Value to indicate a Neumann boundary condition.
integer ixm
the mesh range (within a block with ghost cells)
subroutine, public mg_set_leaves_parents(boxes, level)
Create a list of leaves and a list of parents for a level.
character(len=std_len), dimension(:, :), allocatable typeboundary
Array indicating the type of boundary condition per variable and per physical boundary.
double precision, dimension(:,:), allocatable dx
integer, parameter, public mg_bc_dirichlet
Value to indicate a Dirichlet boundary condition.
subroutine, public mg_set_neighbors_lvl(mg, lvl)
integer, dimension(:), allocatable, parameter d
integer, parameter, public mg_ndim
Problem dimension.
integer mype
The rank of the current MPI task.
procedure(after_new_tree), pointer mg_after_new_tree
If defined, this routine is called after a new multigrid tree is constructed.
integer block_nx
number of cells for each dimension in grid block excluding ghostcells
subroutine, public mg_load_balance_parents(mg)
Load balance the parents (non-leafs). Assign them to the rank that has most children.
integer, parameter, public mg_iphi
Index of solution.
subroutine mg_copy_from_tree_gc(iw_from, iw_to, state_to)
Copy from multigrid tree with one layer of ghost cells. Corner ghost cells are not used/set...
integer, parameter, public mg_cylindrical
Cylindrical coordinate system.
subroutine mg_copy_to_tree(iw_from, iw_to, restrict, restrict_gc, factor, state_from)
Copy a variable to the multigrid tree, including a layer of ghost cells.
subroutine, public mg_build_rectangle(mg, domain_size, box_size, dx, r_min, periodic, n_finer)
subroutine mg_setup_multigrid()
Setup multigrid for usage.
subroutine, public mg_allocate_storage(mg)
Allocate communication buffers and local boxes for a tree that has already been created.
subroutine, public mg_set_next_level_ids(mg, lvl)
integer, parameter ndim
Number of spatial dimensions for grid variables.
subroutine, public mg_set_refinement_boundaries(boxes, level)
Create a list of refinement boundaries (from the coarse side)
subroutine, public mg_fill_ghost_cells(mg, iv)
Fill ghost cells at all grid levels.
integer, save nparents
Number of parent blocks.
Definition: mod_forest.t:73
logical, dimension(ndim) periodb
True for dimensions with periodic boundaries.
integer, save nleafs
Number of leaf block.
Definition: mod_forest.t:76