30 #define SACADO_VIEW_CUDA_HIERARCHICAL_DFAD 1
31 #define SACADO_KOKKOS_USE_MEMORY_POOL 1
37 #include "Kokkos_Timer.hpp"
39 template <
typename ViewTypeA,
typename ViewTypeB,
typename ViewTypeC>
42 typedef typename ViewTypeC::value_type scalar_type;
43 typedef typename ViewTypeC::execution_space execution_space;
45 #if defined (KOKKOS_ENABLE_CUDA)
47 const unsigned vector_size = is_cuda ? 32 : 1;
48 const unsigned team_size = is_cuda ? 128 / vector_size : 1;
49 #elif defined (KOKKOS_ENABLE_HIP)
51 const unsigned vector_size = is_hip ? 64 : 1;
52 const unsigned team_size = is_hip ? 128 / vector_size : 1;
54 const unsigned vector_size = 1;
55 const unsigned team_size = 1;
58 const int m = A.extent(0);
59 const int n = A.extent(1);
60 const int range = (m+team_size-1)/team_size;
62 typedef Kokkos::TeamPolicy<execution_space> Policy;
64 Policy( range,team_size,vector_size ),
65 KOKKOS_LAMBDA (
const typename Policy::member_type& team) {
66 const int i = team.league_rank()*team.team_size() + team.team_rank();
71 for (
int j=0; j<n; ++j)
78 template <
typename ViewTypeA,
typename ViewTypeB,
typename ViewTypeC>
80 const ViewTypeA&
A,
const ViewTypeB& b,
const ViewTypeC&
c) {
81 typedef typename ViewTypeC::value_type scalar_type;
82 typedef typename ViewTypeC::execution_space execution_space;
83 typedef Kokkos::TeamPolicy<execution_space> Policy;
84 typedef typename Policy::member_type team_member;
85 typedef Kokkos::View<scalar_type*,Kokkos::LayoutLeft, typename execution_space::scratch_memory_space, Kokkos::MemoryUnmanaged> TmpScratchSpace;
87 #if defined (KOKKOS_ENABLE_CUDA)
89 const unsigned VectorSize = is_cuda ? 32 : 1;
90 const unsigned TeamSize = is_cuda ? 128 / VectorSize : 1;
91 #elif defined (KOKKOS_ENABLE_HIP)
93 const unsigned VectorSize = is_hip ? 64 : 1;
94 const unsigned TeamSize = is_hip ? 128 / VectorSize : 1;
96 const unsigned VectorSize = 1;
97 const unsigned TeamSize = 1;
100 const int m = A.extent(0);
101 const int n = A.extent(1);
102 const int p = dimension_scalar(A);
103 const int N = (m+TeamSize-1)/TeamSize;
105 Policy policy(N, TeamSize, VectorSize);
106 const size_t bytes = TmpScratchSpace::shmem_size(TeamSize,p);
107 Kokkos::parallel_for(
108 policy.set_scratch_size(0, Kokkos::PerTeam(bytes)),
109 KOKKOS_LAMBDA (
const team_member& team) {
110 const int team_rank = team.team_rank();
111 const int team_size = team.team_size();
112 TmpScratchSpace t(team.team_scratch(0), team_size,
p);
113 const int i = team.league_rank()*team_size + team_rank;
116 for (
int j=0; j<n; ++j)
117 t(team_rank) +=
A(i,j)*b(j);
124 template <
typename ViewTypeA,
typename ViewTypeB,
typename ViewTypeC>
128 const double tol = 1.0e-14;
129 typedef typename ViewTypeC::value_type value_type;
130 typename ViewTypeC::HostMirror h_c = Kokkos::create_mirror_view(c);
131 Kokkos::deep_copy(h_c, c);
132 const size_t m = A.extent(0);
133 const size_t n = A.extent(1);
134 const size_t p = Kokkos::dimension_scalar(A);
135 for (
size_t i=0;
i<m; ++
i) {
136 for (
size_t j=0; j<
p; ++j) {
137 value_type t = (j == p-1 ? n : 2*n);
139 std::cout <<
"Comparison failed! " <<
i <<
"," << j <<
" : "
140 << h_c(
i).fastAccessDx(j) <<
" , " << t << std::endl;
146 template <
typename FadType,
typename ... ViewArgs>
149 const size_t nloop,
const bool check)
151 typedef Kokkos::View<
FadType**, ViewArgs...> ViewTypeA;
152 typedef Kokkos::View<FadType*, ViewArgs...> ViewTypeB;
153 typedef Kokkos::View<FadType*, ViewArgs...> ViewTypeC;
154 typedef typename ViewTypeA::execution_space execution_space;
158 typedef Kokkos::View<FadType**, ConLayoutA, execution_space> ConViewTypeA;
159 typedef Kokkos::View<FadType*, ConLayoutB, execution_space> ConViewTypeB;
160 typedef Kokkos::View<FadType*, ConLayoutC, execution_space> ConViewTypeC;
162 ConViewTypeA
A(
"A",m,n,p+1);
163 ConViewTypeB b(
"B",n,p+1);
164 ConViewTypeC
c(
"c",m,p+1);
169 Kokkos::deep_copy(
typename ConViewTypeA::array_type(A), 1.0);
170 Kokkos::deep_copy(
typename ConViewTypeB::array_type(b), 1.0);
172 Kokkos::Timer wall_clock;
175 #if defined (KOKKOS_ENABLE_CUDA)
177 const size_t warp_dim = is_cuda ? 32 : 1;
178 #elif defined (KOKKOS_ENABLE_HIP)
180 const size_t warp_dim = is_hip ? 64 : 1;
182 const size_t warp_dim = 1;
185 const size_t concurrency = execution_space().concurrency();
186 const size_t block_size = p*
sizeof(double);
187 const size_t nkernels = concurrency / warp_dim;
188 const size_t mem_pool_size =
189 static_cast<size_t>(1.2*nkernels*block_size);
190 const size_t superblock_size = std::max<size_t>(nkernels / 100, 1) * block_size;
191 execution_space space;
200 execution_space().fence();
203 for (
size_t l=0; l<nloop; l++) {
206 execution_space().fence();
208 perf.
time = wall_clock.seconds() / nloop;
221 template <
typename FadType,
typename ... ViewArgs>
224 const size_t m,
const size_t n,
const size_t p,
const size_t nloop,
227 typedef Kokkos::View<
FadType**, ViewArgs...> ViewTypeA;
228 typedef Kokkos::View<FadType*, ViewArgs...> ViewTypeB;
229 typedef Kokkos::View<FadType*, ViewArgs...> ViewTypeC;
230 typedef typename ViewTypeA::execution_space execution_space;
234 typedef Kokkos::View<FadType**, ConLayoutA, execution_space> ConViewTypeA;
235 typedef Kokkos::View<FadType*, ConLayoutB, execution_space> ConViewTypeB;
236 typedef Kokkos::View<FadType*, ConLayoutC, execution_space> ConViewTypeC;
238 ConViewTypeA
A(
"A",m,n,p+1);
239 ConViewTypeB b(
"B",n,p+1);
240 ConViewTypeC
c(
"c",m,p+1);
245 Kokkos::deep_copy(
typename ConViewTypeA::array_type(A), 1.0);
246 Kokkos::deep_copy(
typename ConViewTypeB::array_type(b), 1.0);
248 Kokkos::Timer wall_clock;
253 execution_space().fence();
256 for (
size_t l=0; l<nloop; l++) {
259 execution_space().fence();
261 perf.
time = wall_clock.seconds() / nloop;
274 #define INST_FUNC_FAD_DEV(FAD,DEV) \
275 template Perf do_time_fad_hierarchical_dfad< FAD, Kokkos::LayoutLeft, DEV > ( const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check ); \
276 template Perf do_time_fad_hierarchical_dfad< FAD, Kokkos::LayoutRight, DEV > ( const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check ); \
277 template Perf do_time_fad_hierarchical_dfad< FAD, DEV > ( const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check ); \
278 template Perf do_time_fad_hierarchical_dfad_scratch< FAD, Kokkos::LayoutLeft, DEV > ( const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check ); \
279 template Perf do_time_fad_hierarchical_dfad_scratch< FAD, Kokkos::LayoutRight, DEV > ( const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check ); \
280 template Perf do_time_fad_hierarchical_dfad_scratch< FAD, DEV > ( const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check );
282 #define INST_FUNC_DEV(DEV) \
283 INST_FUNC_FAD_DEV( DFad_type, DEV )
285 #ifdef KOKKOS_ENABLE_SERIAL
289 #ifdef KOKKOS_ENABLE_OPENMP
293 #ifdef KOKKOS_ENABLE_THREADS
297 #ifdef KOKKOS_ENABLE_CUDA
301 #ifdef KOKKOS_ENABLE_HIP
Sacado::Fad::DFad< double > DFad_type
std::enable_if< !Kokkos::is_view_fad< View2 >::value, bool >::type check(const View1 &v_gold, const View2 &v, const double tol)
void createGlobalMemoryPool(const ExecSpace &space, const size_t min_total_alloc_size, const uint32_t min_block_alloc_size, const uint32_t max_block_alloc_size, const uint32_t min_superblock_size)
void check_deriv_hierarchical_dfad(const ViewTypeA &A, const ViewTypeB &b, const ViewTypeC &c)
Sacado::Fad::DFad< double > FadType
Perf do_time_fad_hierarchical_dfad_scratch(const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check)
expr expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c *expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr1 c expr2 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 expr2 expr1 expr2 expr1 expr1 expr1 c
Perf do_time_fad_hierarchical_dfad(const size_t m, const size_t n, const size_t p, const size_t nloop, const bool check)
#define INST_FUNC_DEV(DEV)
expr expr expr fastAccessDx(i)) FAD_UNARYOP_MACRO(exp
void run_mat_vec_hierarchical_dfad_scratch(const ViewTypeA &A, const ViewTypeB &b, const ViewTypeC &c)
void destroyGlobalMemoryPool(const ExecSpace &space)
void run_mat_vec_hierarchical_dfad(const ViewTypeA &A, const ViewTypeB &b, const ViewTypeC &c)