Commit d1d2a590 authored by Jannis Teunissen's avatar Jannis Teunissen

Simplify m_af_restrict module

parent 2dd041c5
......@@ -13,8 +13,7 @@ PROGS_XD := random_refinement poisson_basic poisson_benchmark advection \
poisson_dielectric poisson_helmholtz reaction_diffusion dielectric_surface advection_v2
PROGS_2D := $(PROGS_XD:%=%_2d) poisson_cyl poisson_cyl_dielectric \
simple_streamer poisson_cyl_analytic poisson_helmholtz_cyl h_domain_expts h_bc \
h_advection solid_body_rotation amr_solid_body_rotation euler_gas_dynamics
simple_streamer poisson_cyl_analytic poisson_helmholtz_cyl solid_body_rotation amr_solid_body_rotation euler_gas_dynamics
PROGS_3D := $(PROGS_XD:%=%_3d) poisson_div_cleaning
......
......@@ -86,7 +86,7 @@ program advection
call af_print_info(tree)
! Restrict the initial conditions
call af_restrict_tree(tree, i_phi)
call af_restrict_tree(tree, [i_phi])
! Fill ghost cells for variables i_phi on the sides of all boxes
call af_gc_tree(tree, [i_phi])
......@@ -146,7 +146,7 @@ program advection
call af_loop_box_arg(tree, update_solution, [dt])
! Restrict variables i_phi to all parent boxes
call af_restrict_tree(tree, i_phi)
call af_restrict_tree(tree, [i_phi])
end do
! Take average of phi_old and phi
......
......@@ -73,7 +73,7 @@ program advection
call af_print_info(tree)
! Restrict the initial conditions
call af_restrict_tree(tree, i_phi)
call af_restrict_tree(tree, [i_phi])
! Fill ghost cells for variables i_phi on the sides of all boxes
call af_gc_tree(tree, [i_phi])
......
......@@ -44,7 +44,7 @@ program solid_body_rotation
call af_gc_tree(tree, [iq])
!One level of refinement
call af_adjust_refinement(tree, refine_routine, refine_info, 1)
call af_restrict_tree(tree, iq)
call af_restrict_tree(tree, [iq])
call af_gc_tree(tree, [iq])
!call af_write_silo(tree, "testsbr_amr", dir="output")
......@@ -70,7 +70,7 @@ program solid_body_rotation
call af_loop_box_arg(tree, updateSoln, [dt])
!print *, 'finished time stepping'
call af_restrict_tree(tree, iq)
call af_restrict_tree(tree, [iq])
call af_gc_tree(tree, [iq])
!print *, 'Updating ghost cells'
......
......@@ -117,10 +117,7 @@ program KT_euler
call af_adjust_refinement(tree, ref_rout, refine_info, 1)
if (refine_info%n_add == 0) exit
end do
call af_restrict_tree(tree, i_rho)
call af_restrict_tree(tree, i_mom(1))
call af_restrict_tree(tree, i_mom(2))
call af_restrict_tree(tree, i_e)
call af_restrict_tree(tree, variables)
call af_gc_tree(tree, variables)
else
call af_loop_box(tree, set_init_conds)
......
......@@ -62,7 +62,7 @@ program dielectric_test
! Average epsilon on coarse grids. In the future, it could be better to define
! epsilon on cell faces, and to perform this restriction in a matrix fashion:
! A_coarse = M_restrict * A_fine * M_prolong (A = matrix operator, M = matrix)
call af_restrict_tree(tree, i_eps)
call af_restrict_tree(tree, [i_eps])
write(*,"(A,Es10.3,A)") " Wall-clock time generating AMR grid: ", &
(t_end-t_start) / real(count_rate,dp), " seconds"
......
......@@ -77,8 +77,6 @@ program reaction_diffusion
call af_set_cc_methods(tree, i_u, af_bc_neumann_zero)
call af_set_cc_methods(tree, i_v, af_bc_neumann_zero)
call af_set_cc_methods(tree, i_u+1, af_bc_neumann_zero)
call af_set_cc_methods(tree, i_v+1, af_bc_neumann_zero)
! Initialize tree
call af_init(tree, & ! Tree to initialize
......
......@@ -727,7 +727,7 @@ contains
do i_ch = 1, af_num_children
ch_id = tree%boxes(id)%children(i_ch)
call tree%cc_methods(iv)%restrict(tree%boxes(ch_id), &
tree%boxes(id), iv)
tree%boxes(id), [iv])
end do
end if
end do
......
......@@ -1085,10 +1085,10 @@ contains
if (iv == mg%i_phi) then
! Don't use geometry for restriction, since this is inconsistent with the
! filling of ghost cells near refinement boundaries
call af_restrict_box(box_c, box_p, iv, use_geometry=.false.)
call af_restrict_box(box_c, box_p, [iv], use_geometry=.false.)
else
! For the right-hand side, use the geometry
call af_restrict_box(box_c, box_p, iv, use_geometry=.true.)
call af_restrict_box(box_c, box_p, [iv], use_geometry=.true.)
end if
end subroutine mg_box_rstr_lpl
......
......@@ -126,7 +126,7 @@ contains
error stop "af_particles_to_grid: Invalid interpolation order"
end select
call af_restrict_tree(tree, iv)
call af_restrict_tree(tree, [iv])
if (fill_gc_at_end) then
if (.not. tree%has_cc_method(iv)) then
......
......@@ -11,7 +11,6 @@ module m_af_restrict
public :: af_restrict_to_boxes
public :: af_restrict_tree
public :: af_restrict_box
public :: af_restrict_box_vars
public :: af_restrict_ref_boundary
! public :: af_restrict_box_face
......@@ -19,57 +18,53 @@ contains
!> Restrict the children of a box to the box (e.g., in 2D, average the values
!> at the four children to get the value for the parent)
subroutine af_restrict_to_box(boxes, id, iv, iv_to)
type(box_t), intent(inout) :: boxes(:) !< List of all the boxes
integer, intent(in) :: id !< Box whose children will be restricted to it
integer, intent(in) :: iv !< Variable to restrict
integer, intent(in), optional :: iv_to !< Destination (if /= iv)
integer :: nc, i_c, c_id
subroutine af_restrict_to_box(boxes, id, ivs)
type(box_t), intent(inout) :: boxes(:) !< List of all the boxes
integer, intent(in) :: id !< Box whose children will be restricted to it
integer, intent(in) :: ivs(:) !< Variable to restrict
integer :: nc, i_c, c_id
nc = boxes(id)%n_cell
do i_c = 1, af_num_children
c_id = boxes(id)%children(i_c)
if (c_id == af_no_box) cycle
call af_restrict_box(boxes(c_id), boxes(id), iv, iv_to)
call af_restrict_box(boxes(c_id), boxes(id), ivs)
end do
end subroutine af_restrict_to_box
!> Restrict the children of boxes ids(:) to them.
subroutine af_restrict_to_boxes(boxes, ids, iv, iv_to)
type(box_t), intent(inout) :: boxes(:) !< List of all the boxes
integer, intent(in) :: ids(:) !< Boxes whose children will be restricted to it
integer, intent(in) :: iv !< Variable to restrict
integer, intent(in), optional :: iv_to !< Destination (if /= iv)
integer :: i
subroutine af_restrict_to_boxes(boxes, ids, ivs)
type(box_t), intent(inout) :: boxes(:) !< List of all the boxes
integer, intent(in) :: ids(:) !< Boxes whose children will be restricted to it
integer, intent(in) :: ivs(:) !< Variables to restrict
integer :: i
!$omp parallel do
do i = 1, size(ids)
call af_restrict_to_box(boxes, ids(i), iv, iv_to)
call af_restrict_to_box(boxes, ids(i), ivs)
end do
!$omp end parallel do
end subroutine af_restrict_to_boxes
!> Restrict variables iv to all parent boxes, from the highest to the lowest level
subroutine af_restrict_tree(tree, iv, iv_to)
type(af_t), intent(inout) :: tree !< Tree to restrict on
integer, intent(in) :: iv !< Variable to restrict
integer, intent(in), optional :: iv_to !< Destination (if /= iv)
integer :: lvl
subroutine af_restrict_tree(tree, ivs)
type(af_t), intent(inout) :: tree !< Tree to restrict on
integer, intent(in) :: ivs(:) !< Variables to restrict
integer :: lvl
if (.not. tree%ready) stop "Tree not ready"
do lvl = tree%highest_lvl-1, 1, -1
call af_restrict_to_boxes(tree%boxes, tree%lvls(lvl)%parents, iv, iv_to)
call af_restrict_to_boxes(tree%boxes, tree%lvls(lvl)%parents, ivs)
end do
end subroutine af_restrict_tree
!> Restriction of child box (box_c) to its parent (box_p)
subroutine af_restrict_box(box_c, box_p, iv, iv_to, use_geometry)
type(box_t), intent(in) :: box_c !< Child box to restrict
type(box_t), intent(inout) :: box_p !< Parent box to restrict to
integer, intent(in) :: iv !< Variable to restrict
integer, intent(in), optional :: iv_to !< Destination (if /= iv)
subroutine af_restrict_box(box_c, box_p, ivs, use_geometry)
type(box_t), intent(in) :: box_c !< Child box to restrict
type(box_t), intent(inout) :: box_p !< Parent box to restrict to
integer, intent(in) :: ivs(:) !< Variable to restrict
logical, intent(in), optional :: use_geometry !< If set to false, don't use geometry
integer :: i, j, i_f, j_f, i_c, j_c, i_dest
integer :: i, j, i_f, j_f, i_c, j_c, n, iv
integer :: hnc, ix_offset(NDIM)
logical :: use_geom
#if NDIM == 2
......@@ -81,76 +76,59 @@ contains
hnc = ishft(box_c%n_cell, -1) ! n_cell / 2
ix_offset = af_get_child_offset(box_c)
if (present(iv_to)) then
i_dest = iv_to
else
i_dest = iv
end if
if (present(use_geometry)) then
use_geom = use_geometry
else
use_geom = .true.
end if
do n = 1, size(ivs)
iv = ivs(n)
#if NDIM == 2
if (box_p%coord_t == af_cyl .and. use_geom) then
do j = 1, hnc
j_c = ix_offset(2) + j
j_f = 2 * j - 1
do i = 1, hnc
i_c = ix_offset(1) + i
i_f = 2 * i - 1
call af_cyl_child_weights(box_p, i_c, w1, w2)
box_p%cc(i_c, j_c, i_dest) = 0.25_dp * (&
w1 * sum(box_c%cc(i_f, j_f:j_f+1, iv)) + &
w2 * sum(box_c%cc(i_f+1, j_f:j_f+1, iv)))
if (box_p%coord_t == af_cyl .and. use_geom) then
do j = 1, hnc
j_c = ix_offset(2) + j
j_f = 2 * j - 1
do i = 1, hnc
i_c = ix_offset(1) + i
i_f = 2 * i - 1
call af_cyl_child_weights(box_p, i_c, w1, w2)
box_p%cc(i_c, j_c, iv) = 0.25_dp * (&
w1 * sum(box_c%cc(i_f, j_f:j_f+1, iv)) + &
w2 * sum(box_c%cc(i_f+1, j_f:j_f+1, iv)))
end do
end do
end do
else
do j = 1, hnc
j_c = ix_offset(2) + j
j_f = 2 * j - 1
do i = 1, hnc
i_c = ix_offset(1) + i
i_f = 2 * i - 1
box_p%cc(i_c, j_c, i_dest) = 0.25_dp * &
sum(box_c%cc(i_f:i_f+1, j_f:j_f+1, iv))
else
do j = 1, hnc
j_c = ix_offset(2) + j
j_f = 2 * j - 1
do i = 1, hnc
i_c = ix_offset(1) + i
i_f = 2 * i - 1
box_p%cc(i_c, j_c, iv) = 0.25_dp * &
sum(box_c%cc(i_f:i_f+1, j_f:j_f+1, iv))
end do
end do
end do
endif
endif
#elif NDIM == 3
do k = 1, hnc
k_c = ix_offset(3) + k
k_f = 2 * k - 1
do j = 1, hnc
j_c = ix_offset(2) + j
j_f = 2 * j - 1
do i = 1, hnc
i_c = ix_offset(1) + i
i_f = 2 * i - 1
box_p%cc(i_c, j_c, k_c, i_dest) = 0.125_dp * &
sum(box_c%cc(i_f:i_f+1, j_f:j_f+1, k_f:k_f+1, iv))
do k = 1, hnc
k_c = ix_offset(3) + k
k_f = 2 * k - 1
do j = 1, hnc
j_c = ix_offset(2) + j
j_f = 2 * j - 1
do i = 1, hnc
i_c = ix_offset(1) + i
i_f = 2 * i - 1
box_p%cc(i_c, j_c, k_c, iv) = 0.125_dp * &
sum(box_c%cc(i_f:i_f+1, j_f:j_f+1, k_f:k_f+1, iv))
end do
end do
end do
end do
#endif
end subroutine af_restrict_box
!> Restriction of child box (box_c) to its parent (box_p)
subroutine af_restrict_box_vars(box_c, box_p, ivs, use_geometry)
type(box_t), intent(in) :: box_c !< Child box to restrict
type(box_t), intent(inout) :: box_p !< Parent box to restrict to
integer, intent(in) :: ivs(:) !< Variables to restrict
logical, intent(in), optional :: use_geometry !< If set to false, don't use geometry
integer :: i, iv
do i = 1, size(ivs)
iv = ivs(i)
call af_restrict_box(box_c, box_p, iv, use_geometry=use_geometry)
end do
end subroutine af_restrict_box_vars
end subroutine af_restrict_box
!> Restrict only next to refinement boundaries, which which can be required
!> for filling coarse-grid ghost cells
......@@ -169,8 +147,7 @@ contains
! Only restrict near refinement boundaries
if (p_id > af_no_box .and. &
any(tree%boxes(id)%neighbors == af_no_box)) then
call af_restrict_box_vars(tree%boxes(id), tree%boxes(p_id), &
ivs)
call af_restrict_box(tree%boxes(id), tree%boxes(p_id), ivs)
end if
end do
!$omp end do
......
......@@ -422,12 +422,11 @@ module m_af_types
end subroutine af_subr_prolong
!> Subroutine for restriction
subroutine af_subr_restrict(box_c, box_p, iv, iv_to, use_geometry)
subroutine af_subr_restrict(box_c, box_p, ivs, use_geometry)
import
type(box_t), intent(in) :: box_c !< Child box to restrict
type(box_t), intent(inout) :: box_p !< Parent box to restrict to
integer, intent(in) :: iv !< Variable to restrict
integer, intent(in), optional :: iv_to !< Destination (if /= iv)
type(box_t), intent(in) :: box_c !< Child box to restrict
type(box_t), intent(inout) :: box_p !< Parent box to restrict to
integer, intent(in) :: ivs(:) !< Variables to restrict
logical, intent(in), optional :: use_geometry !< If set to false, don't use geometry
end subroutine af_subr_restrict
end interface
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment