Reference documentation for deal.II version GIT relicensing-422-gb369f187d8 2024-04-17 18:10: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\}}\)
Loading...
Searching...
No Matches
portable_tensor_product_kernels.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 2017 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15
16#ifndef dealii__tensor_product_kernels_h
17#define dealii__tensor_product_kernels_h
18
19#include <deal.II/base/config.h>
20
22
23#include <deal.II/matrix_free/cuda_matrix_free.templates.h>
24
26
27
28namespace Portable
29{
30 namespace internal
31 {
39 // TODO: for now only the general variant is implemented
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);
195 Kokkos::parallel_for(
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
221 Kokkos::parallel_for(
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];
275 Kokkos::parallel_for(
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
307 Kokkos::parallel_for(
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
338 template <EvaluatorVariant variant,
339 int dim,
340 int fe_degree,
341 int n_q_points_1d,
342 typename Number>
345
346
347
355 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
357 dim,
358 fe_degree,
359 n_q_points_1d,
360 Number>
361 {
362 public:
363 using TeamHandle = Kokkos::TeamPolicy<
364 MemorySpace::Default::kokkos_space::execution_space>::member_type;
365
368 const TeamHandle &team_member,
369 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
370 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
371 shape_gradients,
372 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
373 co_shape_gradients);
374
378 template <typename ViewType>
380 evaluate_values(ViewType u);
381
386 template <typename ViewTypeIn, typename ViewTypeOut>
388 evaluate_gradients(const ViewTypeIn u, ViewTypeOut grad_u);
389
394 template <typename ViewType1, typename ViewType2>
396 evaluate_values_and_gradients(ViewType1 u, ViewType2 grad_u);
397
401 template <typename ViewType>
403 integrate_values(ViewType u);
404
409 template <bool add, typename ViewType1, typename ViewType2>
411 integrate_gradients(ViewType1 u, ViewType2 grad_u);
412
417 template <typename ViewType1, typename ViewType2>
419 integrate_values_and_gradients(ViewType1 u, ViewType2 grad_u);
420
425 template <int direction,
426 bool dof_to_quad,
427 bool add,
428 bool in_place,
429 typename ViewTypeIn,
430 typename ViewTypeOut>
432 values(const ViewTypeIn in, ViewTypeOut out) const;
433
438 template <int direction,
439 bool dof_to_quad,
440 bool add,
441 bool in_place,
442 typename ViewTypeIn,
443 typename ViewTypeOut>
445 gradients(const ViewTypeIn in, ViewTypeOut out) const;
446
447 public:
452 template <int direction,
453 bool dof_to_quad,
454 bool add,
455 bool in_place,
456 typename ViewTypeIn,
457 typename ViewTypeOut>
459 co_gradients(const ViewTypeIn in, ViewTypeOut out) const;
460
465
469 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values;
470
474 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
476
480 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
482 };
483
484
485
486 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
489 dim,
490 fe_degree,
491 n_q_points_1d,
492 Number>::
493 EvaluatorTensorProduct(
494 const TeamHandle &team_member,
495 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
496 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
497 shape_gradients,
498 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
499 co_shape_gradients)
500 : team_member(team_member)
501 , shape_values(shape_values)
502 , shape_gradients(shape_gradients)
503 , co_shape_gradients(co_shape_gradients)
504 {}
505
506
507
508 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
509 template <int direction,
510 bool dof_to_quad,
511 bool add,
512 bool in_place,
513 typename ViewTypeIn,
514 typename ViewTypeOut>
517 dim,
518 fe_degree,
519 n_q_points_1d,
520 Number>::values(const ViewTypeIn in,
521 ViewTypeOut out) const
522 {
523 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
524 team_member, shape_values, in, out);
525 }
526
527
528
529 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
530 template <int direction,
531 bool dof_to_quad,
532 bool add,
533 bool in_place,
534 typename ViewTypeIn,
535 typename ViewTypeOut>
538 dim,
539 fe_degree,
540 n_q_points_1d,
541 Number>::gradients(const ViewTypeIn in,
542 ViewTypeOut out) const
543 {
544 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
545 team_member, shape_gradients, in, out);
546 }
547
548
549
550 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
551 template <int direction,
552 bool dof_to_quad,
553 bool add,
554 bool in_place,
555 typename ViewTypeIn,
556 typename ViewTypeOut>
559 dim,
560 fe_degree,
561 n_q_points_1d,
562 Number>::co_gradients(const ViewTypeIn in,
563 ViewTypeOut out) const
564 {
565 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
566 team_member, co_shape_gradients, in, out);
567 }
568
569
570
571 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
572 template <typename ViewType>
573 DEAL_II_HOST_DEVICE inline void
575 dim,
576 fe_degree,
577 n_q_points_1d,
578 Number>::evaluate_values(ViewType u)
579 {
580 if constexpr (dim == 1)
581 values<0, true, false, true>(u, u);
582 else if constexpr (dim == 2)
583 {
584 values<0, true, false, true>(u, u);
585 team_member.team_barrier();
586 values<1, true, false, true>(u, u);
587 }
588 else if constexpr (dim == 3)
589 {
590 values<0, true, false, true>(u, u);
591 team_member.team_barrier();
592 values<1, true, false, true>(u, u);
593 team_member.team_barrier();
594 values<2, true, false, true>(u, u);
595 }
596 else
597 Kokkos::abort("dim must not exceed 3!");
598 }
599
600
601
602 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
603 template <typename ViewType>
604 DEAL_II_HOST_DEVICE inline void
606 dim,
607 fe_degree,
608 n_q_points_1d,
609 Number>::integrate_values(ViewType u)
610 {
611 if constexpr (dim == 1)
612 values<0, false, false, true>(u, u);
613 else if constexpr (dim == 2)
614 {
615 values<0, false, false, true>(u, u);
616 team_member.team_barrier();
617 values<1, false, false, true>(u, u);
618 }
619 else if constexpr (dim == 3)
620 {
621 values<0, false, false, true>(u, u);
622 team_member.team_barrier();
623 values<1, false, false, true>(u, u);
624 team_member.team_barrier();
625 values<2, false, false, true>(u, u);
626 }
627 else
628 Kokkos::abort("dim must not exceed 3!");
629 }
630
631
632
633 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
634 template <typename ViewTypeIn, typename ViewTypeOut>
635 DEAL_II_HOST_DEVICE inline void
637 dim,
638 fe_degree,
639 n_q_points_1d,
640 Number>::evaluate_gradients(const ViewTypeIn u,
641 ViewTypeOut grad_u)
642 {
643 if constexpr (dim == 1)
644 {
645 gradients<0, true, false, false>(
646 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
647 }
648 else if constexpr (dim == 2)
649 {
650 gradients<0, true, false, false>(
651 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
652 values<0, true, false, false>(
653 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
654
655 team_member.team_barrier();
656
657 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
658 Kokkos::subview(grad_u, Kokkos::ALL, 0));
659 gradients<1, true, false, true>(
660 Kokkos::subview(grad_u, Kokkos::ALL, 1),
661 Kokkos::subview(grad_u, Kokkos::ALL, 1));
662 }
663 else if constexpr (dim == 3)
664 {
665 gradients<0, true, false, false>(
666 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
667 values<0, true, false, false>(
668 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
669 values<0, true, false, false>(
670 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
671
672 team_member.team_barrier();
673
674 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
675 Kokkos::subview(grad_u, Kokkos::ALL, 0));
676 gradients<1, true, false, true>(
677 Kokkos::subview(grad_u, Kokkos::ALL, 1),
678 Kokkos::subview(grad_u, Kokkos::ALL, 1));
679 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
680 Kokkos::subview(grad_u, Kokkos::ALL, 2));
681
682 team_member.team_barrier();
683
684 values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
685 Kokkos::subview(grad_u, Kokkos::ALL, 0));
686 values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
687 Kokkos::subview(grad_u, Kokkos::ALL, 1));
688 gradients<2, true, false, true>(
689 Kokkos::subview(grad_u, Kokkos::ALL, 2),
690 Kokkos::subview(grad_u, Kokkos::ALL, 2));
691 }
692 else
693 Kokkos::abort("dim must not exceed 3!");
694 }
695
696
697
698 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
699 template <typename ViewType1, typename ViewType2>
700 DEAL_II_HOST_DEVICE inline void
702 dim,
703 fe_degree,
704 n_q_points_1d,
705 Number>::evaluate_values_and_gradients(ViewType1 u,
706 ViewType2
707 grad_u)
708 {
709 if constexpr (dim == 1)
710 {
711 values<0, true, false, true>(u, u);
712 team_member.team_barrier();
713
714 co_gradients<0, true, false, false>(
715 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
716 }
717 else if constexpr (dim == 2)
718 {
719 values<0, true, false, true>(u, u);
720 team_member.team_barrier();
721 values<1, true, false, true>(u, u);
722 team_member.team_barrier();
723
724 co_gradients<0, true, false, false>(
725 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
726 co_gradients<1, true, false, false>(
727 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
728 }
729 else if constexpr (dim == 3)
730 {
731 values<0, true, false, true>(u, u);
732 team_member.team_barrier();
733 values<1, true, false, true>(u, u);
734 team_member.team_barrier();
735 values<2, true, false, true>(u, u);
736 team_member.team_barrier();
737
738 co_gradients<0, true, false, false>(
739 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
740 co_gradients<1, true, false, false>(
741 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
742 co_gradients<2, true, false, false>(
743 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
744 }
745 else
746 Kokkos::abort("dim must not exceed 3!");
747 }
748
749
750
751 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
752 template <bool add, typename ViewType1, typename ViewType2>
753 DEAL_II_HOST_DEVICE inline void
755 dim,
756 fe_degree,
757 n_q_points_1d,
758 Number>::integrate_gradients(ViewType1 u,
759 ViewType2 grad_u)
760 {
761 if constexpr (dim == 1)
762 {
763 gradients<0, false, add, false>(
764 Kokkos::subview(grad_u, Kokkos::ALL, dim), u);
765 }
766 else if constexpr (dim == 2)
767 {
768 gradients<0, false, false, true>(
769 Kokkos::subview(grad_u, Kokkos::ALL, 0),
770 Kokkos::subview(grad_u, Kokkos::ALL, 0));
771 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
772 Kokkos::subview(grad_u,
773 Kokkos::ALL,
774 1));
775
776 team_member.team_barrier();
777
778 values<1, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
779 u);
780 team_member.team_barrier();
781 gradients<1, false, true, false>(
782 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
783 }
784 else if constexpr (dim == 3)
785 {
786 gradients<0, false, false, true>(
787 Kokkos::subview(grad_u, Kokkos::ALL, 0),
788 Kokkos::subview(grad_u, Kokkos::ALL, 0));
789 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
790 Kokkos::subview(grad_u,
791 Kokkos::ALL,
792 1));
793 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
794 Kokkos::subview(grad_u,
795 Kokkos::ALL,
796 2));
797
798 team_member.team_barrier();
799
800 values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
801 Kokkos::subview(grad_u,
802 Kokkos::ALL,
803 0));
804 gradients<1, false, false, true>(
805 Kokkos::subview(grad_u, Kokkos::ALL, 1),
806 Kokkos::subview(grad_u, Kokkos::ALL, 1));
807 values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
808 Kokkos::subview(grad_u,
809 Kokkos::ALL,
810 2));
811
812 team_member.team_barrier();
813
814 values<2, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
815 u);
816 team_member.team_barrier();
817 values<2, false, true, false>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
818 u);
819 team_member.team_barrier();
820 gradients<2, false, true, false>(
821 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
822 }
823 else
824 Kokkos::abort("dim must not exceed 3!");
825 }
826
827
828
829 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
830 template <typename ViewType1, typename ViewType2>
831 DEAL_II_HOST_DEVICE inline void
833 dim,
834 fe_degree,
835 n_q_points_1d,
836 Number>::integrate_values_and_gradients(ViewType1 u,
837 ViewType2
838 grad_u)
839 {
840 if constexpr (dim == 1)
841 {
842 co_gradients<0, false, true, false>(
843 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
844 team_member.team_barrier();
845
846 values<0, false, false, true>(u, u);
847 }
848 else if constexpr (dim == 2)
849 {
850 co_gradients<1, false, true, false>(
851 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
852 team_member.team_barrier();
853 co_gradients<0, false, true, false>(
854 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
855 team_member.team_barrier();
856
857 values<1, false, false, true>(u, u);
858 team_member.team_barrier();
859 values<0, false, false, true>(u, u);
860 team_member.team_barrier();
861 }
862 else if constexpr (dim == 3)
863 {
864 co_gradients<2, false, true, false>(
865 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
866 team_member.team_barrier();
867 co_gradients<1, false, true, false>(
868 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
869 team_member.team_barrier();
870 co_gradients<0, false, true, false>(
871 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
872 team_member.team_barrier();
873
874 values<2, false, false, true>(u, u);
875 team_member.team_barrier();
876 values<1, false, false, true>(u, u);
877 team_member.team_barrier();
878 values<0, false, false, true>(u, u);
879 team_member.team_barrier();
880 }
881 else
882 Kokkos::abort("dim must not exceed 3!");
883 }
884 } // namespace internal
885} // namespace Portable
886
888
889#endif
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
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:966
const InputIterator OutputIterator out
Definition parallel.h:167
#define DEAL_II_HOST_DEVICE
Definition numbers.h:34
Kokkos::TeamPolicy< MemorySpace::Default::kokkos_space::execution_space >::member_type TeamHandle