deal.II version GIT relicensing-2167-g9622207b8f 2024-11-21 12:40:00+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
21#include <deal.II/base/utilities.h>
22
23DEAL_II_NAMESPACE_OPEN
24
25
26namespace Portable
27{
28 namespace internal
29 {
34 // TODO: for now only the general variant is implemented
41
42
43
44#if KOKKOS_VERSION >= 40000
48 template <int n_q_points_1d,
49 typename Number,
50 int direction,
51 bool dof_to_quad,
52 bool add,
53 bool in_place,
54 typename ViewTypeIn,
55 typename ViewTypeOut>
56 DEAL_II_HOST_DEVICE void
57 apply_1d(const Kokkos::TeamPolicy<
58 MemorySpace::Default::kokkos_space::execution_space>::member_type
59 &team_member,
60 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
61 shape_data,
62 const ViewTypeIn in,
63 ViewTypeOut out)
64 {
65 Number t[n_q_points_1d];
66 Kokkos::parallel_for(Kokkos::TeamThreadRange(team_member, n_q_points_1d),
67 [&](const int &q) {
68 t[q] = 0;
69 // This loop simply multiplies the shape function
70 // at the quadrature point by the value finite
71 // element coefficient.
72 // FIXME check why using parallel_reduce
73 // ThreadVector is slower
74 for (int k = 0; k < n_q_points_1d; ++k)
75 {
76 const unsigned int shape_idx =
77 dof_to_quad ? (q + k * n_q_points_1d) :
78 (k + q * n_q_points_1d);
79 const unsigned int source_idx = k;
80 t[q] += shape_data[shape_idx] * in(source_idx);
81 }
82 });
83
84 if constexpr (in_place)
85 team_member.team_barrier();
86
87 Kokkos::parallel_for(Kokkos::TeamThreadRange(team_member, n_q_points_1d),
88 [&](const int &q) {
89 const unsigned int destination_idx = q;
90 if constexpr (add)
91 Kokkos::atomic_add(&out(destination_idx), t[q]);
92 else
93 out(destination_idx) = t[q];
94 });
95 }
96
97
98
102 template <int n_q_points_1d,
103 typename Number,
104 int direction,
105 bool dof_to_quad,
106 bool add,
107 bool in_place,
108 typename ViewTypeIn,
109 typename ViewTypeOut>
110 DEAL_II_HOST_DEVICE void
111 apply_2d(const Kokkos::TeamPolicy<
112 MemorySpace::Default::kokkos_space::execution_space>::member_type
113 &team_member,
114 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
115 shape_data,
116 const ViewTypeIn in,
117 ViewTypeOut out)
118 {
119 using TeamType = Kokkos::TeamPolicy<
120 MemorySpace::Default::kokkos_space::execution_space>::member_type;
121 constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, 2);
122
123 Number t[n_q_points];
124 auto thread_policy =
125 Kokkos::TeamThreadMDRange<Kokkos::Rank<2>, TeamType>(team_member,
126 n_q_points_1d,
127 n_q_points_1d);
128 Kokkos::parallel_for(thread_policy, [&](const int i, const int j) {
129 int q_point = i + j * n_q_points_1d;
130
131 // This loop simply multiplies the shape function at the quadrature
132 // point by the value finite element coefficient.
133 // FIXME check why using parallel_reduce ThreadVector is slower
134 const int base_shape = dof_to_quad ? j : j * n_q_points_1d;
135 const int stride_shape = dof_to_quad ? n_q_points_1d : 1;
136 const int base_in = (direction == 0) ? (n_q_points_1d * i) : i;
137 const int stride_in = Utilities::pow(n_q_points_1d, direction);
138 Number sum = shape_data[base_shape] * in(base_in);
139 for (int k = 1; k < n_q_points_1d; ++k)
140 {
141 sum += shape_data[base_shape + k * stride_shape] *
142 in(base_in + k * stride_in);
143 }
144 t[q_point] = sum;
145 });
146
147 if (in_place)
148 team_member.team_barrier();
149
150 Kokkos::parallel_for(thread_policy, [&](const int i, const int j) {
151 const int q_point = i + j * n_q_points_1d;
152 const unsigned int destination_idx =
153 (direction == 0) ? (j + n_q_points_1d * i) : (i + n_q_points_1d * j);
154
155 if (add)
156 Kokkos::atomic_add(&out(destination_idx), t[q_point]);
157 else
158 out(destination_idx) = t[q_point];
159 });
160 }
161
162
163
167 template <int n_q_points_1d,
168 typename Number,
169 int direction,
170 bool dof_to_quad,
171 bool add,
172 bool in_place,
173 typename ViewTypeIn,
174 typename ViewTypeOut>
175 DEAL_II_HOST_DEVICE void
176 apply_3d(const Kokkos::TeamPolicy<
177 MemorySpace::Default::kokkos_space::execution_space>::member_type
178 &team_member,
179 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
180 shape_data,
181 const ViewTypeIn in,
182 ViewTypeOut out)
183 {
184 using TeamType = Kokkos::TeamPolicy<
185 MemorySpace::Default::kokkos_space::execution_space>::member_type;
186 constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, 3);
187
188 Number t[n_q_points];
189 auto thread_policy = Kokkos::TeamThreadMDRange<Kokkos::Rank<3>, TeamType>(
190 team_member, n_q_points_1d, n_q_points_1d, n_q_points_1d);
191 Kokkos::parallel_for(
192 thread_policy, [&](const int i, const int j, const int q) {
193 const int q_point =
194 i + j * n_q_points_1d + q * n_q_points_1d * n_q_points_1d;
195
196 // This loop simply multiplies the shape function at the quadrature
197 // point by the value finite element coefficient.
198 // FIXME check why using parallel_reduce ThreadVector is slower
199 const int base_shape = dof_to_quad ? q : q * n_q_points_1d;
200 const int stride_shape = dof_to_quad ? n_q_points_1d : 1;
201 const int base_in =
202 (direction == 0 ?
203 (n_q_points_1d * (i + n_q_points_1d * j)) :
204 (direction == 1 ? (i + n_q_points_1d * n_q_points_1d * j) :
205 (i + n_q_points_1d * j)));
206 const int stride_in = Utilities::pow(n_q_points_1d, direction);
207 Number sum = shape_data[base_shape] * in(base_in);
208 for (int k = 1; k < n_q_points_1d; ++k)
209 {
210 sum += shape_data[base_shape + k * stride_shape] *
211 in(base_in + k * stride_in);
212 }
213 t[q_point] = sum;
214 });
215
216 if (in_place)
217 team_member.team_barrier();
218
219 Kokkos::parallel_for(
220 thread_policy, [&](const int i, const int j, const int q) {
221 const int q_point =
222 i + j * n_q_points_1d + q * n_q_points_1d * n_q_points_1d;
223 const unsigned int destination_idx =
224 (direction == 0) ? (q + n_q_points_1d * (i + n_q_points_1d * j)) :
225 (direction == 1) ? (i + n_q_points_1d * (q + n_q_points_1d * j)) :
226 (i + n_q_points_1d * (j + n_q_points_1d * q));
227
228 if (add)
229 Kokkos::atomic_add(&out(destination_idx), t[q_point]);
230 else
231 out(destination_idx) = t[q_point];
232 });
233 }
234#endif
235
236
237
241 template <int dim,
242 int n_q_points_1d,
243 typename Number,
244 int direction,
245 bool dof_to_quad,
246 bool add,
247 bool in_place,
248 typename ViewTypeIn,
249 typename ViewTypeOut>
250 DEAL_II_HOST_DEVICE void
251 apply(const Kokkos::TeamPolicy<
252 MemorySpace::Default::kokkos_space::execution_space>::member_type
253 &team_member,
254 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
255 shape_data,
256 const ViewTypeIn in,
257 ViewTypeOut out)
258 {
259#if KOKKOS_VERSION >= 40000
260 if constexpr (dim == 1)
261 apply_1d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
262 team_member, shape_data, in, out);
263 if constexpr (dim == 2)
264 apply_2d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
265 team_member, shape_data, in, out);
266 if constexpr (dim == 3)
267 apply_3d<n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
268 team_member, shape_data, in, out);
269#else
270 constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim);
271
272 Number t[n_q_points];
273 Kokkos::parallel_for(
274 Kokkos::TeamThreadRange(team_member, n_q_points),
275 [&](const int &q_point) {
276 const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
277 const unsigned int j =
278 (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
279 const unsigned int q =
280 (dim == 1) ? q_point :
281 (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
282 q_point / (n_q_points_1d * n_q_points_1d);
283
284 // This loop simply multiplies the shape function at the quadrature
285 // point by the value finite element coefficient.
286 const int stride_shape = dof_to_quad ? n_q_points_1d : 1;
287 const int stride = Utilities::pow(n_q_points_1d, direction);
288 const int base_shape = dof_to_quad ? q : (q * n_q_points_1d);
289 const int base =
290 (direction == 0) ? (n_q_points_1d * (i + n_q_points_1d * j)) :
291 (direction == 1) ? (i + n_q_points_1d * (n_q_points_1d * j)) :
292 (i + n_q_points_1d * j);
293 Number sum =
294 shape_data[base_shape] * (in_place ? out(base) : in(base));
295 for (int k = 1; k < n_q_points_1d; ++k)
296 {
297 sum +=
298 shape_data[base_shape + k * stride_shape] *
299 (in_place ? out(base + k * stride) : in(base + k * stride));
300 }
301 t[q_point] = sum;
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
335 template <EvaluatorVariant variant,
336 int dim,
337 int fe_degree,
338 int n_q_points_1d,
339 typename Number>
342
343
344
349 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
350 struct EvaluatorTensorProduct<evaluate_general,
351 dim,
352 fe_degree,
353 n_q_points_1d,
354 Number>
355 {
356 public:
357 using TeamHandle = Kokkos::TeamPolicy<
358 MemorySpace::Default::kokkos_space::execution_space>::member_type;
359
360 DEAL_II_HOST_DEVICE
361 EvaluatorTensorProduct(
362 const TeamHandle &team_member,
363 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
364 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
365 shape_gradients,
366 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
367 co_shape_gradients);
368
372 template <typename ViewType>
373 DEAL_II_HOST_DEVICE void
374 evaluate_values(ViewType u);
375
380 template <typename ViewTypeIn, typename ViewTypeOut>
381 DEAL_II_HOST_DEVICE void
382 evaluate_gradients(const ViewTypeIn u, ViewTypeOut grad_u);
383
388 template <typename ViewType1, typename ViewType2>
389 DEAL_II_HOST_DEVICE void
390 evaluate_values_and_gradients(ViewType1 u, ViewType2 grad_u);
391
395 template <typename ViewType>
396 DEAL_II_HOST_DEVICE void
397 integrate_values(ViewType u);
398
403 template <bool add, typename ViewType1, typename ViewType2>
404 DEAL_II_HOST_DEVICE void
405 integrate_gradients(ViewType1 u, ViewType2 grad_u);
406
411 template <typename ViewType1, typename ViewType2>
412 DEAL_II_HOST_DEVICE void
413 integrate_values_and_gradients(ViewType1 u, ViewType2 grad_u);
414
419 template <int direction,
420 bool dof_to_quad,
421 bool add,
422 bool in_place,
423 typename ViewTypeIn,
424 typename ViewTypeOut>
425 DEAL_II_HOST_DEVICE void
426 values(const ViewTypeIn in, ViewTypeOut out) const;
427
432 template <int direction,
433 bool dof_to_quad,
434 bool add,
435 bool in_place,
436 typename ViewTypeIn,
437 typename ViewTypeOut>
438 DEAL_II_HOST_DEVICE void
439 gradients(const ViewTypeIn in, ViewTypeOut out) const;
440
441 public:
446 template <int direction,
447 bool dof_to_quad,
448 bool add,
449 bool in_place,
450 typename ViewTypeIn,
451 typename ViewTypeOut>
452 DEAL_II_HOST_DEVICE void
453 co_gradients(const ViewTypeIn in, ViewTypeOut out) const;
454
458 const TeamHandle &team_member;
459
463 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values;
464
468 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
469 shape_gradients;
470
474 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
475 co_shape_gradients;
476 };
477
478
479
480 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
481 DEAL_II_HOST_DEVICE
482 EvaluatorTensorProduct<evaluate_general,
483 dim,
484 fe_degree,
485 n_q_points_1d,
486 Number>::
487 EvaluatorTensorProduct(
488 const TeamHandle &team_member,
489 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
490 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
491 shape_gradients,
492 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
493 co_shape_gradients)
494 : team_member(team_member)
495 , shape_values(shape_values)
496 , shape_gradients(shape_gradients)
497 , co_shape_gradients(co_shape_gradients)
498 {}
499
500
501
502 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
503 template <int direction,
504 bool dof_to_quad,
505 bool add,
506 bool in_place,
507 typename ViewTypeIn,
508 typename ViewTypeOut>
509 DEAL_II_HOST_DEVICE void
510 EvaluatorTensorProduct<evaluate_general,
511 dim,
512 fe_degree,
513 n_q_points_1d,
514 Number>::values(const ViewTypeIn in,
515 ViewTypeOut out) const
516 {
517 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
518 team_member, shape_values, in, out);
519 }
520
521
522
523 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
524 template <int direction,
525 bool dof_to_quad,
526 bool add,
527 bool in_place,
528 typename ViewTypeIn,
529 typename ViewTypeOut>
530 DEAL_II_HOST_DEVICE void
531 EvaluatorTensorProduct<evaluate_general,
532 dim,
533 fe_degree,
534 n_q_points_1d,
535 Number>::gradients(const ViewTypeIn in,
536 ViewTypeOut out) const
537 {
538 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
539 team_member, shape_gradients, in, out);
540 }
541
542
543
544 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
545 template <int direction,
546 bool dof_to_quad,
547 bool add,
548 bool in_place,
549 typename ViewTypeIn,
550 typename ViewTypeOut>
551 DEAL_II_HOST_DEVICE void
552 EvaluatorTensorProduct<evaluate_general,
553 dim,
554 fe_degree,
555 n_q_points_1d,
556 Number>::co_gradients(const ViewTypeIn in,
557 ViewTypeOut out) const
558 {
559 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
560 team_member, co_shape_gradients, in, out);
561 }
562
563
564
565 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
566 template <typename ViewType>
567 DEAL_II_HOST_DEVICE inline void
568 EvaluatorTensorProduct<evaluate_general,
569 dim,
570 fe_degree,
571 n_q_points_1d,
572 Number>::evaluate_values(ViewType u)
573 {
574 if constexpr (dim == 1)
575 values<0, true, false, true>(u, u);
576 else if constexpr (dim == 2)
577 {
578 values<0, true, false, true>(u, u);
579 team_member.team_barrier();
580 values<1, true, false, true>(u, u);
581 }
582 else if constexpr (dim == 3)
583 {
584 values<0, true, false, true>(u, u);
585 team_member.team_barrier();
586 values<1, true, false, true>(u, u);
587 team_member.team_barrier();
588 values<2, true, false, true>(u, u);
589 }
590 else
591 Kokkos::abort("dim must not exceed 3!");
592 }
593
594
595
596 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
597 template <typename ViewType>
598 DEAL_II_HOST_DEVICE inline void
599 EvaluatorTensorProduct<evaluate_general,
600 dim,
601 fe_degree,
602 n_q_points_1d,
603 Number>::integrate_values(ViewType u)
604 {
605 if constexpr (dim == 1)
606 values<0, false, false, true>(u, u);
607 else if constexpr (dim == 2)
608 {
609 values<0, false, false, true>(u, u);
610 team_member.team_barrier();
611 values<1, false, false, true>(u, u);
612 }
613 else if constexpr (dim == 3)
614 {
615 values<0, false, false, true>(u, u);
616 team_member.team_barrier();
617 values<1, false, false, true>(u, u);
618 team_member.team_barrier();
619 values<2, false, false, true>(u, u);
620 }
621 else
622 Kokkos::abort("dim must not exceed 3!");
623 }
624
625
626
627 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
628 template <typename ViewTypeIn, typename ViewTypeOut>
629 DEAL_II_HOST_DEVICE inline void
630 EvaluatorTensorProduct<evaluate_general,
631 dim,
632 fe_degree,
633 n_q_points_1d,
634 Number>::evaluate_gradients(const ViewTypeIn u,
635 ViewTypeOut grad_u)
636 {
637 if constexpr (dim == 1)
638 {
639 gradients<0, true, false, false>(
640 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
641 }
642 else if constexpr (dim == 2)
643 {
644 gradients<0, true, false, false>(
645 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
646 values<0, true, false, false>(
647 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
648
649 team_member.team_barrier();
650
651 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
652 Kokkos::subview(grad_u, Kokkos::ALL, 0));
653 gradients<1, true, false, true>(
654 Kokkos::subview(grad_u, Kokkos::ALL, 1),
655 Kokkos::subview(grad_u, Kokkos::ALL, 1));
656 }
657 else if constexpr (dim == 3)
658 {
659 gradients<0, true, false, false>(
660 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
661 values<0, true, false, false>(
662 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
663 values<0, true, false, false>(
664 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
665
666 team_member.team_barrier();
667
668 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
669 Kokkos::subview(grad_u, Kokkos::ALL, 0));
670 gradients<1, true, false, true>(
671 Kokkos::subview(grad_u, Kokkos::ALL, 1),
672 Kokkos::subview(grad_u, Kokkos::ALL, 1));
673 values<1, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
674 Kokkos::subview(grad_u, Kokkos::ALL, 2));
675
676 team_member.team_barrier();
677
678 values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
679 Kokkos::subview(grad_u, Kokkos::ALL, 0));
680 values<2, true, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
681 Kokkos::subview(grad_u, Kokkos::ALL, 1));
682 gradients<2, true, false, true>(
683 Kokkos::subview(grad_u, Kokkos::ALL, 2),
684 Kokkos::subview(grad_u, Kokkos::ALL, 2));
685 }
686 else
687 Kokkos::abort("dim must not exceed 3!");
688 }
689
690
691
692 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
693 template <typename ViewType1, typename ViewType2>
694 DEAL_II_HOST_DEVICE inline void
695 EvaluatorTensorProduct<evaluate_general,
696 dim,
697 fe_degree,
698 n_q_points_1d,
699 Number>::evaluate_values_and_gradients(ViewType1 u,
700 ViewType2
701 grad_u)
702 {
703 if constexpr (dim == 1)
704 {
705 values<0, true, false, true>(u, u);
706 team_member.team_barrier();
707
708 co_gradients<0, true, false, false>(
709 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
710 }
711 else if constexpr (dim == 2)
712 {
713 values<0, true, false, true>(u, u);
714 team_member.team_barrier();
715 values<1, true, false, true>(u, u);
716 team_member.team_barrier();
717
718 co_gradients<0, true, false, false>(
719 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
720 co_gradients<1, true, false, false>(
721 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
722 }
723 else if constexpr (dim == 3)
724 {
725 values<0, true, false, true>(u, u);
726 team_member.team_barrier();
727 values<1, true, false, true>(u, u);
728 team_member.team_barrier();
729 values<2, true, false, true>(u, u);
730 team_member.team_barrier();
731
732 co_gradients<0, true, false, false>(
733 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
734 co_gradients<1, true, false, false>(
735 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
736 co_gradients<2, true, false, false>(
737 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
738 }
739 else
740 Kokkos::abort("dim must not exceed 3!");
741 }
742
743
744
745 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
746 template <bool add, typename ViewType1, typename ViewType2>
747 DEAL_II_HOST_DEVICE inline void
748 EvaluatorTensorProduct<evaluate_general,
749 dim,
750 fe_degree,
751 n_q_points_1d,
752 Number>::integrate_gradients(ViewType1 u,
753 ViewType2 grad_u)
754 {
755 if constexpr (dim == 1)
756 {
757 gradients<0, false, add, false>(
758 Kokkos::subview(grad_u, Kokkos::ALL, dim), u);
759 }
760 else if constexpr (dim == 2)
761 {
762 gradients<0, false, false, true>(
763 Kokkos::subview(grad_u, Kokkos::ALL, 0),
764 Kokkos::subview(grad_u, Kokkos::ALL, 0));
765 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
766 Kokkos::subview(grad_u,
767 Kokkos::ALL,
768 1));
769
770 team_member.team_barrier();
771
772 values<1, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
773 u);
774 team_member.team_barrier();
775 gradients<1, false, true, false>(
776 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
777 }
778 else if constexpr (dim == 3)
779 {
780 gradients<0, false, false, true>(
781 Kokkos::subview(grad_u, Kokkos::ALL, 0),
782 Kokkos::subview(grad_u, Kokkos::ALL, 0));
783 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
784 Kokkos::subview(grad_u,
785 Kokkos::ALL,
786 1));
787 values<0, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
788 Kokkos::subview(grad_u,
789 Kokkos::ALL,
790 2));
791
792 team_member.team_barrier();
793
794 values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
795 Kokkos::subview(grad_u,
796 Kokkos::ALL,
797 0));
798 gradients<1, false, false, true>(
799 Kokkos::subview(grad_u, Kokkos::ALL, 1),
800 Kokkos::subview(grad_u, Kokkos::ALL, 1));
801 values<1, false, false, true>(Kokkos::subview(grad_u, Kokkos::ALL, 2),
802 Kokkos::subview(grad_u,
803 Kokkos::ALL,
804 2));
805
806 team_member.team_barrier();
807
808 values<2, false, add, false>(Kokkos::subview(grad_u, Kokkos::ALL, 0),
809 u);
810 team_member.team_barrier();
811 values<2, false, true, false>(Kokkos::subview(grad_u, Kokkos::ALL, 1),
812 u);
813 team_member.team_barrier();
814 gradients<2, false, true, false>(
815 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
816 }
817 else
818 Kokkos::abort("dim must not exceed 3!");
819 }
820
821
822
823 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
824 template <typename ViewType1, typename ViewType2>
825 DEAL_II_HOST_DEVICE inline void
826 EvaluatorTensorProduct<evaluate_general,
827 dim,
828 fe_degree,
829 n_q_points_1d,
830 Number>::integrate_values_and_gradients(ViewType1 u,
831 ViewType2
832 grad_u)
833 {
834 if constexpr (dim == 1)
835 {
836 co_gradients<0, false, true, false>(
837 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
838 team_member.team_barrier();
839
840 values<0, false, false, true>(u, u);
841 }
842 else if constexpr (dim == 2)
843 {
844 co_gradients<1, false, true, false>(
845 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
846 team_member.team_barrier();
847 co_gradients<0, false, true, false>(
848 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
849 team_member.team_barrier();
850
851 values<1, false, false, true>(u, u);
852 team_member.team_barrier();
853 values<0, false, false, true>(u, u);
854 team_member.team_barrier();
855 }
856 else if constexpr (dim == 3)
857 {
858 co_gradients<2, false, true, false>(
859 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
860 team_member.team_barrier();
861 co_gradients<1, false, true, false>(
862 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
863 team_member.team_barrier();
864 co_gradients<0, false, true, false>(
865 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
866 team_member.team_barrier();
867
868 values<2, false, false, true>(u, u);
869 team_member.team_barrier();
870 values<1, false, false, true>(u, u);
871 team_member.team_barrier();
872 values<0, false, false, true>(u, u);
873 team_member.team_barrier();
874 }
875 else
876 Kokkos::abort("dim must not exceed 3!");
877 }
878 } // namespace internal
879} // namespace Portable
880
881DEAL_II_NAMESPACE_CLOSE
882
883#endif
DEAL_II_NAMESPACE_OPEN
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
DEAL_II_NAMESPACE_CLOSE
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
Portable::internal::apply
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)
Definition portable_tensor_product_kernels.h:251
Utilities::MPI::sum
T sum(const T &t, const MPI_Comm mpi_communicator)
Utilities::pow
constexpr T pow(const T base, const int iexp)
Definition utilities.h:966
DEAL_II_HOST_DEVICE
#define DEAL_II_HOST_DEVICE
Definition numbers.h:30
Portable::internal::EvaluatorTensorProduct< evaluate_general, dim, fe_degree, n_q_points_1d, Number >::team_member
const TeamHandle & team_member
Definition portable_tensor_product_kernels.h:458
Portable::internal::EvaluatorTensorProduct< evaluate_general, dim, fe_degree, n_q_points_1d, Number >::co_shape_gradients
Kokkos::View< Number *, MemorySpace::Default::kokkos_space > co_shape_gradients
Definition portable_tensor_product_kernels.h:475
Portable::internal::EvaluatorTensorProduct< evaluate_general, dim, fe_degree, n_q_points_1d, Number >::shape_values
Kokkos::View< Number *, MemorySpace::Default::kokkos_space > shape_values
Definition portable_tensor_product_kernels.h:463
Portable::internal::EvaluatorTensorProduct< evaluate_general, dim, fe_degree, n_q_points_1d, Number >::shape_gradients
Kokkos::View< Number *, MemorySpace::Default::kokkos_space > shape_gradients
Definition portable_tensor_product_kernels.h:469
Portable::internal::EvaluatorTensorProduct< evaluate_general, dim, fe_degree, n_q_points_1d, Number >::TeamHandle
Kokkos::TeamPolicy< MemorySpace::Default::kokkos_space::execution_space >::member_type TeamHandle
Definition portable_tensor_product_kernels.h:358
Generated by doxygen 1.9.8