Reference documentation for deal.II version GIT 9042b9283b 2023-12-02 14:50:02+00:00
\(\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}}\) \(\newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=}\) \(\newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]}\) \(\newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}\)
cuda_tensor_product_kernels.h
Go to the documentation of this file.
1 // ---------------------------------------------------------------------
2 //
3 // Copyright (C) 2017 - 2023 by the deal.II authors
4 //
5 // This file is part of the deal.II library.
6 //
7 // The deal.II library is free software; you can use it, redistribute
8 // it, and/or modify it under the terms of the GNU Lesser General
9 // Public License as published by the Free Software Foundation; either
10 // version 2.1 of the License, or (at your option) any later version.
11 // The full text of the license can be found in the file LICENSE.md at
12 // the top level directory of deal.II.
13 //
14 // ---------------------------------------------------------------------
15 
16 
17 #ifndef dealii_cuda_tensor_product_kernels_h
18 #define dealii_cuda_tensor_product_kernels_h
19 
20 #include <deal.II/base/config.h>
21 
22 #include <deal.II/base/utilities.h>
23 
24 #include <deal.II/matrix_free/cuda_matrix_free.templates.h>
25 
27 
28 
29 namespace CUDAWrappers
30 {
31  namespace internal
32  {
39  // TODO: for now only the general variant is implemented
41  {
45  };
46 
47 
48 
49 #if KOKKOS_VERSION >= 40000
53  template <int n_q_points_1d,
54  typename Number,
55  int direction,
56  bool dof_to_quad,
57  bool add,
58  bool in_place,
59  typename ViewTypeIn,
60  typename ViewTypeOut>
62  apply_1d(const Kokkos::TeamPolicy<
63  MemorySpace::Default::kokkos_space::execution_space>::member_type
64  &team_member,
65  const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
66  shape_data,
67  const ViewTypeIn in,
68  ViewTypeOut out)
69  {
70  Number t[n_q_points_1d];
71  Kokkos::parallel_for(Kokkos::TeamThreadRange(team_member, n_q_points_1d),
72  [&](const int &q) {
73  t[q] = 0;
74  // This loop simply multiplies the shape function
75  // at the quadrature point by the value finite
76  // element coefficient.
77  // FIXME check why using parallel_reduce
78  // ThreadVector is slower
79  for (int k = 0; k < n_q_points_1d; ++k)
80  {
81  const unsigned int shape_idx =
82  dof_to_quad ? (q + k * n_q_points_1d) :
83  (k + q * n_q_points_1d);
84  const unsigned int source_idx = k;
85  t[q] += shape_data[shape_idx] * in(source_idx);
86  }
87  });
88 
89  if constexpr (in_place)
90  team_member.team_barrier();
91 
92  Kokkos::parallel_for(Kokkos::TeamThreadRange(team_member, n_q_points_1d),
93  [&](const int &q) {
94  const unsigned int destination_idx = q;
95  if constexpr (add)
96  Kokkos::atomic_add(&out(destination_idx), t[q]);
97  else
98  out(destination_idx) = t[q];
99  });
100  }
101 
102 
103 
107  template <int n_q_points_1d,
108  typename Number,
109  int direction,
110  bool dof_to_quad,
111  bool add,
112  bool in_place,
113  typename ViewTypeIn,
114  typename ViewTypeOut>
116  apply_2d(const Kokkos::TeamPolicy<
117  MemorySpace::Default::kokkos_space::execution_space>::member_type
118  &team_member,
119  const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
120  shape_data,
121  const ViewTypeIn in,
122  ViewTypeOut out)
123  {
124  using TeamType = Kokkos::TeamPolicy<
125  MemorySpace::Default::kokkos_space::execution_space>::member_type;
126  constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, 2);
127 
128  Number t[n_q_points];
129  auto thread_policy =
130  Kokkos::TeamThreadMDRange<Kokkos::Rank<2>, TeamType>(team_member,
131  n_q_points_1d,
132  n_q_points_1d);
133  Kokkos::parallel_for(thread_policy, [&](const int i, const int j) {
134  int q_point = i + j * n_q_points_1d;
135  t[q_point] = 0;
136 
137  // This loop simply multiplies the shape function at the quadrature
138  // point by the value finite element coefficient.
139  // FIXME check why using parallel_reduce ThreadVector is slower
140  for (int k = 0; k < n_q_points_1d; ++k)
141  {
142  const unsigned int shape_idx =
143  dof_to_quad ? (j + k * n_q_points_1d) : (k + j * n_q_points_1d);
144  const unsigned int source_idx = (direction == 0) ?
145  (k + n_q_points_1d * i) :
146  (i + n_q_points_1d * k);
147  t[q_point] += shape_data[shape_idx] * in(source_idx);
148  }
149  });
150 
151  if (in_place)
152  team_member.team_barrier();
153 
154  Kokkos::parallel_for(thread_policy, [&](const int i, const int j) {
155  const int q_point = i + j * n_q_points_1d;
156  const unsigned int destination_idx =
157  (direction == 0) ? (j + n_q_points_1d * i) : (i + n_q_points_1d * j);
158 
159  if (add)
160  Kokkos::atomic_add(&out(destination_idx), t[q_point]);
161  else
162  out(destination_idx) = t[q_point];
163  });
164  }
165 
166 
167 
171  template <int n_q_points_1d,
172  typename Number,
173  int direction,
174  bool dof_to_quad,
175  bool add,
176  bool in_place,
177  typename ViewTypeIn,
178  typename ViewTypeOut>
180  apply_3d(const Kokkos::TeamPolicy<
181  MemorySpace::Default::kokkos_space::execution_space>::member_type
182  &team_member,
183  const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
184  shape_data,
185  const ViewTypeIn in,
186  ViewTypeOut out)
187  {
188  using TeamType = Kokkos::TeamPolicy<
189  MemorySpace::Default::kokkos_space::execution_space>::member_type;
190  constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, 3);
191 
192  Number t[n_q_points];
193  auto thread_policy = Kokkos::TeamThreadMDRange<Kokkos::Rank<3>, TeamType>(
194  team_member, n_q_points_1d, n_q_points_1d, n_q_points_1d);
196  thread_policy, [&](const int i, const int j, const int q) {
197  const int q_point =
198  i + j * n_q_points_1d + q * n_q_points_1d * n_q_points_1d;
199  t[q_point] = 0;
200 
201  // This loop simply multiplies the shape function at the quadrature
202  // point by the value finite element coefficient.
203  // FIXME check why using parallel_reduce ThreadVector is slower
204  for (int k = 0; k < n_q_points_1d; ++k)
205  {
206  const unsigned int shape_idx =
207  dof_to_quad ? (q + k * n_q_points_1d) : (k + q * n_q_points_1d);
208  const unsigned int source_idx =
209  (direction == 0) ?
210  (k + n_q_points_1d * (i + n_q_points_1d * j)) :
211  (direction == 1) ?
212  (i + n_q_points_1d * (k + n_q_points_1d * j)) :
213  (i + n_q_points_1d * (j + n_q_points_1d * k));
214  t[q_point] += shape_data[shape_idx] * in(source_idx);
215  }
216  });
217 
218  if (in_place)
219  team_member.team_barrier();
220 
222  thread_policy, [&](const int i, const int j, const int q) {
223  const int q_point =
224  i + j * n_q_points_1d + q * n_q_points_1d * n_q_points_1d;
225  const unsigned int destination_idx =
226  (direction == 0) ? (q + n_q_points_1d * (i + n_q_points_1d * j)) :
227  (direction == 1) ? (i + n_q_points_1d * (q + n_q_points_1d * j)) :
228  (i + n_q_points_1d * (j + n_q_points_1d * q));
229 
230  if (add)
231  Kokkos::atomic_add(&out(destination_idx), t[q_point]);
232  else
233  out(destination_idx) = t[q_point];
234  });
235  }
236 #endif
237 
238 
239 
243  template <int dim,
244  int n_q_points_1d,
245  typename Number,
246  int direction,
247  bool dof_to_quad,
248  bool add,
249  bool in_place,
250  typename ViewTypeIn,
251  typename ViewTypeOut>
253  apply(const Kokkos::TeamPolicy<
254  MemorySpace::Default::kokkos_space::execution_space>::member_type
255  &team_member,
256  const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
257  shape_data,
258  const ViewTypeIn in,
259  ViewTypeOut out)
260  {
261 #if KOKKOS_VERSION >= 40000
262  if constexpr (dim == 1)
263  apply_1d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
264  team_member, shape_data, in, out);
265  if constexpr (dim == 2)
266  apply_2d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
267  team_member, shape_data, in, out);
268  if constexpr (dim == 3)
269  apply_3d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
270  team_member, shape_data, in, out);
271 #else
272  constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim);
273 
274  Number t[n_q_points];
276  Kokkos::TeamThreadRange(team_member, n_q_points),
277  [&](const int &q_point) {
278  const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
279  const unsigned int j =
280  (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
281  const unsigned int q =
282  (dim == 1) ? q_point :
283  (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
284  q_point / (n_q_points_1d * n_q_points_1d);
285 
286  // This loop simply multiplies the shape function at the quadrature
287  // point by the value finite element coefficient.
288  t[q_point] = 0;
289  for (int k = 0; k < n_q_points_1d; ++k)
290  {
291  const unsigned int shape_idx =
292  dof_to_quad ? (q + k * n_q_points_1d) : (k + q * n_q_points_1d);
293  const unsigned int source_idx =
294  (direction == 0) ?
295  (k + n_q_points_1d * (i + n_q_points_1d * j)) :
296  (direction == 1) ?
297  (i + n_q_points_1d * (k + n_q_points_1d * j)) :
298  (i + n_q_points_1d * (j + n_q_points_1d * k));
299  t[q_point] += shape_data[shape_idx] *
300  (in_place ? out(source_idx) : in(source_idx));
301  }
302  });
303 
304  if (in_place)
305  team_member.team_barrier();
306 
308  Kokkos::TeamThreadRange(team_member, n_q_points),
309  [&](const int &q_point) {
310  const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
311  const unsigned int j =
312  (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
313  const unsigned int q =
314  (dim == 1) ? q_point :
315  (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
316  q_point / (n_q_points_1d * n_q_points_1d);
317 
318  const unsigned int destination_idx =
319  (direction == 0) ? (q + n_q_points_1d * (i + n_q_points_1d * j)) :
320  (direction == 1) ? (i + n_q_points_1d * (q + n_q_points_1d * j)) :
321  (i + n_q_points_1d * (j + n_q_points_1d * q));
322 
323  if (add)
324  Kokkos::atomic_add(&out(destination_idx), t[q_point]);
325  else
326  out(destination_idx) = t[q_point];
327  });
328 #endif
329  }
330 
331 
337  template <EvaluatorVariant variant,
338  int dim,
339  int fe_degree,
340  int n_q_points_1d,
341  typename Number>
343  {};
344 
345 
346 
353  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
355  dim,
356  fe_degree,
357  n_q_points_1d,
358  Number>
359  {
360  using TeamHandle = Kokkos::TeamPolicy<
361  MemorySpace::Default::kokkos_space::execution_space>::member_type;
362 
365  const TeamHandle &team_member,
366  Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
367  Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
368  shape_gradients,
369  Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
370  co_shape_gradients);
371 
376  template <int direction,
377  bool dof_to_quad,
378  bool add,
379  bool in_place,
380  typename ViewTypeIn,
381  typename ViewTypeOut>
383  values(const ViewTypeIn in, ViewTypeOut out) const;
384 
389  template <int direction,
390  bool dof_to_quad,
391  bool add,
392  bool in_place,
393  typename ViewTypeIn,
394  typename ViewTypeOut>
396  gradients(const ViewTypeIn in, ViewTypeOut out) const;
397 
402  template <int direction,
403  bool dof_to_quad,
404  bool add,
405  bool in_place,
406  typename ViewTypeIn,
407  typename ViewTypeOut>
409  co_gradients(const ViewTypeIn in, ViewTypeOut out) const;
410 
414  template <typename ViewType>
416  value_at_quad_pts(ViewType u);
417 
421  template <typename ViewType>
423  integrate_value(ViewType u);
424 
429  template <typename ViewTypeIn, typename ViewTypeOut>
431  gradient_at_quad_pts(const ViewTypeIn u, ViewTypeOut grad_u);
432 
437  template <typename ViewType1, typename ViewType2>
439  value_and_gradient_at_quad_pts(ViewType1 u, ViewType2 grad_u);
440 
445  template <bool add, typename ViewType1, typename ViewType2>
447  integrate_gradient(ViewType1 u, ViewType2 grad_u);
448 
453  template <typename ViewType1, typename ViewType2>
455  integrate_value_and_gradient(ViewType1 u, ViewType2 grad_u);
456 
461 
465  Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values;
466 
470  Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
472 
476  Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
478  };
479 
480 
481 
482  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
485  dim,
486  fe_degree,
487  n_q_points_1d,
488  Number>::
489  EvaluatorTensorProduct(
490  const TeamHandle &team_member,
491  Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
492  Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
493  shape_gradients,
494  Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
495  co_shape_gradients)
496  : team_member(team_member)
497  , shape_values(shape_values)
498  , shape_gradients(shape_gradients)
499  , co_shape_gradients(co_shape_gradients)
500  {}
501 
502 
503 
504  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
505  template <int direction,
506  bool dof_to_quad,
507  bool add,
508  bool in_place,
509  typename ViewTypeIn,
510  typename ViewTypeOut>
513  dim,
514  fe_degree,
515  n_q_points_1d,
516  Number>::values(const ViewTypeIn in,
517  ViewTypeOut out) const
518  {
519  apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
520  team_member, shape_values, in, out);
521  }
522 
523 
524 
525  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
526  template <int direction,
527  bool dof_to_quad,
528  bool add,
529  bool in_place,
530  typename ViewTypeIn,
531  typename ViewTypeOut>
534  dim,
535  fe_degree,
536  n_q_points_1d,
537  Number>::gradients(const ViewTypeIn in,
538  ViewTypeOut out) const
539  {
540  apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
541  team_member, shape_gradients, in, out);
542  }
543 
544 
545 
546  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
547  template <int direction,
548  bool dof_to_quad,
549  bool add,
550  bool in_place,
551  typename ViewTypeIn,
552  typename ViewTypeOut>
555  dim,
556  fe_degree,
557  n_q_points_1d,
558  Number>::co_gradients(const ViewTypeIn in,
559  ViewTypeOut out) const
560  {
561  apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
562  team_member, co_shape_gradients, in, out);
563  }
564 
565 
566 
567  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
568  template <typename ViewType>
569  DEAL_II_HOST_DEVICE inline void
571  dim,
572  fe_degree,
573  n_q_points_1d,
574  Number>::value_at_quad_pts(ViewType u)
575  {
576  if constexpr (dim == 1)
577  values<0, true, false, true>(u, u);
578  else if constexpr (dim == 2)
579  {
580  values<0, true, false, true>(u, u);
581  team_member.team_barrier();
582  values<1, true, false, true>(u, u);
583  }
584  else if constexpr (dim == 3)
585  {
586  values<0, true, false, true>(u, u);
587  team_member.team_barrier();
588  values<1, true, false, true>(u, u);
589  team_member.team_barrier();
590  values<2, true, false, true>(u, u);
591  }
592  else
593  Kokkos::abort("dim must not exceed 3!");
594  }
595 
596 
597 
598  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
599  template <typename ViewType>
600  DEAL_II_HOST_DEVICE inline void
602  dim,
603  fe_degree,
604  n_q_points_1d,
605  Number>::integrate_value(ViewType u)
606  {
607  if constexpr (dim == 1)
608  values<0, false, false, true>(u, u);
609  else if constexpr (dim == 2)
610  {
611  values<0, false, false, true>(u, u);
612  team_member.team_barrier();
613  values<1, false, false, true>(u, u);
614  }
615  else if constexpr (dim == 3)
616  {
617  values<0, false, false, true>(u, u);
618  team_member.team_barrier();
619  values<1, false, false, true>(u, u);
620  team_member.team_barrier();
621  values<2, false, false, true>(u, u);
622  }
623  else
624  Kokkos::abort("dim must not exceed 3!");
625  }
626 
627 
628 
629  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
630  template <typename ViewTypeIn, typename ViewTypeOut>
631  DEAL_II_HOST_DEVICE inline void
633  dim,
634  fe_degree,
635  n_q_points_1d,
636  Number>::gradient_at_quad_pts(const ViewTypeIn u,
637  ViewTypeOut grad_u)
638  {
639  if constexpr (dim == 1)
640  {
641  gradients<0, true, false, false>(
642  u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
643  }
644  else if constexpr (dim == 2)
645  {
646  gradients<0, true, false, false>(
647  u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
648  values<0, true, false, false>(
649  u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
650 
651  team_member.team_barrier();
652 
653  values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
654  Kokkos::subview(grad_u, Kokkos::ALL, 0));
655  gradients<1, true, false, true>(
656  Kokkos::subview(grad_u, Kokkos::ALL, 1),
657  Kokkos::subview(grad_u, Kokkos::ALL, 1));
658  }
659  else if constexpr (dim == 3)
660  {
661  gradients<0, true, false, false>(
662  u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
663  values<0, true, false, false>(
664  u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
665  values<0, true, false, false>(
666  u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
667 
668  team_member.team_barrier();
669 
670  values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
671  Kokkos::subview(grad_u, Kokkos::ALL, 0));
672  gradients<1, true, false, true>(
673  Kokkos::subview(grad_u, Kokkos::ALL, 1),
674  Kokkos::subview(grad_u, Kokkos::ALL, 1));
675  values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
676  Kokkos::subview(grad_u, Kokkos::ALL, 2));
677 
678  team_member.team_barrier();
679 
680  values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
681  Kokkos::subview(grad_u, Kokkos::ALL, 0));
682  values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
683  Kokkos::subview(grad_u, Kokkos::ALL, 1));
684  gradients<2, true, false, true>(
685  Kokkos::subview(grad_u, Kokkos::ALL, 2),
686  Kokkos::subview(grad_u, Kokkos::ALL, 2));
687  }
688  else
689  Kokkos::abort("dim must not exceed 3!");
690  }
691 
692 
693 
694  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
695  template <typename ViewType1, typename ViewType2>
696  DEAL_II_HOST_DEVICE inline void
698  dim,
699  fe_degree,
700  n_q_points_1d,
701  Number>::value_and_gradient_at_quad_pts(ViewType1 u,
702  ViewType2
703  grad_u)
704  {
705  if constexpr (dim == 1)
706  {
707  values<0, true, false, true>(u, u);
708  team_member.team_barrier();
709 
710  co_gradients<0, true, false, false>(
711  u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
712  }
713  else if constexpr (dim == 2)
714  {
715  values<0, true, false, true>(u, u);
716  team_member.team_barrier();
717  values<1, true, false, true>(u, u);
718  team_member.team_barrier();
719 
720  co_gradients<0, true, false, false>(
721  u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
722  co_gradients<1, true, false, false>(
723  u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
724  }
725  else if constexpr (dim == 3)
726  {
727  values<0, true, false, true>(u, u);
728  team_member.team_barrier();
729  values<1, true, false, true>(u, u);
730  team_member.team_barrier();
731  values<2, true, false, true>(u, u);
732  team_member.team_barrier();
733 
734  co_gradients<0, true, false, false>(
735  u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
736  co_gradients<1, true, false, false>(
737  u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
738  co_gradients<2, true, false, false>(
739  u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
740  }
741  else
742  Kokkos::abort("dim must not exceed 3!");
743  }
744 
745 
746 
747  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
748  template <bool add, typename ViewType1, typename ViewType2>
749  DEAL_II_HOST_DEVICE inline void
751  dim,
752  fe_degree,
753  n_q_points_1d,
754  Number>::integrate_gradient(ViewType1 u,
755  ViewType2 grad_u)
756  {
757  if constexpr (dim == 1)
758  {
759  gradients<0, false, add, false>(
760  Kokkos::subview(grad_u, Kokkos::ALL, dim), u);
761  }
762  else if constexpr (dim == 2)
763  {
764  gradients<0, false, false, true>(
765  Kokkos::subview(grad_u, Kokkos::ALL, 0),
766  Kokkos::subview(grad_u, Kokkos::ALL, 0));
767  values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
768  Kokkos::subview(grad_u,
769  Kokkos::ALL,
770  1));
771 
772  team_member.team_barrier();
773 
774  values<1, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
775  u);
776  team_member.team_barrier();
777  gradients<1, false, true, false>(
778  Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
779  }
780  else if constexpr (dim == 3)
781  {
782  gradients<0, false, false, true>(
783  Kokkos::subview(grad_u, Kokkos::ALL, 0),
784  Kokkos::subview(grad_u, Kokkos::ALL, 0));
785  values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
786  Kokkos::subview(grad_u,
787  Kokkos::ALL,
788  1));
789  values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
790  Kokkos::subview(grad_u,
791  Kokkos::ALL,
792  2));
793 
794  team_member.team_barrier();
795 
796  values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
797  Kokkos::subview(grad_u,
798  Kokkos::ALL,
799  0));
800  gradients<1, false, false, true>(
801  Kokkos::subview(grad_u, Kokkos::ALL, 1),
802  Kokkos::subview(grad_u, Kokkos::ALL, 1));
803  values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
804  Kokkos::subview(grad_u,
805  Kokkos::ALL,
806  2));
807 
808  team_member.team_barrier();
809 
810  values<2, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
811  u);
812  team_member.team_barrier();
813  values<2, false, true, false>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
814  u);
815  team_member.team_barrier();
816  gradients<2, false, true, false>(
817  Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
818  }
819  else
820  Kokkos::abort("dim must not exceed 3!");
821  }
822 
823 
824 
825  template <int dim, int fe_degree, int n_q_points_1d, typename Number>
826  template <typename ViewType1, typename ViewType2>
827  DEAL_II_HOST_DEVICE inline void
829  dim,
830  fe_degree,
831  n_q_points_1d,
832  Number>::integrate_value_and_gradient(ViewType1 u,
833  ViewType2
834  grad_u)
835  {
836  if constexpr (dim == 1)
837  {
838  co_gradients<0, false, true, false>(
839  Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
840  team_member.team_barrier();
841 
842  values<0, false, false, true>(u, u);
843  }
844  else if constexpr (dim == 2)
845  {
846  co_gradients<1, false, true, false>(
847  Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
848  team_member.team_barrier();
849  co_gradients<0, false, true, false>(
850  Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
851  team_member.team_barrier();
852 
853  values<1, false, false, true>(u, u);
854  team_member.team_barrier();
855  values<0, false, false, true>(u, u);
856  team_member.team_barrier();
857  }
858  else if constexpr (dim == 3)
859  {
860  co_gradients<2, false, true, false>(
861  Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
862  team_member.team_barrier();
863  co_gradients<1, false, true, false>(
864  Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
865  team_member.team_barrier();
866  co_gradients<0, false, true, false>(
867  Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
868  team_member.team_barrier();
869 
870  values<2, false, false, true>(u, u);
871  team_member.team_barrier();
872  values<1, false, false, true>(u, u);
873  team_member.team_barrier();
874  values<0, false, false, true>(u, u);
875  team_member.team_barrier();
876  }
877  else
878  Kokkos::abort("dim must not exceed 3!");
879  }
880  } // namespace internal
881 } // namespace CUDAWrappers
882 
884 
885 #endif
#define DEAL_II_NAMESPACE_OPEN
Definition: config.h:477
#define DEAL_II_NAMESPACE_CLOSE
Definition: config.h:478
void apply(const Kokkos::TeamPolicy< MemorySpace::Default::kokkos_space::execution_space >::member_type &team_member, const Kokkos::View< Number *, MemorySpace::Default::kokkos_space > shape_data, const ViewTypeIn in, ViewTypeOut out)
constexpr T pow(const T base, const int iexp)
Definition: utilities.h:447
void abort(const ExceptionBase &exc) noexcept
Definition: exceptions.cc:460
void parallel_for(Iterator x_begin, Iterator x_end, const Functor &functor, const unsigned int grainsize)
Definition: parallel.h:84
#define DEAL_II_HOST_DEVICE
Definition: numbers.h:35
Kokkos::TeamPolicy< MemorySpace::Default::kokkos_space::execution_space >::member_type TeamHandle