Reference documentation for deal.II version GIT relicensing-487-ge9eb5ab491 2024-04-25 07:20: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
136 // This loop simply multiplies the shape function at the quadrature
137 // point by the value finite element coefficient.
138 // FIXME check why using parallel_reduce ThreadVector is slower
139 const int base_shape = dof_to_quad ? j : j * n_q_points_1d;
140 const int stride_shape = dof_to_quad ? n_q_points_1d : 1;
141 const int base_in = (direction == 0) ? (n_q_points_1d * i) : i;
142 const int stride_in = Utilities::pow(n_q_points_1d, direction);
143 Number sum = shape_data[base_shape] * in(base_in);
144 for (int k = 1; k < n_q_points_1d; ++k)
145 {
146 sum += shape_data[base_shape + k * stride_shape] *
147 in(base_in + k * stride_in);
148 }
149 t[q_point] = sum;
150 });
151
152 if (in_place)
153 team_member.team_barrier();
154
155 Kokkos::parallel_for(thread_policy, [&](const int i, const int j) {
156 const int q_point = i + j * n_q_points_1d;
157 const unsigned int destination_idx =
158 (direction == 0) ? (j + n_q_points_1d * i) : (i + n_q_points_1d * j);
159
160 if (add)
161 Kokkos::atomic_add(&out(destination_idx), t[q_point]);
162 else
163 out(destination_idx) = t[q_point];
164 });
165 }
166
167
168
172 template <int n_q_points_1d,
173 typename Number,
174 int direction,
175 bool dof_to_quad,
176 bool add,
177 bool in_place,
178 typename ViewTypeIn,
179 typename ViewTypeOut>
181 apply_3d(const Kokkos::TeamPolicy<
182 MemorySpace::Default::kokkos_space::execution_space>::member_type
183 &team_member,
184 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
185 shape_data,
186 const ViewTypeIn in,
187 ViewTypeOut out)
188 {
189 using TeamType = Kokkos::TeamPolicy<
190 MemorySpace::Default::kokkos_space::execution_space>::member_type;
191 constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, 3);
192
193 Number t[n_q_points];
194 auto thread_policy = Kokkos::TeamThreadMDRange<Kokkos::Rank<3>, TeamType>(
195 team_member, n_q_points_1d, n_q_points_1d, n_q_points_1d);
196 Kokkos::parallel_for(
197 thread_policy, [&](const int i, const int j, const int q) {
198 const int q_point =
199 i + j * n_q_points_1d + q * n_q_points_1d * n_q_points_1d;
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 const int base_shape = dof_to_quad ? q : q * n_q_points_1d;
205 const int stride_shape = dof_to_quad ? n_q_points_1d : 1;
206 const int base_in =
207 (direction == 0 ?
208 (n_q_points_1d * (i + n_q_points_1d * j)) :
209 (direction == 1 ? (i + n_q_points_1d * n_q_points_1d * j) :
210 (i + n_q_points_1d * j)));
211 const int stride_in = Utilities::pow(n_q_points_1d, direction);
212 Number sum = shape_data[base_shape] * in(base_in);
213 for (int k = 1; k < n_q_points_1d; ++k)
214 {
215 sum += shape_data[base_shape + k * stride_shape] *
216 in(base_in + k * stride_in);
217 }
218 t[q_point] = sum;
219 });
220
221 if (in_place)
222 team_member.team_barrier();
223
224 Kokkos::parallel_for(
225 thread_policy, [&](const int i, const int j, const int q) {
226 const int q_point =
227 i + j * n_q_points_1d + q * n_q_points_1d * n_q_points_1d;
228 const unsigned int destination_idx =
229 (direction == 0) ? (q + n_q_points_1d * (i + n_q_points_1d * j)) :
230 (direction == 1) ? (i + n_q_points_1d * (q + n_q_points_1d * j)) :
231 (i + n_q_points_1d * (j + n_q_points_1d * q));
232
233 if (add)
234 Kokkos::atomic_add(&out(destination_idx), t[q_point]);
235 else
236 out(destination_idx) = t[q_point];
237 });
238 }
239#endif
240
241
242
246 template <int dim,
247 int n_q_points_1d,
248 typename Number,
249 int direction,
250 bool dof_to_quad,
251 bool add,
252 bool in_place,
253 typename ViewTypeIn,
254 typename ViewTypeOut>
256 apply(const Kokkos::TeamPolicy<
257 MemorySpace::Default::kokkos_space::execution_space>::member_type
258 &team_member,
259 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
260 shape_data,
261 const ViewTypeIn in,
262 ViewTypeOut out)
263 {
264#if KOKKOS_VERSION >= 40000
265 if constexpr (dim == 1)
266 apply_1d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
267 team_member, shape_data, in, out);
268 if constexpr (dim == 2)
269 apply_2d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
270 team_member, shape_data, in, out);
271 if constexpr (dim == 3)
272 apply_3d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
273 team_member, shape_data, in, out);
274#else
275 constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim);
276
277 Number t[n_q_points];
278 Kokkos::parallel_for(
279 Kokkos::TeamThreadRange(team_member, n_q_points),
280 [&](const int &q_point) {
281 const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
282 const unsigned int j =
283 (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
284 const unsigned int q =
285 (dim == 1) ? q_point :
286 (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
287 q_point / (n_q_points_1d * n_q_points_1d);
288
289 // This loop simply multiplies the shape function at the quadrature
290 // point by the value finite element coefficient.
291 const int stride_shape = dof_to_quad ? n_q_points_1d : 1;
292 const int stride = Utilities::pow(n_q_points_1d, direction);
293 const int base_shape = dof_to_quad ? q : (q * n_q_points_1d);
294 const int base =
295 (direction == 0) ? (n_q_points_1d * (i + n_q_points_1d * j)) :
296 (direction == 1) ? (i + n_q_points_1d * (n_q_points_1d * j)) :
297 (i + n_q_points_1d * j);
298 Number sum =
299 shape_data[base_shape] * (in_place ? out(base) : in(base));
300 for (int k = 1; k < n_q_points_1d; ++k)
301 {
302 sum +=
303 shape_data[base_shape + k * stride_shape] *
304 (in_place ? out(base + k * stride) : in(base + k * stride));
305 }
306 t[q_point] = sum;
307 });
308
309 if (in_place)
310 team_member.team_barrier();
311
312 Kokkos::parallel_for(
313 Kokkos::TeamThreadRange(team_member, n_q_points),
314 [&](const int &q_point) {
315 const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
316 const unsigned int j =
317 (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
318 const unsigned int q =
319 (dim == 1) ? q_point :
320 (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
321 q_point / (n_q_points_1d * n_q_points_1d);
322
323 const unsigned int destination_idx =
324 (direction == 0) ? (q + n_q_points_1d * (i + n_q_points_1d * j)) :
325 (direction == 1) ? (i + n_q_points_1d * (q + n_q_points_1d * j)) :
326 (i + n_q_points_1d * (j + n_q_points_1d * q));
327
328 if (add)
329 Kokkos::atomic_add(&out(destination_idx), t[q_point]);
330 else
331 out(destination_idx) = t[q_point];
332 });
333#endif
334 }
335
336
343 template <EvaluatorVariant variant,
344 int dim,
345 int fe_degree,
346 int n_q_points_1d,
347 typename Number>
350
351
352
360 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
362 dim,
363 fe_degree,
364 n_q_points_1d,
365 Number>
366 {
367 public:
368 using TeamHandle = Kokkos::TeamPolicy<
369 MemorySpace::Default::kokkos_space::execution_space>::member_type;
370
373 const TeamHandle &team_member,
374 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
375 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
376 shape_gradients,
377 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
378 co_shape_gradients);
379
383 template <typename ViewType>
385 evaluate_values(ViewType u);
386
391 template <typename ViewTypeIn, typename ViewTypeOut>
393 evaluate_gradients(const ViewTypeIn u, ViewTypeOut grad_u);
394
399 template <typename ViewType1, typename ViewType2>
401 evaluate_values_and_gradients(ViewType1 u, ViewType2 grad_u);
402
406 template <typename ViewType>
408 integrate_values(ViewType u);
409
414 template <bool add, typename ViewType1, typename ViewType2>
416 integrate_gradients(ViewType1 u, ViewType2 grad_u);
417
422 template <typename ViewType1, typename ViewType2>
424 integrate_values_and_gradients(ViewType1 u, ViewType2 grad_u);
425
430 template <int direction,
431 bool dof_to_quad,
432 bool add,
433 bool in_place,
434 typename ViewTypeIn,
435 typename ViewTypeOut>
437 values(const ViewTypeIn in, ViewTypeOut out) const;
438
443 template <int direction,
444 bool dof_to_quad,
445 bool add,
446 bool in_place,
447 typename ViewTypeIn,
448 typename ViewTypeOut>
450 gradients(const ViewTypeIn in, ViewTypeOut out) const;
451
452 public:
457 template <int direction,
458 bool dof_to_quad,
459 bool add,
460 bool in_place,
461 typename ViewTypeIn,
462 typename ViewTypeOut>
464 co_gradients(const ViewTypeIn in, ViewTypeOut out) const;
465
470
474 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values;
475
479 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
481
485 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
487 };
488
489
490
491 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
494 dim,
495 fe_degree,
496 n_q_points_1d,
497 Number>::
498 EvaluatorTensorProduct(
499 const TeamHandle &team_member,
500 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
501 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
502 shape_gradients,
503 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
504 co_shape_gradients)
505 : team_member(team_member)
506 , shape_values(shape_values)
507 , shape_gradients(shape_gradients)
508 , co_shape_gradients(co_shape_gradients)
509 {}
510
511
512
513 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
514 template <int direction,
515 bool dof_to_quad,
516 bool add,
517 bool in_place,
518 typename ViewTypeIn,
519 typename ViewTypeOut>
522 dim,
523 fe_degree,
524 n_q_points_1d,
525 Number>::values(const ViewTypeIn in,
526 ViewTypeOut out) const
527 {
528 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
529 team_member, shape_values, in, out);
530 }
531
532
533
534 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
535 template <int direction,
536 bool dof_to_quad,
537 bool add,
538 bool in_place,
539 typename ViewTypeIn,
540 typename ViewTypeOut>
543 dim,
544 fe_degree,
545 n_q_points_1d,
546 Number>::gradients(const ViewTypeIn in,
547 ViewTypeOut out) const
548 {
549 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
550 team_member, shape_gradients, in, out);
551 }
552
553
554
555 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
556 template <int direction,
557 bool dof_to_quad,
558 bool add,
559 bool in_place,
560 typename ViewTypeIn,
561 typename ViewTypeOut>
564 dim,
565 fe_degree,
566 n_q_points_1d,
567 Number>::co_gradients(const ViewTypeIn in,
568 ViewTypeOut out) const
569 {
570 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
571 team_member, co_shape_gradients, in, out);
572 }
573
574
575
576 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
577 template <typename ViewType>
578 DEAL_II_HOST_DEVICE inline void
580 dim,
581 fe_degree,
582 n_q_points_1d,
583 Number>::evaluate_values(ViewType u)
584 {
585 if constexpr (dim == 1)
586 values<0, true, false, true>(u, u);
587 else if constexpr (dim == 2)
588 {
589 values<0, true, false, true>(u, u);
590 team_member.team_barrier();
591 values<1, true, false, true>(u, u);
592 }
593 else if constexpr (dim == 3)
594 {
595 values<0, true, false, true>(u, u);
596 team_member.team_barrier();
597 values<1, true, false, true>(u, u);
598 team_member.team_barrier();
599 values<2, true, false, true>(u, u);
600 }
601 else
602 Kokkos::abort("dim must not exceed 3!");
603 }
604
605
606
607 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
608 template <typename ViewType>
609 DEAL_II_HOST_DEVICE inline void
611 dim,
612 fe_degree,
613 n_q_points_1d,
614 Number>::integrate_values(ViewType u)
615 {
616 if constexpr (dim == 1)
617 values<0, false, false, true>(u, u);
618 else if constexpr (dim == 2)
619 {
620 values<0, false, false, true>(u, u);
621 team_member.team_barrier();
622 values<1, false, false, true>(u, u);
623 }
624 else if constexpr (dim == 3)
625 {
626 values<0, false, false, true>(u, u);
627 team_member.team_barrier();
628 values<1, false, false, true>(u, u);
629 team_member.team_barrier();
630 values<2, false, false, true>(u, u);
631 }
632 else
633 Kokkos::abort("dim must not exceed 3!");
634 }
635
636
637
638 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
639 template <typename ViewTypeIn, typename ViewTypeOut>
640 DEAL_II_HOST_DEVICE inline void
642 dim,
643 fe_degree,
644 n_q_points_1d,
645 Number>::evaluate_gradients(const ViewTypeIn u,
646 ViewTypeOut grad_u)
647 {
648 if constexpr (dim == 1)
649 {
650 gradients<0, true, false, false>(
651 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
652 }
653 else if constexpr (dim == 2)
654 {
655 gradients<0, true, false, false>(
656 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
657 values<0, true, false, false>(
658 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
659
660 team_member.team_barrier();
661
662 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
663 Kokkos::subview(grad_u, Kokkos::ALL, 0));
664 gradients<1, true, false, true>(
665 Kokkos::subview(grad_u, Kokkos::ALL, 1),
666 Kokkos::subview(grad_u, Kokkos::ALL, 1));
667 }
668 else if constexpr (dim == 3)
669 {
670 gradients<0, true, false, false>(
671 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
672 values<0, true, false, false>(
673 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
674 values<0, true, false, false>(
675 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
676
677 team_member.team_barrier();
678
679 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
680 Kokkos::subview(grad_u, Kokkos::ALL, 0));
681 gradients<1, true, false, true>(
682 Kokkos::subview(grad_u, Kokkos::ALL, 1),
683 Kokkos::subview(grad_u, Kokkos::ALL, 1));
684 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
685 Kokkos::subview(grad_u, Kokkos::ALL, 2));
686
687 team_member.team_barrier();
688
689 values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
690 Kokkos::subview(grad_u, Kokkos::ALL, 0));
691 values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
692 Kokkos::subview(grad_u, Kokkos::ALL, 1));
693 gradients<2, true, false, true>(
694 Kokkos::subview(grad_u, Kokkos::ALL, 2),
695 Kokkos::subview(grad_u, Kokkos::ALL, 2));
696 }
697 else
698 Kokkos::abort("dim must not exceed 3!");
699 }
700
701
702
703 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
704 template <typename ViewType1, typename ViewType2>
705 DEAL_II_HOST_DEVICE inline void
707 dim,
708 fe_degree,
709 n_q_points_1d,
710 Number>::evaluate_values_and_gradients(ViewType1 u,
711 ViewType2
712 grad_u)
713 {
714 if constexpr (dim == 1)
715 {
716 values<0, true, false, true>(u, u);
717 team_member.team_barrier();
718
719 co_gradients<0, true, false, false>(
720 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
721 }
722 else if constexpr (dim == 2)
723 {
724 values<0, true, false, true>(u, u);
725 team_member.team_barrier();
726 values<1, true, false, true>(u, u);
727 team_member.team_barrier();
728
729 co_gradients<0, true, false, false>(
730 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
731 co_gradients<1, true, false, false>(
732 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
733 }
734 else if constexpr (dim == 3)
735 {
736 values<0, true, false, true>(u, u);
737 team_member.team_barrier();
738 values<1, true, false, true>(u, u);
739 team_member.team_barrier();
740 values<2, true, false, true>(u, u);
741 team_member.team_barrier();
742
743 co_gradients<0, true, false, false>(
744 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
745 co_gradients<1, true, false, false>(
746 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
747 co_gradients<2, true, false, false>(
748 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
749 }
750 else
751 Kokkos::abort("dim must not exceed 3!");
752 }
753
754
755
756 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
757 template <bool add, typename ViewType1, typename ViewType2>
758 DEAL_II_HOST_DEVICE inline void
760 dim,
761 fe_degree,
762 n_q_points_1d,
763 Number>::integrate_gradients(ViewType1 u,
764 ViewType2 grad_u)
765 {
766 if constexpr (dim == 1)
767 {
768 gradients<0, false, add, false>(
769 Kokkos::subview(grad_u, Kokkos::ALL, dim), u);
770 }
771 else if constexpr (dim == 2)
772 {
773 gradients<0, false, false, true>(
774 Kokkos::subview(grad_u, Kokkos::ALL, 0),
775 Kokkos::subview(grad_u, Kokkos::ALL, 0));
776 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
777 Kokkos::subview(grad_u,
778 Kokkos::ALL,
779 1));
780
781 team_member.team_barrier();
782
783 values<1, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
784 u);
785 team_member.team_barrier();
786 gradients<1, false, true, false>(
787 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
788 }
789 else if constexpr (dim == 3)
790 {
791 gradients<0, false, false, true>(
792 Kokkos::subview(grad_u, Kokkos::ALL, 0),
793 Kokkos::subview(grad_u, Kokkos::ALL, 0));
794 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
795 Kokkos::subview(grad_u,
796 Kokkos::ALL,
797 1));
798 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
799 Kokkos::subview(grad_u,
800 Kokkos::ALL,
801 2));
802
803 team_member.team_barrier();
804
805 values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
806 Kokkos::subview(grad_u,
807 Kokkos::ALL,
808 0));
809 gradients<1, false, false, true>(
810 Kokkos::subview(grad_u, Kokkos::ALL, 1),
811 Kokkos::subview(grad_u, Kokkos::ALL, 1));
812 values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
813 Kokkos::subview(grad_u,
814 Kokkos::ALL,
815 2));
816
817 team_member.team_barrier();
818
819 values<2, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
820 u);
821 team_member.team_barrier();
822 values<2, false, true, false>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
823 u);
824 team_member.team_barrier();
825 gradients<2, false, true, false>(
826 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
827 }
828 else
829 Kokkos::abort("dim must not exceed 3!");
830 }
831
832
833
834 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
835 template <typename ViewType1, typename ViewType2>
836 DEAL_II_HOST_DEVICE inline void
838 dim,
839 fe_degree,
840 n_q_points_1d,
841 Number>::integrate_values_and_gradients(ViewType1 u,
842 ViewType2
843 grad_u)
844 {
845 if constexpr (dim == 1)
846 {
847 co_gradients<0, false, true, false>(
848 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
849 team_member.team_barrier();
850
851 values<0, false, false, true>(u, u);
852 }
853 else if constexpr (dim == 2)
854 {
855 co_gradients<1, false, true, false>(
856 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
857 team_member.team_barrier();
858 co_gradients<0, false, true, false>(
859 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
860 team_member.team_barrier();
861
862 values<1, false, false, true>(u, u);
863 team_member.team_barrier();
864 values<0, false, false, true>(u, u);
865 team_member.team_barrier();
866 }
867 else if constexpr (dim == 3)
868 {
869 co_gradients<2, false, true, false>(
870 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
871 team_member.team_barrier();
872 co_gradients<1, false, true, false>(
873 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
874 team_member.team_barrier();
875 co_gradients<0, false, true, false>(
876 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
877 team_member.team_barrier();
878
879 values<2, false, false, true>(u, u);
880 team_member.team_barrier();
881 values<1, false, false, true>(u, u);
882 team_member.team_barrier();
883 values<0, false, false, true>(u, u);
884 team_member.team_barrier();
885 }
886 else
887 Kokkos::abort("dim must not exceed 3!");
888 }
889 } // namespace internal
890} // namespace Portable
891
893
894#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)
T sum(const T &t, const MPI_Comm mpi_communicator)
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