This program can be used to benchmark the ghostcells routines.
1#include "../src/cpp_macros.h"
2
3
4
5program ghostcell_benchmark
7
8 implicit none
9
10 integer :: i_phi
11
12 type(af_t) :: tree
13 type(ref_info_t) :: ref_info
14 integer :: n_args
15 integer :: n_cell, it, n_iterations, max_ref_lvl
16 real(dp) :: time
17 character(len=100) :: arg_string
18 integer :: count_rate, t_start, t_end
19
20 print *, "Running ghostcell_benchmark_" // dimname // ""
21 print *, "Number of threads", af_get_max_threads()
22
23
24 n_args = command_argument_count()
25
26 if (n_args >= 1) then
27 call get_command_argument(1, arg_string)
28 read(arg_string, *) n_cell
29 else
30 print *, "No arguments specified, using default values"
31 print *, "Usage: ./ghostcell_benchmark_" // dimname // " n_cell max_ref_lvl n_iterations"
32 print *, ""
33 n_cell = 16
34 end if
35
36 if (n_args >= 2) then
37 call get_command_argument(2, arg_string)
38 read(arg_string, *) max_ref_lvl
39 else
40 max_ref_lvl = 4
41 end if
42
43 if (n_args >= 3) then
44 call get_command_argument(3, arg_string)
45 read(arg_string, *) n_iterations
46 else
47 n_iterations = 100
48 end if
49
50 print *, "Box size: ", n_cell
51 print *, "Max refinement lvl: ", max_ref_lvl
52 print *, "Num iterations: ", n_iterations
53
54 call af_add_cc_variable(tree, "phi", ix=i_phi)
55 call af_set_cc_methods(tree, i_phi, af_bc_dirichlet_zero)
56
57
58 call af_init(tree, &
59 n_cell, &
60 [dtimes(1.0_dp)], &
61 [dtimes(n_cell)])
62
63 call system_clock(t_start, count_rate)
64 do
65
66 call af_loop_box(tree, set_init_cond)
67
68
69
70
71
72 call af_adjust_refinement(tree, ref_routine, ref_info)
73
74
75 if (ref_info%n_add == 0) exit
76 end do
77 call system_clock(t_end, count_rate)
78
79 write(*,"(A,Es10.3,A)") " Wall-clock time generating AMR grid: ", &
80 (t_end-t_start) / real(count_rate,dp), " seconds"
81
82 call af_print_info(tree)
83
84
85 call system_clock(t_start, count_rate)
86 do it = 1, n_iterations
87 call af_gc_tree(tree, [i_phi], .false.)
88 end do
89 call system_clock(t_end, count_rate)
90
91 time = (t_end-t_start) / real(count_rate, dp)
92 write(*, "(A,I0,A,E10.3,A)") &
93 " Wall-clock time after ", n_iterations, &
94 " iterations: ", time, " seconds"
95 write(*, "(A,E10.3,A)") " Per iteration: ", time/n_iterations, " seconds"
96
97contains
98
99
100 subroutine ref_routine(box, cell_flags)
101 type(box_t), intent(in) :: box
102 integer, intent(out) :: cell_flags(DTIMES(box%n_cell))
103
104
105 if (box%lvl < max_ref_lvl) then
106 cell_flags = af_do_ref
107 else
108 cell_flags = af_keep_ref
109 end if
110 end subroutine ref_routine
111
112
113 subroutine set_init_cond(box)
114 type(box_t), intent(inout) :: box
115 integer :: nc
116
117 nc = box%n_cell
118 box%cc(dtimes(1:nc), i_phi) = 1.0_dp
119 end subroutine set_init_cond
120
121end program ghostcell_benchmark
Module which contains all Afivo modules, so that a user does not have to include them separately.