Afivo 0.3
All Classes Namespaces Functions Variables Pages
random_refinement.f90

This example shows how to create an AMR tree, perform random refinement, write output files, and how to fill ghost cells.

This example shows how to create an AMR tree, perform random refinement, write output files, and how to fill ghost cells.

1#include "../src/cpp_macros.h"
2!> \example random_refinement.f90
3!>
4!> This example shows how to create an AMR tree, perform random refinement,
5!> write output files, and how to fill ghost cells.
6program random_refinement
7 use m_af_all
8
9 implicit none
10
11 type(af_t) :: tree
12 integer :: grid_size(NDIM)
13 real(dp) :: domain_size(NDIM)
14 logical :: periodic(NDIM) = .true.
15 integer :: iter, boxes_used
16 integer, parameter :: coord_type = af_xyz ! af_xyz or af_cyl
17
18 integer, parameter :: box_size = 8
19 integer, parameter :: i_phi = 1
20 type(ref_info_t) :: ref_info
21 real(dp) :: sum_phi_t0, sum_phi
22 character(len=100) :: fname
23 integer :: count_rate,t_start,t_end
24
25 write(*,'(A,I0,A)') 'program random_refinement_', ndim, "d"
26
27 print *, "Number of threads", af_get_max_threads()
28
29 call af_add_cc_variable(tree, "phi")
30
31 call af_set_cc_methods(tree, 1, af_bc_dirichlet_zero, &
32 prolong=af_prolong_linear)
33
34 ! Initialize tree
35 grid_size(1) = 2 * box_size
36 grid_size(2:) = box_size
37 domain_size(1) = 2 * acos(-1.0_dp)
38 domain_size(2:) = acos(-1.0_dp)
39
40 if (command_argument_count() == 1) then
41 ! Load tree from file
42 call get_command_argument(1, fname)
43 call af_read_tree(tree, fname, read_other_data=read_string)
44 else
45 call af_init(tree, & ! Tree to initialize
46 box_size, & ! A box contains box_size**DIM cells
47 domain_size, &
48 grid_size, &
49 periodic=periodic, &
50 coord=coord_type)
51
52 call af_print_info(tree)
53
54 ! Set variables on base by using the helper functions af_loop_box(tree, sub)
55 ! and af_loop_boxes(tree, sub). These functions call the subroutine sub for
56 ! each box in the tree, with a slightly different syntax.
57 call af_loop_box(tree, set_init_cond)
58
59 ! Fill ghost cells for phi. The third argument is a subroutine that fills
60 ! ghost cells near refinement boundaries, and the fourth argument fill ghost
61 ! cells near physical boundaries.
62 call af_gc_tree(tree, [i_phi])
63 end if
64
65 call af_tree_sum_cc(tree, i_phi, sum_phi_t0)
66
67 call system_clock(t_start, count_rate)
68 boxes_used = 1
69 do iter = 1, 50
70 ! This writes a VTK output file containing the cell-centered values of the
71 ! leaves of the tree (the boxes not covered by refinement).
72 ! Variables are the names given as the third argument.
73 write(fname, "(A,I0)") "output/random_refinement_" // dimname // "_", iter
74 call af_write_silo(tree, trim(fname), n_cycle=iter)
75 call af_write_tree(tree, fname, write_other_data=write_string)
76
77 ! This updates the refinement of the tree, by at most one level per call.
78 ! The second argument is a subroutine that is called for each box that can
79 ! be refined or derefined, and it should set refinement flags. Information
80 ! about the changes in refinement are returned in the third argument.
81 !
82 ! Within each (de)refinement step the subroutine af_tidy_up is called. This
83 ! means that every now and then whether too much holes appear in the tree
84 ! box list. Therefore reordering of tree is not necessary at the end of the
85 ! iteration process
86 call af_adjust_refinement(tree, ref_routine, ref_info, ref_buffer=0)
87
88 call af_tree_sum_cc(tree, i_phi, sum_phi)
89 write(*, "(A,E10.2)") " conservation error: ", sum_phi - sum_phi_t0
90 boxes_used = boxes_used + ref_info%n_add - ref_info%n_rm
91 write(*,'(4(3x,A,1x,i6))') "# new boxes", ref_info%n_add, &
92 "# removed boxes", ref_info%n_rm, &
93 "# boxes used ", boxes_used, &
94 " highest level ", tree%highest_lvl
95 print *, "boxes / in use / not used: ", tree%highest_id, &
96 count(tree%boxes(1:tree%highest_id)%in_use), &
97 tree%n_removed_ids
98 end do
99 call system_clock(t_end, count_rate)
100
101 write(*, '(A,i3,1x,A,f8.2,1x,A,/)') &
102 ' Wall-clock time after ',iter, &
103 ' iterations: ', (t_end-t_start) / real(count_rate, dp), &
104 ' seconds'
105
106 call af_print_info(tree)
107
108 ! This call is not really necessary here, but cleaning up the data in a tree
109 ! is important if your program continues with other tasks.
110 call af_destroy(tree)
111
112contains
113
114 ! Return the refinement flag for boxes(id)
115 subroutine ref_routine(box, cell_flags)
116 type(box_t), intent(in) :: box ! A list of all boxes in the tree
117 integer, intent(out) :: cell_flags(DTIMES(box%n_cell))
118 real(dp) :: rr
119
120 ! Draw a [0, 1) random number
121 call random_number(rr)
122
123 if (rr < 0.5_dp**ndim .and. box%lvl < 5) then
124 cell_flags = af_do_ref ! Add refinement
125 else
126 cell_flags = af_rm_ref ! Ask to remove this box, which will not always
127 ! happen (see documentation)
128 end if
129 end subroutine ref_routine
130
131 ! This routine sets the initial conditions for each box
132 subroutine set_init_cond(box)
133 type(box_t), intent(inout) :: box
134 integer :: IJK, nc
135 real(dp) :: rr(NDIM)
136
137 nc = box%n_cell
138 do kji_do(1,nc)
139 ! Get the coordinate of the cell center at i,j
140 rr = af_r_cc(box, [ijk])
141
142 ! Set the values at each cell according to some function
143#if NDIM > 1
144 box%cc(ijk, i_phi) = sin(0.5_dp * rr(1))**2 * cos(rr(2))**2
145#else
146 box%cc(ijk, i_phi) = sin(0.5_dp * rr(1))**2
147#endif
148 end do; close_do
149 end subroutine set_init_cond
150
151 subroutine write_string(my_unit)
152 integer, intent(in) :: my_unit
153 character(len=10) :: str
154
155 str = "Hello"
156 write(my_unit) str
157 end subroutine write_string
158
159 subroutine read_string(my_unit)
160 integer, intent(in) :: my_unit
161 character(len=10) :: str
162 read(my_unit) str
163
164 if (str /= "Hello") error stop "error reading other data"
165 end subroutine read_string
166
167end program random_refinement
Module which contains all Afivo modules, so that a user does not have to include them separately.
Definition: m_af_all.f90:3