Reference documentation for deal.II version 9.5.0
\(\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
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
23
24#include <deal.II/matrix_free/cuda_matrix_free.templates.h>
25
27
28
29namespace CUDAWrappers
30{
31 namespace internal
32 {
39 // TODO: for now only the general variant is implemented
41 {
45 };
46
50 template <int dim,
51 int n_q_points_1d,
52 typename Number,
53 int direction,
54 bool dof_to_quad,
55 bool add,
56 bool in_place,
57 typename ViewTypeIn,
58 typename ViewTypeOut>
60 apply(const Kokkos::TeamPolicy<
61 MemorySpace::Default::kokkos_space::execution_space>::member_type
62 &team_member,
63 const Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
64 shape_data,
65 const ViewTypeIn in,
66 ViewTypeOut out)
67 {
68 constexpr unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim);
69
70 Number t[n_q_points];
71 Kokkos::parallel_for(
72 Kokkos::TeamThreadRange(team_member, n_q_points),
73 [&](const int &q_point) {
74 const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
75 const unsigned int j =
76 (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
77 const unsigned int q =
78 (dim == 1) ? q_point :
79 (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
80 q_point / (n_q_points_1d * n_q_points_1d);
81
82 // This loop simply multiplies the shape function at the quadrature
83 // point by the value finite element coefficient.
84 t[q_point] = 0;
85 for (int k = 0; k < n_q_points_1d; ++k)
86 {
87 const unsigned int shape_idx =
88 dof_to_quad ? (q + k * n_q_points_1d) : (k + q * n_q_points_1d);
89 const unsigned int source_idx =
90 (direction == 0) ?
91 (k + n_q_points_1d * (i + n_q_points_1d * j)) :
92 (direction == 1) ?
93 (i + n_q_points_1d * (k + n_q_points_1d * j)) :
94 (i + n_q_points_1d * (j + n_q_points_1d * k));
95 t[q_point] += shape_data[shape_idx] *
96 (in_place ? out(source_idx) : in(source_idx));
97 }
98 });
99
100 if (in_place)
101 team_member.team_barrier();
102
103 Kokkos::parallel_for(
104 Kokkos::TeamThreadRange(team_member, n_q_points),
105 [&](const int &q_point) {
106 const unsigned int i = (dim == 1) ? 0 : q_point % n_q_points_1d;
107 const unsigned int j =
108 (dim == 3) ? (q_point / n_q_points_1d) % n_q_points_1d : 0;
109 const unsigned int q =
110 (dim == 1) ? q_point :
111 (dim == 2) ? (q_point / n_q_points_1d) % n_q_points_1d :
112 q_point / (n_q_points_1d * n_q_points_1d);
113
114 const unsigned int destination_idx =
115 (direction == 0) ? (q + n_q_points_1d * (i + n_q_points_1d * j)) :
116 (direction == 1) ? (i + n_q_points_1d * (q + n_q_points_1d * j)) :
117 (i + n_q_points_1d * (j + n_q_points_1d * q));
118
119 if (add)
120 Kokkos::atomic_add(&out(destination_idx), t[q_point]);
121 else
122 out(destination_idx) = t[q_point];
123 });
124 }
125
126
132 template <EvaluatorVariant variant,
133 int dim,
134 int fe_degree,
135 int n_q_points_1d,
136 typename Number>
138 {};
139
140
141
148 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
150 dim,
151 fe_degree,
152 n_q_points_1d,
153 Number>
154 {
155 using TeamHandle = Kokkos::TeamPolicy<
156 MemorySpace::Default::kokkos_space::execution_space>::member_type;
157
160 const TeamHandle & team_member,
161 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
162 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
163 shape_gradients,
164 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
165 co_shape_gradients);
166
171 template <int direction,
172 bool dof_to_quad,
173 bool add,
174 bool in_place,
175 typename ViewTypeIn,
176 typename ViewTypeOut>
178 values(const ViewTypeIn in, ViewTypeOut out) const;
179
184 template <int direction,
185 bool dof_to_quad,
186 bool add,
187 bool in_place,
188 typename ViewTypeIn,
189 typename ViewTypeOut>
191 gradients(const ViewTypeIn in, ViewTypeOut out) const;
192
197 template <int direction,
198 bool dof_to_quad,
199 bool add,
200 bool in_place,
201 typename ViewTypeIn,
202 typename ViewTypeOut>
204 co_gradients(const ViewTypeIn in, ViewTypeOut out) const;
205
209 template <typename ViewType>
211 value_at_quad_pts(ViewType u);
212
216 template <typename ViewType>
218 integrate_value(ViewType u);
219
224 template <typename ViewTypeIn, typename ViewTypeOut>
226 gradient_at_quad_pts(const ViewTypeIn u, ViewTypeOut grad_u);
227
232 template <typename ViewType1, typename ViewType2>
234 value_and_gradient_at_quad_pts(ViewType1 u, ViewType2 grad_u);
235
240 template <bool add, typename ViewType1, typename ViewType2>
242 integrate_gradient(ViewType1 u, ViewType2 grad_u);
243
248 template <typename ViewType1, typename ViewType2>
250 integrate_value_and_gradient(ViewType1 u, ViewType2 grad_u);
251
256
260 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values;
261
265 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
267
271 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
273 };
274
275
276
277 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
280 dim,
281 fe_degree,
282 n_q_points_1d,
283 Number>::
284 EvaluatorTensorProduct(
285 const TeamHandle & team_member,
286 Kokkos::View<Number *, MemorySpace::Default::kokkos_space> shape_values,
287 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
288 shape_gradients,
289 Kokkos::View<Number *, MemorySpace::Default::kokkos_space>
290 co_shape_gradients)
291 : team_member(team_member)
292 , shape_values(shape_values)
293 , shape_gradients(shape_gradients)
294 , co_shape_gradients(co_shape_gradients)
295 {}
296
297
298
299 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
300 template <int direction,
301 bool dof_to_quad,
302 bool add,
303 bool in_place,
304 typename ViewTypeIn,
305 typename ViewTypeOut>
308 dim,
309 fe_degree,
310 n_q_points_1d,
311 Number>::values(const ViewTypeIn in,
312 ViewTypeOut out) const
313 {
314 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
315 team_member, shape_values, in, out);
316 }
317
318
319
320 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
321 template <int direction,
322 bool dof_to_quad,
323 bool add,
324 bool in_place,
325 typename ViewTypeIn,
326 typename ViewTypeOut>
329 dim,
330 fe_degree,
331 n_q_points_1d,
332 Number>::gradients(const ViewTypeIn in,
333 ViewTypeOut out) const
334 {
335 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
336 team_member, shape_gradients, in, out);
337 }
338
339
340
341 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
342 template <int direction,
343 bool dof_to_quad,
344 bool add,
345 bool in_place,
346 typename ViewTypeIn,
347 typename ViewTypeOut>
350 dim,
351 fe_degree,
352 n_q_points_1d,
353 Number>::co_gradients(const ViewTypeIn in,
354 ViewTypeOut out) const
355 {
356 apply<dim, n_q_points_1d, Number, direction, dof_to_quad, add, in_place>(
357 team_member, co_shape_gradients, in, out);
358 }
359
360
361
362 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
363 template <typename ViewType>
364 DEAL_II_HOST_DEVICE inline void
366 dim,
367 fe_degree,
368 n_q_points_1d,
369 Number>::value_at_quad_pts(ViewType u)
370 {
371 switch (dim)
372 {
373 case 1:
374 {
375 values<0, true, false, true>(u, u);
376
377 break;
378 }
379 case 2:
380 {
381 values<0, true, false, true>(u, u);
382 team_member.team_barrier();
383 values<1, true, false, true>(u, u);
384
385 break;
386 }
387 case 3:
388 {
389 values<0, true, false, true>(u, u);
390 team_member.team_barrier();
391 values<1, true, false, true>(u, u);
392 team_member.team_barrier();
393 values<2, true, false, true>(u, u);
394
395 break;
396 }
397 default:
398 {
399 // Do nothing. We should throw but we can't from a __device__
400 // function.
401 }
402 }
403 }
404
405
406
407 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
408 template <typename ViewType>
409 DEAL_II_HOST_DEVICE inline void
411 dim,
412 fe_degree,
413 n_q_points_1d,
414 Number>::integrate_value(ViewType u)
415 {
416 switch (dim)
417 {
418 case 1:
419 {
420 values<0, false, false, true>(u, u);
421
422 break;
423 }
424 case 2:
425 {
426 values<0, false, false, true>(u, u);
427 team_member.team_barrier();
428 values<1, false, false, true>(u, u);
429
430 break;
431 }
432 case 3:
433 {
434 values<0, false, false, true>(u, u);
435 team_member.team_barrier();
436 values<1, false, false, true>(u, u);
437 team_member.team_barrier();
438 values<2, false, false, true>(u, u);
439
440 break;
441 }
442 default:
443 {
444 // Do nothing. We should throw but we can't from a __device__
445 // function.
446 }
447 }
448 }
449
450
451
452 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
453 template <typename ViewTypeIn, typename ViewTypeOut>
454 DEAL_II_HOST_DEVICE inline void
456 dim,
457 fe_degree,
458 n_q_points_1d,
459 Number>::gradient_at_quad_pts(const ViewTypeIn u,
460 ViewTypeOut grad_u)
461 {
462 switch (dim)
463 {
464 case 1:
465 {
466 gradients<0, true, false, false>(
467 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
468
469 break;
470 }
471 case 2:
472 {
473 gradients<0, true, false, false>(
474 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
475 values<0, true, false, false>(
476 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
477
478 team_member.team_barrier();
479
480 values<1, true, false, true>(
481 Kokkos::subview(grad_u, Kokkos::ALL, 0),
482 Kokkos::subview(grad_u, Kokkos::ALL, 0));
483 gradients<1, true, false, true>(
484 Kokkos::subview(grad_u, Kokkos::ALL, 1),
485 Kokkos::subview(grad_u, Kokkos::ALL, 1));
486
487 break;
488 }
489 case 3:
490 {
491 gradients<0, true, false, false>(
492 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
493 values<0, true, false, false>(
494 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
495 values<0, true, false, false>(
496 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
497
498 team_member.team_barrier();
499
500 values<1, true, false, true>(
501 Kokkos::subview(grad_u, Kokkos::ALL, 0),
502 Kokkos::subview(grad_u, Kokkos::ALL, 0));
503 gradients<1, true, false, true>(
504 Kokkos::subview(grad_u, Kokkos::ALL, 1),
505 Kokkos::subview(grad_u, Kokkos::ALL, 1));
506 values<1, true, false, true>(
507 Kokkos::subview(grad_u, Kokkos::ALL, 2),
508 Kokkos::subview(grad_u, Kokkos::ALL, 2));
509
510 team_member.team_barrier();
511
512 values<2, true, false, true>(
513 Kokkos::subview(grad_u, Kokkos::ALL, 0),
514 Kokkos::subview(grad_u, Kokkos::ALL, 0));
515 values<2, true, false, true>(
516 Kokkos::subview(grad_u, Kokkos::ALL, 1),
517 Kokkos::subview(grad_u, Kokkos::ALL, 1));
518 gradients<2, true, false, true>(
519 Kokkos::subview(grad_u, Kokkos::ALL, 2),
520 Kokkos::subview(grad_u, Kokkos::ALL, 2));
521
522 break;
523 }
524 default:
525 {
526 // Do nothing. We should throw but we can't from a __device__
527 // function.
528 }
529 }
530 }
531
532
533
534 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
535 template <typename ViewType1, typename ViewType2>
536 DEAL_II_HOST_DEVICE inline void
538 dim,
539 fe_degree,
540 n_q_points_1d,
541 Number>::value_and_gradient_at_quad_pts(ViewType1 u,
542 ViewType2
543 grad_u)
544 {
545 switch (dim)
546 {
547 case 1:
548 {
549 values<0, true, false, true>(u, u);
550 team_member.team_barrier();
551
552 co_gradients<0, true, false, false>(
553 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
554
555 break;
556 }
557 case 2:
558 {
559 values<0, true, false, true>(u, u);
560 team_member.team_barrier();
561 values<1, true, false, true>(u, u);
562 team_member.team_barrier();
563
564 co_gradients<0, true, false, false>(
565 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
566 co_gradients<1, true, false, false>(
567 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
568
569 break;
570 }
571 case 3:
572 {
573 values<0, true, false, true>(u, u);
574 team_member.team_barrier();
575 values<1, true, false, true>(u, u);
576 team_member.team_barrier();
577 values<2, true, false, true>(u, u);
578 team_member.team_barrier();
579
580 co_gradients<0, true, false, false>(
581 u, Kokkos::subview(grad_u, Kokkos::ALL, 0));
582 co_gradients<1, true, false, false>(
583 u, Kokkos::subview(grad_u, Kokkos::ALL, 1));
584 co_gradients<2, true, false, false>(
585 u, Kokkos::subview(grad_u, Kokkos::ALL, 2));
586
587 break;
588 }
589 default:
590 {
591 // Do nothing. We should throw but we can't from a __device__
592 // function.
593 }
594 }
595 }
596
597
598
599 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
600 template <bool add, typename ViewType1, typename ViewType2>
601 DEAL_II_HOST_DEVICE inline void
603 dim,
604 fe_degree,
605 n_q_points_1d,
606 Number>::integrate_gradient(ViewType1 u,
607 ViewType2 grad_u)
608 {
609 switch (dim)
610 {
611 case 1:
612 {
613 gradients<0, false, add, false>(
614 Kokkos::subview(grad_u, Kokkos::ALL, dim), u);
615
616 break;
617 }
618 case 2:
619 {
620 gradients<0, false, false, true>(
621 Kokkos::subview(grad_u, Kokkos::ALL, 0),
622 Kokkos::subview(grad_u, Kokkos::ALL, 0));
623 values<0, false, false, true>(
624 Kokkos::subview(grad_u, Kokkos::ALL, 1),
625 Kokkos::subview(grad_u, Kokkos::ALL, 1));
626
627 team_member.team_barrier();
628
629 values<1, false, add, false>(
630 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
631 team_member.team_barrier();
632 gradients<1, false, true, false>(
633 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
634
635 break;
636 }
637 case 3:
638 {
639 gradients<0, false, false, true>(
640 Kokkos::subview(grad_u, Kokkos::ALL, 0),
641 Kokkos::subview(grad_u, Kokkos::ALL, 0));
642 values<0, false, false, true>(
643 Kokkos::subview(grad_u, Kokkos::ALL, 1),
644 Kokkos::subview(grad_u, Kokkos::ALL, 1));
645 values<0, false, false, true>(
646 Kokkos::subview(grad_u, Kokkos::ALL, 2),
647 Kokkos::subview(grad_u, Kokkos::ALL, 2));
648
649 team_member.team_barrier();
650
651 values<1, false, false, true>(
652 Kokkos::subview(grad_u, Kokkos::ALL, 0),
653 Kokkos::subview(grad_u, Kokkos::ALL, 0));
654 gradients<1, false, false, true>(
655 Kokkos::subview(grad_u, Kokkos::ALL, 1),
656 Kokkos::subview(grad_u, Kokkos::ALL, 1));
657 values<1, false, false, true>(
658 Kokkos::subview(grad_u, Kokkos::ALL, 2),
659 Kokkos::subview(grad_u, Kokkos::ALL, 2));
660
661 team_member.team_barrier();
662
663 values<2, false, add, false>(
664 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
665 team_member.team_barrier();
666 values<2, false, true, false>(
667 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
668 team_member.team_barrier();
669 gradients<2, false, true, false>(
670 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
671
672 break;
673 }
674 default:
675 {
676 // Do nothing. We should throw but we can't from a __device__
677 // function.
678 }
679 }
680 }
681
682
683
684 template <int dim, int fe_degree, int n_q_points_1d, typename Number>
685 template <typename ViewType1, typename ViewType2>
686 DEAL_II_HOST_DEVICE inline void
688 dim,
689 fe_degree,
690 n_q_points_1d,
691 Number>::integrate_value_and_gradient(ViewType1 u,
692 ViewType2
693 grad_u)
694 {
695 switch (dim)
696 {
697 case 1:
698 {
699 co_gradients<0, false, true, false>(
700 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
701 team_member.team_barrier();
702
703 values<0, false, false, true>(u, u);
704
705 break;
706 }
707 case 2:
708 {
709 co_gradients<1, false, true, false>(
710 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
711 team_member.team_barrier();
712 co_gradients<0, false, true, false>(
713 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
714 team_member.team_barrier();
715
716 values<1, false, false, true>(u, u);
717 team_member.team_barrier();
718 values<0, false, false, true>(u, u);
719 team_member.team_barrier();
720
721 break;
722 }
723 case 3:
724 {
725 co_gradients<2, false, true, false>(
726 Kokkos::subview(grad_u, Kokkos::ALL, 2), u);
727 team_member.team_barrier();
728 co_gradients<1, false, true, false>(
729 Kokkos::subview(grad_u, Kokkos::ALL, 1), u);
730 team_member.team_barrier();
731 co_gradients<0, false, true, false>(
732 Kokkos::subview(grad_u, Kokkos::ALL, 0), u);
733 team_member.team_barrier();
734
735 values<2, false, false, true>(u, u);
736 team_member.team_barrier();
737 values<1, false, false, true>(u, u);
738 team_member.team_barrier();
739 values<0, false, false, true>(u, u);
740 team_member.team_barrier();
741
742 break;
743 }
744 default:
745 {
746 // Do nothing. We should throw but we can't from a __device__
747 // function.
748 }
749 }
750 }
751 } // namespace internal
752} // namespace CUDAWrappers
753
755
756#endif
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
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
#define DEAL_II_HOST_DEVICE
Definition numbers.h:35
Kokkos::TeamPolicy< MemorySpace::Default::kokkos_space::execution_space >::member_type TeamHandle