Kokkos Core Kernels Package Version of the Day
Loading...
Searching...
No Matches
Kokkos_Parallel_Reduce.hpp
1//@HEADER
2// ************************************************************************
3//
4// Kokkos v. 4.0
5// Copyright (2022) National Technology & Engineering
6// Solutions of Sandia, LLC (NTESS).
7//
8// Under the terms of Contract DE-NA0003525 with NTESS,
9// the U.S. Government retains certain rights in this software.
10//
11// Part of Kokkos, under the Apache License v2.0 with LLVM Exceptions.
12// See https://kokkos.org/LICENSE for license information.
13// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
14//
15//@HEADER
16
17#ifndef KOKKOS_IMPL_PUBLIC_INCLUDE
18#include <Kokkos_Macros.hpp>
19static_assert(false,
20 "Including non-public Kokkos header files is not allowed.");
21#endif
22#ifndef KOKKOS_PARALLEL_REDUCE_HPP
23#define KOKKOS_PARALLEL_REDUCE_HPP
24
25#include <Kokkos_ReductionIdentity.hpp>
26#include <Kokkos_View.hpp>
27#include <impl/Kokkos_FunctorAnalysis.hpp>
28#include <impl/Kokkos_Tools_Generic.hpp>
29#include <type_traits>
30#include <iostream>
31
32namespace Kokkos {
33
34template <class Scalar, class Space>
35struct Sum {
36 public:
37 // Required
38 using reducer = Sum<Scalar, Space>;
39 using value_type = std::remove_cv_t<Scalar>;
40
41 using result_view_type = Kokkos::View<value_type, Space>;
42
43 private:
44 result_view_type value;
45 bool references_scalar_v;
46
47 public:
48 KOKKOS_INLINE_FUNCTION
49 Sum(value_type& value_) : value(&value_), references_scalar_v(true) {}
50
51 KOKKOS_INLINE_FUNCTION
52 Sum(const result_view_type& value_)
53 : value(value_), references_scalar_v(false) {}
54
55 // Required
56 KOKKOS_INLINE_FUNCTION
57 void join(value_type& dest, const value_type& src) const { dest += src; }
58
59 KOKKOS_INLINE_FUNCTION
60 void init(value_type& val) const {
61 val = reduction_identity<value_type>::sum();
62 }
63
64 KOKKOS_INLINE_FUNCTION
65 value_type& reference() const { return *value.data(); }
66
67 KOKKOS_INLINE_FUNCTION
68 result_view_type view() const { return value; }
69
70 KOKKOS_INLINE_FUNCTION
71 bool references_scalar() const { return references_scalar_v; }
72};
73
74template <typename Scalar, typename... Properties>
75Sum(View<Scalar, Properties...> const&)
76 ->Sum<Scalar, typename View<Scalar, Properties...>::memory_space>;
77
78template <class Scalar, class Space>
79struct Prod {
80 public:
81 // Required
82 using reducer = Prod<Scalar, Space>;
83 using value_type = std::remove_cv_t<Scalar>;
84
85 using result_view_type = Kokkos::View<value_type, Space>;
86
87 private:
88 result_view_type value;
89 bool references_scalar_v;
90
91 public:
92 KOKKOS_INLINE_FUNCTION
93 Prod(value_type& value_) : value(&value_), references_scalar_v(true) {}
94
95 KOKKOS_INLINE_FUNCTION
96 Prod(const result_view_type& value_)
97 : value(value_), references_scalar_v(false) {}
98
99 // Required
100 KOKKOS_INLINE_FUNCTION
101 void join(value_type& dest, const value_type& src) const { dest *= src; }
102
103 KOKKOS_INLINE_FUNCTION
104 void init(value_type& val) const {
105 val = reduction_identity<value_type>::prod();
106 }
107
108 KOKKOS_INLINE_FUNCTION
109 value_type& reference() const { return *value.data(); }
110
111 KOKKOS_INLINE_FUNCTION
112 result_view_type view() const { return value; }
113
114 KOKKOS_INLINE_FUNCTION
115 bool references_scalar() const { return references_scalar_v; }
116};
117
118template <typename Scalar, typename... Properties>
119Prod(View<Scalar, Properties...> const&)
120 ->Prod<Scalar, typename View<Scalar, Properties...>::memory_space>;
121
122template <class Scalar, class Space>
123struct Min {
124 public:
125 // Required
126 using reducer = Min<Scalar, Space>;
127 using value_type = std::remove_cv_t<Scalar>;
128
129 using result_view_type = Kokkos::View<value_type, Space>;
130
131 private:
132 result_view_type value;
133 bool references_scalar_v;
134
135 public:
136 KOKKOS_INLINE_FUNCTION
137 Min(value_type& value_) : value(&value_), references_scalar_v(true) {}
138
139 KOKKOS_INLINE_FUNCTION
140 Min(const result_view_type& value_)
141 : value(value_), references_scalar_v(false) {}
142
143 // Required
144 KOKKOS_INLINE_FUNCTION
145 void join(value_type& dest, const value_type& src) const {
146 if (src < dest) dest = src;
147 }
148
149 KOKKOS_INLINE_FUNCTION
150 void init(value_type& val) const {
151 val = reduction_identity<value_type>::min();
152 }
153
154 KOKKOS_INLINE_FUNCTION
155 value_type& reference() const { return *value.data(); }
156
157 KOKKOS_INLINE_FUNCTION
158 result_view_type view() const { return value; }
159
160 KOKKOS_INLINE_FUNCTION
161 bool references_scalar() const { return references_scalar_v; }
162};
163
164template <typename Scalar, typename... Properties>
165Min(View<Scalar, Properties...> const&)
166 ->Min<Scalar, typename View<Scalar, Properties...>::memory_space>;
167
168template <class Scalar, class Space>
169struct Max {
170 public:
171 // Required
172 using reducer = Max<Scalar, Space>;
173 using value_type = std::remove_cv_t<Scalar>;
174
175 using result_view_type = Kokkos::View<value_type, Space>;
176
177 private:
178 result_view_type value;
179 bool references_scalar_v;
180
181 public:
182 KOKKOS_INLINE_FUNCTION
183 Max(value_type& value_) : value(&value_), references_scalar_v(true) {}
184
185 KOKKOS_INLINE_FUNCTION
186 Max(const result_view_type& value_)
187 : value(value_), references_scalar_v(false) {}
188
189 // Required
190 KOKKOS_INLINE_FUNCTION
191 void join(value_type& dest, const value_type& src) const {
192 if (src > dest) dest = src;
193 }
194
195 // Required
196 KOKKOS_INLINE_FUNCTION
197 void init(value_type& val) const {
198 val = reduction_identity<value_type>::max();
199 }
200
201 KOKKOS_INLINE_FUNCTION
202 value_type& reference() const { return *value.data(); }
203
204 KOKKOS_INLINE_FUNCTION
205 result_view_type view() const { return value; }
206
207 KOKKOS_INLINE_FUNCTION
208 bool references_scalar() const { return references_scalar_v; }
209};
210
211template <typename Scalar, typename... Properties>
212Max(View<Scalar, Properties...> const&)
213 ->Max<Scalar, typename View<Scalar, Properties...>::memory_space>;
214
215template <class Scalar, class Space>
216struct LAnd {
217 public:
218 // Required
219 using reducer = LAnd<Scalar, Space>;
220 using value_type = std::remove_cv_t<Scalar>;
221
222 using result_view_type = Kokkos::View<value_type, Space>;
223
224 private:
225 result_view_type value;
226 bool references_scalar_v;
227
228 public:
229 KOKKOS_INLINE_FUNCTION
230 LAnd(value_type& value_) : value(&value_), references_scalar_v(true) {}
231
232 KOKKOS_INLINE_FUNCTION
233 LAnd(const result_view_type& value_)
234 : value(value_), references_scalar_v(false) {}
235
236 KOKKOS_INLINE_FUNCTION
237 void join(value_type& dest, const value_type& src) const {
238 dest = dest && src;
239 }
240
241 KOKKOS_INLINE_FUNCTION
242 void init(value_type& val) const {
243 val = reduction_identity<value_type>::land();
244 }
245
246 KOKKOS_INLINE_FUNCTION
247 value_type& reference() const { return *value.data(); }
248
249 KOKKOS_INLINE_FUNCTION
250 result_view_type view() const { return value; }
251
252 KOKKOS_INLINE_FUNCTION
253 bool references_scalar() const { return references_scalar_v; }
254};
255
256template <typename Scalar, typename... Properties>
257LAnd(View<Scalar, Properties...> const&)
258 ->LAnd<Scalar, typename View<Scalar, Properties...>::memory_space>;
259
260template <class Scalar, class Space>
261struct LOr {
262 public:
263 // Required
264 using reducer = LOr<Scalar, Space>;
265 using value_type = std::remove_cv_t<Scalar>;
266
267 using result_view_type = Kokkos::View<value_type, Space>;
268
269 private:
270 result_view_type value;
271 bool references_scalar_v;
272
273 public:
274 KOKKOS_INLINE_FUNCTION
275 LOr(value_type& value_) : value(&value_), references_scalar_v(true) {}
276
277 KOKKOS_INLINE_FUNCTION
278 LOr(const result_view_type& value_)
279 : value(value_), references_scalar_v(false) {}
280
281 // Required
282 KOKKOS_INLINE_FUNCTION
283 void join(value_type& dest, const value_type& src) const {
284 dest = dest || src;
285 }
286
287 KOKKOS_INLINE_FUNCTION
288 void init(value_type& val) const {
289 val = reduction_identity<value_type>::lor();
290 }
291
292 KOKKOS_INLINE_FUNCTION
293 value_type& reference() const { return *value.data(); }
294
295 KOKKOS_INLINE_FUNCTION
296 result_view_type view() const { return value; }
297
298 KOKKOS_INLINE_FUNCTION
299 bool references_scalar() const { return references_scalar_v; }
300};
301
302template <typename Scalar, typename... Properties>
303LOr(View<Scalar, Properties...> const&)
304 ->LOr<Scalar, typename View<Scalar, Properties...>::memory_space>;
305
306template <class Scalar, class Space>
307struct BAnd {
308 public:
309 // Required
310 using reducer = BAnd<Scalar, Space>;
311 using value_type = std::remove_cv_t<Scalar>;
312
313 using result_view_type = Kokkos::View<value_type, Space>;
314
315 private:
316 result_view_type value;
317 bool references_scalar_v;
318
319 public:
320 KOKKOS_INLINE_FUNCTION
321 BAnd(value_type& value_) : value(&value_), references_scalar_v(true) {}
322
323 KOKKOS_INLINE_FUNCTION
324 BAnd(const result_view_type& value_)
325 : value(value_), references_scalar_v(false) {}
326
327 // Required
328 KOKKOS_INLINE_FUNCTION
329 void join(value_type& dest, const value_type& src) const {
330 dest = dest & src;
331 }
332
333 KOKKOS_INLINE_FUNCTION
334 void init(value_type& val) const {
335 val = reduction_identity<value_type>::band();
336 }
337
338 KOKKOS_INLINE_FUNCTION
339 value_type& reference() const { return *value.data(); }
340
341 KOKKOS_INLINE_FUNCTION
342 result_view_type view() const { return value; }
343
344 KOKKOS_INLINE_FUNCTION
345 bool references_scalar() const { return references_scalar_v; }
346};
347
348template <typename Scalar, typename... Properties>
349BAnd(View<Scalar, Properties...> const&)
350 ->BAnd<Scalar, typename View<Scalar, Properties...>::memory_space>;
351
352template <class Scalar, class Space>
353struct BOr {
354 public:
355 // Required
356 using reducer = BOr<Scalar, Space>;
357 using value_type = std::remove_cv_t<Scalar>;
358
359 using result_view_type = Kokkos::View<value_type, Space>;
360
361 private:
362 result_view_type value;
363 bool references_scalar_v;
364
365 public:
366 KOKKOS_INLINE_FUNCTION
367 BOr(value_type& value_) : value(&value_), references_scalar_v(true) {}
368
369 KOKKOS_INLINE_FUNCTION
370 BOr(const result_view_type& value_)
371 : value(value_), references_scalar_v(false) {}
372
373 // Required
374 KOKKOS_INLINE_FUNCTION
375 void join(value_type& dest, const value_type& src) const {
376 dest = dest | src;
377 }
378
379 KOKKOS_INLINE_FUNCTION
380 void init(value_type& val) const {
381 val = reduction_identity<value_type>::bor();
382 }
383
384 KOKKOS_INLINE_FUNCTION
385 value_type& reference() const { return *value.data(); }
386
387 KOKKOS_INLINE_FUNCTION
388 result_view_type view() const { return value; }
389
390 KOKKOS_INLINE_FUNCTION
391 bool references_scalar() const { return references_scalar_v; }
392};
393
394template <typename Scalar, typename... Properties>
395BOr(View<Scalar, Properties...> const&)
396 ->BOr<Scalar, typename View<Scalar, Properties...>::memory_space>;
397
398template <class Scalar, class Index>
399struct ValLocScalar {
400 Scalar val;
401 Index loc;
402
403 KOKKOS_INLINE_FUNCTION
404 void operator=(const ValLocScalar& rhs) {
405 val = rhs.val;
406 loc = rhs.loc;
407 }
408};
409
410template <class Scalar, class Index, class Space>
411struct MinLoc {
412 private:
413 using scalar_type = std::remove_cv_t<Scalar>;
414 using index_type = std::remove_cv_t<Index>;
415
416 public:
417 // Required
418 using reducer = MinLoc<Scalar, Index, Space>;
419 using value_type = ValLocScalar<scalar_type, index_type>;
420
421 using result_view_type = Kokkos::View<value_type, Space>;
422
423 private:
424 result_view_type value;
425 bool references_scalar_v;
426
427 public:
428 KOKKOS_INLINE_FUNCTION
429 MinLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
430
431 KOKKOS_INLINE_FUNCTION
432 MinLoc(const result_view_type& value_)
433 : value(value_), references_scalar_v(false) {}
434
435 // Required
436 KOKKOS_INLINE_FUNCTION
437 void join(value_type& dest, const value_type& src) const {
438 if (src.val < dest.val) dest = src;
439 }
440
441 KOKKOS_INLINE_FUNCTION
442 void init(value_type& val) const {
443 val.val = reduction_identity<scalar_type>::min();
444 val.loc = reduction_identity<index_type>::min();
445 }
446
447 KOKKOS_INLINE_FUNCTION
448 value_type& reference() const { return *value.data(); }
449
450 KOKKOS_INLINE_FUNCTION
451 result_view_type view() const { return value; }
452
453 KOKKOS_INLINE_FUNCTION
454 bool references_scalar() const { return references_scalar_v; }
455};
456
457template <typename Scalar, typename Index, typename... Properties>
458MinLoc(View<ValLocScalar<Scalar, Index>, Properties...> const&)
459 ->MinLoc<Scalar, Index,
460 typename View<ValLocScalar<Scalar, Index>,
461 Properties...>::memory_space>;
462
463template <class Scalar, class Index, class Space>
464struct MaxLoc {
465 private:
466 using scalar_type = std::remove_cv_t<Scalar>;
467 using index_type = std::remove_cv_t<Index>;
468
469 public:
470 // Required
471 using reducer = MaxLoc<Scalar, Index, Space>;
472 using value_type = ValLocScalar<scalar_type, index_type>;
473
474 using result_view_type = Kokkos::View<value_type, Space>;
475
476 private:
477 result_view_type value;
478 bool references_scalar_v;
479
480 public:
481 KOKKOS_INLINE_FUNCTION
482 MaxLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
483
484 KOKKOS_INLINE_FUNCTION
485 MaxLoc(const result_view_type& value_)
486 : value(value_), references_scalar_v(false) {}
487
488 // Required
489 KOKKOS_INLINE_FUNCTION
490 void join(value_type& dest, const value_type& src) const {
491 if (src.val > dest.val) dest = src;
492 }
493
494 KOKKOS_INLINE_FUNCTION
495 void init(value_type& val) const {
496 val.val = reduction_identity<scalar_type>::max();
497 val.loc = reduction_identity<index_type>::min();
498 }
499
500 KOKKOS_INLINE_FUNCTION
501 value_type& reference() const { return *value.data(); }
502
503 KOKKOS_INLINE_FUNCTION
504 result_view_type view() const { return value; }
505
506 KOKKOS_INLINE_FUNCTION
507 bool references_scalar() const { return references_scalar_v; }
508};
509
510template <typename Scalar, typename Index, typename... Properties>
511MaxLoc(View<ValLocScalar<Scalar, Index>, Properties...> const&)
512 ->MaxLoc<Scalar, Index,
513 typename View<ValLocScalar<Scalar, Index>,
514 Properties...>::memory_space>;
515
516template <class Scalar>
517struct MinMaxScalar {
518 Scalar min_val, max_val;
519
520 KOKKOS_INLINE_FUNCTION
521 void operator=(const MinMaxScalar& rhs) {
522 min_val = rhs.min_val;
523 max_val = rhs.max_val;
524 }
525};
526
527template <class Scalar, class Space>
528struct MinMax {
529 private:
530 using scalar_type = std::remove_cv_t<Scalar>;
531
532 public:
533 // Required
534 using reducer = MinMax<Scalar, Space>;
535 using value_type = MinMaxScalar<scalar_type>;
536
537 using result_view_type = Kokkos::View<value_type, Space>;
538
539 private:
540 result_view_type value;
541 bool references_scalar_v;
542
543 public:
544 KOKKOS_INLINE_FUNCTION
545 MinMax(value_type& value_) : value(&value_), references_scalar_v(true) {}
546
547 KOKKOS_INLINE_FUNCTION
548 MinMax(const result_view_type& value_)
549 : value(value_), references_scalar_v(false) {}
550
551 // Required
552 KOKKOS_INLINE_FUNCTION
553 void join(value_type& dest, const value_type& src) const {
554 if (src.min_val < dest.min_val) {
555 dest.min_val = src.min_val;
556 }
557 if (src.max_val > dest.max_val) {
558 dest.max_val = src.max_val;
559 }
560 }
561
562 KOKKOS_INLINE_FUNCTION
563 void init(value_type& val) const {
564 val.max_val = reduction_identity<scalar_type>::max();
565 val.min_val = reduction_identity<scalar_type>::min();
566 }
567
568 KOKKOS_INLINE_FUNCTION
569 value_type& reference() const { return *value.data(); }
570
571 KOKKOS_INLINE_FUNCTION
572 result_view_type view() const { return value; }
573
574 KOKKOS_INLINE_FUNCTION
575 bool references_scalar() const { return references_scalar_v; }
576};
577
578template <typename Scalar, typename... Properties>
579MinMax(View<MinMaxScalar<Scalar>, Properties...> const&)
580 ->MinMax<Scalar,
581 typename View<MinMaxScalar<Scalar>, Properties...>::memory_space>;
582
583template <class Scalar, class Index>
584struct MinMaxLocScalar {
585 Scalar min_val, max_val;
586 Index min_loc, max_loc;
587
588 KOKKOS_INLINE_FUNCTION
589 void operator=(const MinMaxLocScalar& rhs) {
590 min_val = rhs.min_val;
591 min_loc = rhs.min_loc;
592 max_val = rhs.max_val;
593 max_loc = rhs.max_loc;
594 }
595};
596
597template <class Scalar, class Index, class Space>
598struct MinMaxLoc {
599 private:
600 using scalar_type = std::remove_cv_t<Scalar>;
601 using index_type = std::remove_cv_t<Index>;
602
603 public:
604 // Required
605 using reducer = MinMaxLoc<Scalar, Index, Space>;
606 using value_type = MinMaxLocScalar<scalar_type, index_type>;
607
608 using result_view_type = Kokkos::View<value_type, Space>;
609
610 private:
611 result_view_type value;
612 bool references_scalar_v;
613
614 public:
615 KOKKOS_INLINE_FUNCTION
616 MinMaxLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
617
618 KOKKOS_INLINE_FUNCTION
619 MinMaxLoc(const result_view_type& value_)
620 : value(value_), references_scalar_v(false) {}
621
622 // Required
623 KOKKOS_INLINE_FUNCTION
624 void join(value_type& dest, const value_type& src) const {
625 if (src.min_val < dest.min_val) {
626 dest.min_val = src.min_val;
627 dest.min_loc = src.min_loc;
628 }
629 if (src.max_val > dest.max_val) {
630 dest.max_val = src.max_val;
631 dest.max_loc = src.max_loc;
632 }
633 }
634
635 KOKKOS_INLINE_FUNCTION
636 void init(value_type& val) const {
637 val.max_val = reduction_identity<scalar_type>::max();
638 val.min_val = reduction_identity<scalar_type>::min();
639 val.max_loc = reduction_identity<index_type>::min();
640 val.min_loc = reduction_identity<index_type>::min();
641 }
642
643 KOKKOS_INLINE_FUNCTION
644 value_type& reference() const { return *value.data(); }
645
646 KOKKOS_INLINE_FUNCTION
647 result_view_type view() const { return value; }
648
649 KOKKOS_INLINE_FUNCTION
650 bool references_scalar() const { return references_scalar_v; }
651};
652
653template <typename Scalar, typename Index, typename... Properties>
654MinMaxLoc(View<MinMaxLocScalar<Scalar, Index>, Properties...> const&)
655 ->MinMaxLoc<Scalar, Index,
656 typename View<MinMaxLocScalar<Scalar, Index>,
657 Properties...>::memory_space>;
658
659// --------------------------------------------------
660// reducers added to support std algorithms
661// --------------------------------------------------
662
663//
664// MaxFirstLoc
665//
666template <class Scalar, class Index, class Space>
667struct MaxFirstLoc {
668 private:
669 using scalar_type = std::remove_cv_t<Scalar>;
670 using index_type = std::remove_cv_t<Index>;
671
672 public:
673 // Required
674 using reducer = MaxFirstLoc<Scalar, Index, Space>;
675 using value_type = ::Kokkos::ValLocScalar<scalar_type, index_type>;
676
677 using result_view_type = ::Kokkos::View<value_type, Space>;
678
679 private:
680 result_view_type value;
681 bool references_scalar_v;
682
683 public:
684 KOKKOS_INLINE_FUNCTION
685 MaxFirstLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
686
687 KOKKOS_INLINE_FUNCTION
688 MaxFirstLoc(const result_view_type& value_)
689 : value(value_), references_scalar_v(false) {}
690
691 // Required
692 KOKKOS_INLINE_FUNCTION
693 void join(value_type& dest, const value_type& src) const {
694 if (dest.val < src.val) {
695 dest = src;
696 } else if (!(src.val < dest.val)) {
697 dest.loc = (src.loc < dest.loc) ? src.loc : dest.loc;
698 }
699 }
700
701 KOKKOS_INLINE_FUNCTION
702 void init(value_type& val) const {
703 val.val = reduction_identity<scalar_type>::max();
704 val.loc = reduction_identity<index_type>::min();
705 }
706
707 KOKKOS_INLINE_FUNCTION
708 value_type& reference() const { return *value.data(); }
709
710 KOKKOS_INLINE_FUNCTION
711 result_view_type view() const { return value; }
712
713 KOKKOS_INLINE_FUNCTION
714 bool references_scalar() const { return references_scalar_v; }
715};
716
717template <typename Scalar, typename Index, typename... Properties>
718MaxFirstLoc(View<ValLocScalar<Scalar, Index>, Properties...> const&)
719 ->MaxFirstLoc<Scalar, Index,
720 typename View<ValLocScalar<Scalar, Index>,
721 Properties...>::memory_space>;
722
723//
724// MaxFirstLocCustomComparator
725// recall that comp(a,b) returns true is a < b
726//
727template <class Scalar, class Index, class ComparatorType, class Space>
728struct MaxFirstLocCustomComparator {
729 private:
730 using scalar_type = std::remove_cv_t<Scalar>;
731 using index_type = std::remove_cv_t<Index>;
732
733 public:
734 // Required
735 using reducer =
736 MaxFirstLocCustomComparator<Scalar, Index, ComparatorType, Space>;
737 using value_type = ::Kokkos::ValLocScalar<scalar_type, index_type>;
738
739 using result_view_type = ::Kokkos::View<value_type, Space>;
740
741 private:
742 result_view_type value;
743 bool references_scalar_v;
744 ComparatorType m_comp;
745
746 public:
747 KOKKOS_INLINE_FUNCTION
748 MaxFirstLocCustomComparator(value_type& value_, ComparatorType comp_)
749 : value(&value_), references_scalar_v(true), m_comp(comp_) {}
750
751 KOKKOS_INLINE_FUNCTION
752 MaxFirstLocCustomComparator(const result_view_type& value_,
753 ComparatorType comp_)
754 : value(value_), references_scalar_v(false), m_comp(comp_) {}
755
756 // Required
757 KOKKOS_INLINE_FUNCTION
758 void join(value_type& dest, const value_type& src) const {
759 if (m_comp(dest.val, src.val)) {
760 dest = src;
761 } else if (!m_comp(src.val, dest.val)) {
762 dest.loc = (src.loc < dest.loc) ? src.loc : dest.loc;
763 }
764 }
765
766 KOKKOS_INLINE_FUNCTION
767 void init(value_type& val) const {
768 val.val = reduction_identity<scalar_type>::max();
769 val.loc = reduction_identity<index_type>::min();
770 }
771
772 KOKKOS_INLINE_FUNCTION
773 value_type& reference() const { return *value.data(); }
774
775 KOKKOS_INLINE_FUNCTION
776 result_view_type view() const { return value; }
777
778 KOKKOS_INLINE_FUNCTION
779 bool references_scalar() const { return references_scalar_v; }
780};
781
782template <typename Scalar, typename Index, typename ComparatorType,
783 typename... Properties>
784MaxFirstLocCustomComparator(
785 View<ValLocScalar<Scalar, Index>, Properties...> const&, ComparatorType)
786 ->MaxFirstLocCustomComparator<Scalar, Index, ComparatorType,
787 typename View<ValLocScalar<Scalar, Index>,
788 Properties...>::memory_space>;
789
790//
791// MinFirstLoc
792//
793template <class Scalar, class Index, class Space>
794struct MinFirstLoc {
795 private:
796 using scalar_type = std::remove_cv_t<Scalar>;
797 using index_type = std::remove_cv_t<Index>;
798
799 public:
800 // Required
801 using reducer = MinFirstLoc<Scalar, Index, Space>;
802 using value_type = ::Kokkos::ValLocScalar<scalar_type, index_type>;
803
804 using result_view_type = ::Kokkos::View<value_type, Space>;
805
806 private:
807 result_view_type value;
808 bool references_scalar_v;
809
810 public:
811 KOKKOS_INLINE_FUNCTION
812 MinFirstLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
813
814 KOKKOS_INLINE_FUNCTION
815 MinFirstLoc(const result_view_type& value_)
816 : value(value_), references_scalar_v(false) {}
817
818 // Required
819 KOKKOS_INLINE_FUNCTION
820 void join(value_type& dest, const value_type& src) const {
821 if (src.val < dest.val) {
822 dest = src;
823 } else if (!(dest.val < src.val)) {
824 dest.loc = (src.loc < dest.loc) ? src.loc : dest.loc;
825 }
826 }
827
828 KOKKOS_INLINE_FUNCTION
829 void init(value_type& val) const {
830 val.val = reduction_identity<scalar_type>::min();
831 val.loc = reduction_identity<index_type>::min();
832 }
833
834 KOKKOS_INLINE_FUNCTION
835 value_type& reference() const { return *value.data(); }
836
837 KOKKOS_INLINE_FUNCTION
838 result_view_type view() const { return value; }
839
840 KOKKOS_INLINE_FUNCTION
841 bool references_scalar() const { return references_scalar_v; }
842};
843
844template <typename Scalar, typename Index, typename... Properties>
845MinFirstLoc(View<ValLocScalar<Scalar, Index>, Properties...> const&)
846 ->MinFirstLoc<Scalar, Index,
847 typename View<ValLocScalar<Scalar, Index>,
848 Properties...>::memory_space>;
849
850//
851// MinFirstLocCustomComparator
852// recall that comp(a,b) returns true is a < b
853//
854template <class Scalar, class Index, class ComparatorType, class Space>
855struct MinFirstLocCustomComparator {
856 private:
857 using scalar_type = std::remove_cv_t<Scalar>;
858 using index_type = std::remove_cv_t<Index>;
859
860 public:
861 // Required
862 using reducer =
863 MinFirstLocCustomComparator<Scalar, Index, ComparatorType, Space>;
864 using value_type = ::Kokkos::ValLocScalar<scalar_type, index_type>;
865
866 using result_view_type = ::Kokkos::View<value_type, Space>;
867
868 private:
869 result_view_type value;
870 bool references_scalar_v;
871 ComparatorType m_comp;
872
873 public:
874 KOKKOS_INLINE_FUNCTION
875 MinFirstLocCustomComparator(value_type& value_, ComparatorType comp_)
876 : value(&value_), references_scalar_v(true), m_comp(comp_) {}
877
878 KOKKOS_INLINE_FUNCTION
879 MinFirstLocCustomComparator(const result_view_type& value_,
880 ComparatorType comp_)
881 : value(value_), references_scalar_v(false), m_comp(comp_) {}
882
883 // Required
884 KOKKOS_INLINE_FUNCTION
885 void join(value_type& dest, const value_type& src) const {
886 if (m_comp(src.val, dest.val)) {
887 dest = src;
888 } else if (!m_comp(dest.val, src.val)) {
889 dest.loc = (src.loc < dest.loc) ? src.loc : dest.loc;
890 }
891 }
892
893 KOKKOS_INLINE_FUNCTION
894 void init(value_type& val) const {
895 val.val = reduction_identity<scalar_type>::min();
896 val.loc = reduction_identity<index_type>::min();
897 }
898
899 KOKKOS_INLINE_FUNCTION
900 value_type& reference() const { return *value.data(); }
901
902 KOKKOS_INLINE_FUNCTION
903 result_view_type view() const { return value; }
904
905 KOKKOS_INLINE_FUNCTION
906 bool references_scalar() const { return references_scalar_v; }
907};
908
909template <typename Scalar, typename Index, typename ComparatorType,
910 typename... Properties>
911MinFirstLocCustomComparator(
912 View<ValLocScalar<Scalar, Index>, Properties...> const&, ComparatorType)
913 ->MinFirstLocCustomComparator<Scalar, Index, ComparatorType,
914 typename View<ValLocScalar<Scalar, Index>,
915 Properties...>::memory_space>;
916
917//
918// MinMaxFirstLastLoc
919//
920template <class Scalar, class Index, class Space>
921struct MinMaxFirstLastLoc {
922 private:
923 using scalar_type = std::remove_cv_t<Scalar>;
924 using index_type = std::remove_cv_t<Index>;
925
926 public:
927 // Required
928 using reducer = MinMaxFirstLastLoc<Scalar, Index, Space>;
929 using value_type = ::Kokkos::MinMaxLocScalar<scalar_type, index_type>;
930
931 using result_view_type = ::Kokkos::View<value_type, Space>;
932
933 private:
934 result_view_type value;
935 bool references_scalar_v;
936
937 public:
938 KOKKOS_INLINE_FUNCTION
939 MinMaxFirstLastLoc(value_type& value_)
940 : value(&value_), references_scalar_v(true) {}
941
942 KOKKOS_INLINE_FUNCTION
943 MinMaxFirstLastLoc(const result_view_type& value_)
944 : value(value_), references_scalar_v(false) {}
945
946 // Required
947 KOKKOS_INLINE_FUNCTION
948 void join(value_type& dest, const value_type& src) const {
949 if (src.min_val < dest.min_val) {
950 dest.min_val = src.min_val;
951 dest.min_loc = src.min_loc;
952 } else if (!(dest.min_val < src.min_val)) {
953 dest.min_loc = (src.min_loc < dest.min_loc) ? src.min_loc : dest.min_loc;
954 }
955
956 if (dest.max_val < src.max_val) {
957 dest.max_val = src.max_val;
958 dest.max_loc = src.max_loc;
959 } else if (!(src.max_val < dest.max_val)) {
960 dest.max_loc = (src.max_loc > dest.max_loc) ? src.max_loc : dest.max_loc;
961 }
962 }
963
964 KOKKOS_INLINE_FUNCTION
965 void init(value_type& val) const {
966 val.max_val = ::Kokkos::reduction_identity<scalar_type>::max();
967 val.min_val = ::Kokkos::reduction_identity<scalar_type>::min();
968 val.max_loc = ::Kokkos::reduction_identity<index_type>::max();
969 val.min_loc = ::Kokkos::reduction_identity<index_type>::min();
970 }
971
972 KOKKOS_INLINE_FUNCTION
973 value_type& reference() const { return *value.data(); }
974
975 KOKKOS_INLINE_FUNCTION
976 result_view_type view() const { return value; }
977
978 KOKKOS_INLINE_FUNCTION
979 bool references_scalar() const { return references_scalar_v; }
980};
981
982template <typename Scalar, typename Index, typename... Properties>
983MinMaxFirstLastLoc(View<MinMaxLocScalar<Scalar, Index>, Properties...> const&)
984 ->MinMaxFirstLastLoc<Scalar, Index,
985 typename View<MinMaxLocScalar<Scalar, Index>,
986 Properties...>::memory_space>;
987
988//
989// MinMaxFirstLastLocCustomComparator
990// recall that comp(a,b) returns true is a < b
991//
992template <class Scalar, class Index, class ComparatorType, class Space>
993struct MinMaxFirstLastLocCustomComparator {
994 private:
995 using scalar_type = std::remove_cv_t<Scalar>;
996 using index_type = std::remove_cv_t<Index>;
997
998 public:
999 // Required
1000 using reducer =
1001 MinMaxFirstLastLocCustomComparator<Scalar, Index, ComparatorType, Space>;
1002 using value_type = ::Kokkos::MinMaxLocScalar<scalar_type, index_type>;
1003
1004 using result_view_type = ::Kokkos::View<value_type, Space>;
1005
1006 private:
1007 result_view_type value;
1008 bool references_scalar_v;
1009 ComparatorType m_comp;
1010
1011 public:
1012 KOKKOS_INLINE_FUNCTION
1013 MinMaxFirstLastLocCustomComparator(value_type& value_, ComparatorType comp_)
1014 : value(&value_), references_scalar_v(true), m_comp(comp_) {}
1015
1016 KOKKOS_INLINE_FUNCTION
1017 MinMaxFirstLastLocCustomComparator(const result_view_type& value_,
1018 ComparatorType comp_)
1019 : value(value_), references_scalar_v(false), m_comp(comp_) {}
1020
1021 // Required
1022 KOKKOS_INLINE_FUNCTION
1023 void join(value_type& dest, const value_type& src) const {
1024 if (m_comp(src.min_val, dest.min_val)) {
1025 dest.min_val = src.min_val;
1026 dest.min_loc = src.min_loc;
1027 } else if (!m_comp(dest.min_val, src.min_val)) {
1028 dest.min_loc = (src.min_loc < dest.min_loc) ? src.min_loc : dest.min_loc;
1029 }
1030
1031 if (m_comp(dest.max_val, src.max_val)) {
1032 dest.max_val = src.max_val;
1033 dest.max_loc = src.max_loc;
1034 } else if (!m_comp(src.max_val, dest.max_val)) {
1035 dest.max_loc = (src.max_loc > dest.max_loc) ? src.max_loc : dest.max_loc;
1036 }
1037 }
1038
1039 KOKKOS_INLINE_FUNCTION
1040 void init(value_type& val) const {
1041 val.max_val = ::Kokkos::reduction_identity<scalar_type>::max();
1042 val.min_val = ::Kokkos::reduction_identity<scalar_type>::min();
1043 val.max_loc = ::Kokkos::reduction_identity<index_type>::max();
1044 val.min_loc = ::Kokkos::reduction_identity<index_type>::min();
1045 }
1046
1047 KOKKOS_INLINE_FUNCTION
1048 value_type& reference() const { return *value.data(); }
1049
1050 KOKKOS_INLINE_FUNCTION
1051 result_view_type view() const { return value; }
1052
1053 KOKKOS_INLINE_FUNCTION
1054 bool references_scalar() const { return references_scalar_v; }
1055};
1056
1057template <typename Scalar, typename Index, typename ComparatorType,
1058 typename... Properties>
1059MinMaxFirstLastLocCustomComparator(
1060 View<MinMaxLocScalar<Scalar, Index>, Properties...> const&, ComparatorType)
1061 ->MinMaxFirstLastLocCustomComparator<
1062 Scalar, Index, ComparatorType,
1063 typename View<MinMaxLocScalar<Scalar, Index>,
1064 Properties...>::memory_space>;
1065
1066//
1067// FirstLoc
1068//
1069template <class Index>
1070struct FirstLocScalar {
1071 Index min_loc_true;
1072
1073 KOKKOS_INLINE_FUNCTION
1074 void operator=(const FirstLocScalar& rhs) { min_loc_true = rhs.min_loc_true; }
1075};
1076
1077template <class Index, class Space>
1078struct FirstLoc {
1079 private:
1080 using index_type = std::remove_cv_t<Index>;
1081
1082 public:
1083 // Required
1084 using reducer = FirstLoc<Index, Space>;
1085 using value_type = FirstLocScalar<index_type>;
1086
1087 using result_view_type = ::Kokkos::View<value_type, Space>;
1088
1089 private:
1090 result_view_type value;
1091 bool references_scalar_v;
1092
1093 public:
1094 KOKKOS_INLINE_FUNCTION
1095 FirstLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
1096
1097 KOKKOS_INLINE_FUNCTION
1098 FirstLoc(const result_view_type& value_)
1099 : value(value_), references_scalar_v(false) {}
1100
1101 // Required
1102 KOKKOS_INLINE_FUNCTION
1103 void join(value_type& dest, const value_type& src) const {
1104 dest.min_loc_true = (src.min_loc_true < dest.min_loc_true)
1105 ? src.min_loc_true
1106 : dest.min_loc_true;
1107 }
1108
1109 KOKKOS_INLINE_FUNCTION
1110 void init(value_type& val) const {
1111 val.min_loc_true = ::Kokkos::reduction_identity<index_type>::min();
1112 }
1113
1114 KOKKOS_INLINE_FUNCTION
1115 value_type& reference() const { return *value.data(); }
1116
1117 KOKKOS_INLINE_FUNCTION
1118 result_view_type view() const { return value; }
1119
1120 KOKKOS_INLINE_FUNCTION
1121 bool references_scalar() const { return references_scalar_v; }
1122};
1123
1124template <typename Index, typename... Properties>
1125FirstLoc(View<FirstLocScalar<Index>, Properties...> const&)
1126 ->FirstLoc<Index, typename View<FirstLocScalar<Index>,
1127 Properties...>::memory_space>;
1128
1129//
1130// LastLoc
1131//
1132template <class Index>
1133struct LastLocScalar {
1134 Index max_loc_true;
1135
1136 KOKKOS_INLINE_FUNCTION
1137 void operator=(const LastLocScalar& rhs) { max_loc_true = rhs.max_loc_true; }
1138};
1139
1140template <class Index, class Space>
1141struct LastLoc {
1142 private:
1143 using index_type = std::remove_cv_t<Index>;
1144
1145 public:
1146 // Required
1147 using reducer = LastLoc<Index, Space>;
1148 using value_type = LastLocScalar<index_type>;
1149
1150 using result_view_type = ::Kokkos::View<value_type, Space>;
1151
1152 private:
1153 result_view_type value;
1154 bool references_scalar_v;
1155
1156 public:
1157 KOKKOS_INLINE_FUNCTION
1158 LastLoc(value_type& value_) : value(&value_), references_scalar_v(true) {}
1159
1160 KOKKOS_INLINE_FUNCTION
1161 LastLoc(const result_view_type& value_)
1162 : value(value_), references_scalar_v(false) {}
1163
1164 // Required
1165 KOKKOS_INLINE_FUNCTION
1166 void join(value_type& dest, const value_type& src) const {
1167 dest.max_loc_true = (src.max_loc_true > dest.max_loc_true)
1168 ? src.max_loc_true
1169 : dest.max_loc_true;
1170 }
1171
1172 KOKKOS_INLINE_FUNCTION
1173 void init(value_type& val) const {
1174 val.max_loc_true = ::Kokkos::reduction_identity<index_type>::max();
1175 }
1176
1177 KOKKOS_INLINE_FUNCTION
1178 value_type& reference() const { return *value.data(); }
1179
1180 KOKKOS_INLINE_FUNCTION
1181 result_view_type view() const { return value; }
1182
1183 KOKKOS_INLINE_FUNCTION
1184 bool references_scalar() const { return references_scalar_v; }
1185};
1186
1187template <typename Index, typename... Properties>
1188LastLoc(View<LastLocScalar<Index>, Properties...> const&)
1189 ->LastLoc<Index,
1190 typename View<LastLocScalar<Index>, Properties...>::memory_space>;
1191
1192template <class Index>
1193struct StdIsPartScalar {
1194 Index max_loc_true, min_loc_false;
1195
1196 KOKKOS_INLINE_FUNCTION
1197 void operator=(const StdIsPartScalar& rhs) {
1198 min_loc_false = rhs.min_loc_false;
1199 max_loc_true = rhs.max_loc_true;
1200 }
1201};
1202
1203//
1204// StdIsPartitioned
1205//
1206template <class Index, class Space>
1207struct StdIsPartitioned {
1208 private:
1209 using index_type = std::remove_cv_t<Index>;
1210
1211 public:
1212 // Required
1213 using reducer = StdIsPartitioned<Index, Space>;
1214 using value_type = StdIsPartScalar<index_type>;
1215
1216 using result_view_type = ::Kokkos::View<value_type, Space>;
1217
1218 private:
1219 result_view_type value;
1220 bool references_scalar_v;
1221
1222 public:
1223 KOKKOS_INLINE_FUNCTION
1224 StdIsPartitioned(value_type& value_)
1225 : value(&value_), references_scalar_v(true) {}
1226
1227 KOKKOS_INLINE_FUNCTION
1228 StdIsPartitioned(const result_view_type& value_)
1229 : value(value_), references_scalar_v(false) {}
1230
1231 // Required
1232 KOKKOS_INLINE_FUNCTION
1233 void join(value_type& dest, const value_type& src) const {
1234 dest.max_loc_true = (dest.max_loc_true < src.max_loc_true)
1235 ? src.max_loc_true
1236 : dest.max_loc_true;
1237
1238 dest.min_loc_false = (dest.min_loc_false < src.min_loc_false)
1239 ? dest.min_loc_false
1240 : src.min_loc_false;
1241 }
1242
1243 KOKKOS_INLINE_FUNCTION
1244 void init(value_type& val) const {
1245 val.max_loc_true = ::Kokkos::reduction_identity<index_type>::max();
1246 val.min_loc_false = ::Kokkos::reduction_identity<index_type>::min();
1247 }
1248
1249 KOKKOS_INLINE_FUNCTION
1250 value_type& reference() const { return *value.data(); }
1251
1252 KOKKOS_INLINE_FUNCTION
1253 result_view_type view() const { return value; }
1254
1255 KOKKOS_INLINE_FUNCTION
1256 bool references_scalar() const { return references_scalar_v; }
1257};
1258
1259template <typename Index, typename... Properties>
1260StdIsPartitioned(View<StdIsPartScalar<Index>, Properties...> const&)
1261 ->StdIsPartitioned<Index, typename View<StdIsPartScalar<Index>,
1262 Properties...>::memory_space>;
1263
1264template <class Index>
1265struct StdPartPointScalar {
1266 Index min_loc_false;
1267
1268 KOKKOS_INLINE_FUNCTION
1269 void operator=(const StdPartPointScalar& rhs) {
1270 min_loc_false = rhs.min_loc_false;
1271 }
1272};
1273
1274//
1275// StdPartitionPoint
1276//
1277template <class Index, class Space>
1278struct StdPartitionPoint {
1279 private:
1280 using index_type = std::remove_cv_t<Index>;
1281
1282 public:
1283 // Required
1284 using reducer = StdPartitionPoint<Index, Space>;
1285 using value_type = StdPartPointScalar<index_type>;
1286
1287 using result_view_type = ::Kokkos::View<value_type, Space>;
1288
1289 private:
1290 result_view_type value;
1291 bool references_scalar_v;
1292
1293 public:
1294 KOKKOS_INLINE_FUNCTION
1295 StdPartitionPoint(value_type& value_)
1296 : value(&value_), references_scalar_v(true) {}
1297
1298 KOKKOS_INLINE_FUNCTION
1299 StdPartitionPoint(const result_view_type& value_)
1300 : value(value_), references_scalar_v(false) {}
1301
1302 // Required
1303 KOKKOS_INLINE_FUNCTION
1304 void join(value_type& dest, const value_type& src) const {
1305 dest.min_loc_false = (dest.min_loc_false < src.min_loc_false)
1306 ? dest.min_loc_false
1307 : src.min_loc_false;
1308 }
1309
1310 KOKKOS_INLINE_FUNCTION
1311 void init(value_type& val) const {
1312 val.min_loc_false = ::Kokkos::reduction_identity<index_type>::min();
1313 }
1314
1315 KOKKOS_INLINE_FUNCTION
1316 value_type& reference() const { return *value.data(); }
1317
1318 KOKKOS_INLINE_FUNCTION
1319 result_view_type view() const { return value; }
1320
1321 KOKKOS_INLINE_FUNCTION
1322 bool references_scalar() const { return references_scalar_v; }
1323};
1324
1325template <typename Index, typename... Properties>
1326StdPartitionPoint(View<StdPartPointScalar<Index>, Properties...> const&)
1327 ->StdPartitionPoint<Index, typename View<StdPartPointScalar<Index>,
1328 Properties...>::memory_space>;
1329
1330} // namespace Kokkos
1331namespace Kokkos {
1332namespace Impl {
1333
1334template <class T, class ReturnType, class ValueTraits>
1335struct ParallelReduceReturnValue;
1336
1337template <class ReturnType, class FunctorType>
1338struct ParallelReduceReturnValue<
1339 std::enable_if_t<Kokkos::is_view<ReturnType>::value>, ReturnType,
1340 FunctorType> {
1341 using return_type = ReturnType;
1342 using reducer_type = InvalidType;
1343
1344 using value_type_scalar = typename return_type::value_type;
1345 using value_type_array = typename return_type::value_type* const;
1346
1347 using value_type = std::conditional_t<return_type::rank == 0,
1348 value_type_scalar, value_type_array>;
1349
1350 static return_type& return_value(ReturnType& return_val, const FunctorType&) {
1351 return return_val;
1352 }
1353};
1354
1355template <class ReturnType, class FunctorType>
1356struct ParallelReduceReturnValue<
1357 std::enable_if_t<!Kokkos::is_view<ReturnType>::value &&
1358 (!std::is_array<ReturnType>::value &&
1359 !std::is_pointer<ReturnType>::value) &&
1360 !Kokkos::is_reducer<ReturnType>::value>,
1361 ReturnType, FunctorType> {
1362 using return_type =
1364
1365 using reducer_type = InvalidType;
1366
1367 using value_type = typename return_type::value_type;
1368
1369 static return_type return_value(ReturnType& return_val, const FunctorType&) {
1370 return return_type(&return_val);
1371 }
1372};
1373
1374template <class ReturnType, class FunctorType>
1375struct ParallelReduceReturnValue<
1376 std::enable_if_t<(std::is_array<ReturnType>::value ||
1377 std::is_pointer<ReturnType>::value)>,
1378 ReturnType, FunctorType> {
1380 Kokkos::HostSpace, Kokkos::MemoryUnmanaged>;
1381
1382 using reducer_type = InvalidType;
1383
1384 using value_type = typename return_type::value_type[];
1385
1386 static return_type return_value(ReturnType& return_val,
1387 const FunctorType& functor) {
1388 if (std::is_array<ReturnType>::value)
1389 return return_type(return_val);
1390 else
1391 return return_type(return_val, functor.value_count);
1392 }
1393};
1394
1395template <class ReturnType, class FunctorType>
1396struct ParallelReduceReturnValue<
1397 std::enable_if_t<Kokkos::is_reducer<ReturnType>::value>, ReturnType,
1398 FunctorType> {
1399 using return_type = ReturnType;
1400 using reducer_type = ReturnType;
1401 using value_type = typename return_type::value_type;
1402
1403 static return_type return_value(ReturnType& return_val, const FunctorType&) {
1404 return return_val;
1405 }
1406};
1407
1408template <class T, class ReturnType, class FunctorType>
1409struct ParallelReducePolicyType;
1410
1411template <class PolicyType, class FunctorType>
1412struct ParallelReducePolicyType<
1413 std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value>,
1414 PolicyType, FunctorType> {
1415 using policy_type = PolicyType;
1416 static PolicyType policy(const PolicyType& policy_) { return policy_; }
1417};
1418
1419template <class PolicyType, class FunctorType>
1420struct ParallelReducePolicyType<
1421 std::enable_if_t<std::is_integral<PolicyType>::value>, PolicyType,
1422 FunctorType> {
1423 using execution_space =
1424 typename Impl::FunctorPolicyExecutionSpace<FunctorType,
1425 void>::execution_space;
1426
1427 using policy_type = Kokkos::RangePolicy<execution_space>;
1428
1429 static policy_type policy(const PolicyType& policy_) {
1430 return policy_type(0, policy_);
1431 }
1432};
1433
1434template <class FunctorType, class ExecPolicy, class ValueType,
1435 class ExecutionSpace>
1436struct ParallelReduceFunctorType {
1437 using functor_type = FunctorType;
1438 static const functor_type& functor(const functor_type& functor) {
1439 return functor;
1440 }
1441};
1442
1443template <class PolicyType, class FunctorType, class ReturnType>
1444struct ParallelReduceAdaptor {
1445 using return_value_adapter =
1446 Impl::ParallelReduceReturnValue<void, ReturnType, FunctorType>;
1447
1448 static inline void execute_impl(const std::string& label,
1449 const PolicyType& policy,
1450 const FunctorType& functor,
1451 ReturnType& return_value) {
1452 uint64_t kpID = 0;
1453
1454 PolicyType inner_policy = policy;
1455 Kokkos::Tools::Impl::begin_parallel_reduce<
1456 typename return_value_adapter::reducer_type>(inner_policy, functor,
1457 label, kpID);
1458
1459 Kokkos::Impl::shared_allocation_tracking_disable();
1460 Impl::ParallelReduce<FunctorType, PolicyType,
1461 typename return_value_adapter::reducer_type>
1462 closure(functor, inner_policy,
1463 return_value_adapter::return_value(return_value, functor));
1464 Kokkos::Impl::shared_allocation_tracking_enable();
1465 closure.execute();
1466
1467 Kokkos::Tools::Impl::end_parallel_reduce<
1468 typename return_value_adapter::reducer_type>(inner_policy, functor,
1469 label, kpID);
1470 }
1471
1472 static constexpr bool is_array_reduction =
1473 Impl::FunctorAnalysis<Impl::FunctorPatternInterface::REDUCE, PolicyType,
1474 FunctorType>::StaticValueSize == 0;
1475
1476 template <typename Dummy = ReturnType>
1477 static inline std::enable_if_t<!(is_array_reduction &&
1478 std::is_pointer<Dummy>::value)>
1479 execute(const std::string& label, const PolicyType& policy,
1480 const FunctorType& functor, ReturnType& return_value) {
1481 execute_impl(label, policy, functor, return_value);
1482 }
1483};
1484} // namespace Impl
1485
1486//----------------------------------------------------------------------------
1487
1499// Parallel Reduce Blocking behavior
1500
1501namespace Impl {
1502template <typename T>
1503struct ReducerHasTestReferenceFunction {
1504 template <typename E>
1505 static std::true_type test_func(decltype(&E::references_scalar));
1506 template <typename E>
1507 static std::false_type test_func(...);
1508
1509 enum {
1510 value = std::is_same<std::true_type, decltype(test_func<T>(nullptr))>::value
1511 };
1512};
1513
1514template <class ExecutionSpace, class Arg>
1515constexpr std::enable_if_t<
1516 // constraints only necessary because SFINAE lacks subsumption
1517 !ReducerHasTestReferenceFunction<Arg>::value &&
1518 !Kokkos::is_view<Arg>::value,
1519 // return type:
1520 bool>
1521parallel_reduce_needs_fence(ExecutionSpace const&, Arg const&) {
1522 return true;
1523}
1524
1525template <class ExecutionSpace, class Reducer>
1526constexpr std::enable_if_t<
1527 // equivalent to:
1528 // (requires (Reducer const& r) {
1529 // { reducer.references_scalar() } -> std::convertible_to<bool>;
1530 // })
1531 ReducerHasTestReferenceFunction<Reducer>::value,
1532 // return type:
1533 bool>
1534parallel_reduce_needs_fence(ExecutionSpace const&, Reducer const& reducer) {
1535 return reducer.references_scalar();
1536}
1537
1538template <class ExecutionSpace, class ViewLike>
1539constexpr std::enable_if_t<
1540 // requires Kokkos::ViewLike<ViewLike>
1541 Kokkos::is_view<ViewLike>::value,
1542 // return type:
1543 bool>
1544parallel_reduce_needs_fence(ExecutionSpace const&, ViewLike const&) {
1545 return false;
1546}
1547
1548template <class ExecutionSpace, class... Args>
1549struct ParallelReduceFence {
1550 template <class... ArgsDeduced>
1551 static void fence(const ExecutionSpace& ex, const std::string& name,
1552 ArgsDeduced&&... args) {
1553 if (Impl::parallel_reduce_needs_fence(ex, (ArgsDeduced &&) args...)) {
1554 ex.fence(name);
1555 }
1556 }
1557};
1558
1559} // namespace Impl
1560
1599// ReturnValue is scalar or array: take by reference
1600
1601template <class PolicyType, class FunctorType, class ReturnType>
1602inline std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value &&
1603 !(Kokkos::is_view<ReturnType>::value ||
1604 Kokkos::is_reducer<ReturnType>::value ||
1605 std::is_pointer<ReturnType>::value)>
1606parallel_reduce(const std::string& label, const PolicyType& policy,
1607 const FunctorType& functor, ReturnType& return_value) {
1608 static_assert(
1609 !std::is_const<ReturnType>::value,
1610 "A const reduction result type is only allowed for a View, pointer or "
1611 "reducer return type!");
1612
1613 Impl::ParallelReduceAdaptor<PolicyType, FunctorType, ReturnType>::execute(
1614 label, policy, functor, return_value);
1615 Impl::ParallelReduceFence<typename PolicyType::execution_space, ReturnType>::
1616 fence(
1617 policy.space(),
1618 "Kokkos::parallel_reduce: fence due to result being value, not view",
1619 return_value);
1620}
1621
1622template <class PolicyType, class FunctorType, class ReturnType>
1623inline std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value &&
1624 !(Kokkos::is_view<ReturnType>::value ||
1625 Kokkos::is_reducer<ReturnType>::value ||
1626 std::is_pointer<ReturnType>::value)>
1627parallel_reduce(const PolicyType& policy, const FunctorType& functor,
1628 ReturnType& return_value) {
1629 static_assert(
1630 !std::is_const<ReturnType>::value,
1631 "A const reduction result type is only allowed for a View, pointer or "
1632 "reducer return type!");
1633
1634 Impl::ParallelReduceAdaptor<PolicyType, FunctorType, ReturnType>::execute(
1635 "", policy, functor, return_value);
1636 Impl::ParallelReduceFence<typename PolicyType::execution_space, ReturnType>::
1637 fence(
1638 policy.space(),
1639 "Kokkos::parallel_reduce: fence due to result being value, not view",
1640 return_value);
1641}
1642
1643template <class FunctorType, class ReturnType>
1644inline std::enable_if_t<!(Kokkos::is_view<ReturnType>::value ||
1645 Kokkos::is_reducer<ReturnType>::value ||
1646 std::is_pointer<ReturnType>::value)>
1647parallel_reduce(const size_t& policy, const FunctorType& functor,
1648 ReturnType& return_value) {
1649 static_assert(
1650 !std::is_const<ReturnType>::value,
1651 "A const reduction result type is only allowed for a View, pointer or "
1652 "reducer return type!");
1653
1654 using policy_type =
1655 typename Impl::ParallelReducePolicyType<void, size_t,
1656 FunctorType>::policy_type;
1657
1658 Impl::ParallelReduceAdaptor<policy_type, FunctorType, ReturnType>::execute(
1659 "", policy_type(0, policy), functor, return_value);
1660 Impl::ParallelReduceFence<typename policy_type::execution_space, ReturnType>::
1661 fence(
1662 typename policy_type::execution_space(),
1663 "Kokkos::parallel_reduce: fence due to result being value, not view",
1664 return_value);
1665}
1666
1667template <class FunctorType, class ReturnType>
1668inline std::enable_if_t<!(Kokkos::is_view<ReturnType>::value ||
1669 Kokkos::is_reducer<ReturnType>::value ||
1670 std::is_pointer<ReturnType>::value)>
1671parallel_reduce(const std::string& label, const size_t& policy,
1672 const FunctorType& functor, ReturnType& return_value) {
1673 static_assert(
1674 !std::is_const<ReturnType>::value,
1675 "A const reduction result type is only allowed for a View, pointer or "
1676 "reducer return type!");
1677
1678 using policy_type =
1679 typename Impl::ParallelReducePolicyType<void, size_t,
1680 FunctorType>::policy_type;
1681 Impl::ParallelReduceAdaptor<policy_type, FunctorType, ReturnType>::execute(
1682 label, policy_type(0, policy), functor, return_value);
1683 Impl::ParallelReduceFence<typename policy_type::execution_space, ReturnType>::
1684 fence(
1685 typename policy_type::execution_space(),
1686 "Kokkos::parallel_reduce: fence due to result being value, not view",
1687 return_value);
1688}
1689
1690// ReturnValue as View or Reducer: take by copy to allow for inline construction
1691
1692template <class PolicyType, class FunctorType, class ReturnType>
1693inline std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value &&
1694 (Kokkos::is_view<ReturnType>::value ||
1695 Kokkos::is_reducer<ReturnType>::value ||
1696 std::is_pointer<ReturnType>::value)>
1697parallel_reduce(const std::string& label, const PolicyType& policy,
1698 const FunctorType& functor, const ReturnType& return_value) {
1699 ReturnType return_value_impl = return_value;
1700 Impl::ParallelReduceAdaptor<PolicyType, FunctorType, ReturnType>::execute(
1701 label, policy, functor, return_value_impl);
1702 Impl::ParallelReduceFence<typename PolicyType::execution_space, ReturnType>::
1703 fence(
1704 policy.space(),
1705 "Kokkos::parallel_reduce: fence due to result being value, not view",
1706 return_value);
1707}
1708
1709template <class PolicyType, class FunctorType, class ReturnType>
1710inline std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value &&
1711 (Kokkos::is_view<ReturnType>::value ||
1712 Kokkos::is_reducer<ReturnType>::value ||
1713 std::is_pointer<ReturnType>::value)>
1714parallel_reduce(const PolicyType& policy, const FunctorType& functor,
1715 const ReturnType& return_value) {
1716 ReturnType return_value_impl = return_value;
1717 Impl::ParallelReduceAdaptor<PolicyType, FunctorType, ReturnType>::execute(
1718 "", policy, functor, return_value_impl);
1719 Impl::ParallelReduceFence<typename PolicyType::execution_space, ReturnType>::
1720 fence(
1721 policy.space(),
1722 "Kokkos::parallel_reduce: fence due to result being value, not view",
1723 return_value);
1724}
1725
1726template <class FunctorType, class ReturnType>
1727inline std::enable_if_t<Kokkos::is_view<ReturnType>::value ||
1728 Kokkos::is_reducer<ReturnType>::value ||
1729 std::is_pointer<ReturnType>::value>
1730parallel_reduce(const size_t& policy, const FunctorType& functor,
1731 const ReturnType& return_value) {
1732 using policy_type =
1733 typename Impl::ParallelReducePolicyType<void, size_t,
1734 FunctorType>::policy_type;
1735 ReturnType return_value_impl = return_value;
1736 Impl::ParallelReduceAdaptor<policy_type, FunctorType, ReturnType>::execute(
1737 "", policy_type(0, policy), functor, return_value_impl);
1738 Impl::ParallelReduceFence<typename policy_type::execution_space, ReturnType>::
1739 fence(
1740 typename policy_type::execution_space(),
1741 "Kokkos::parallel_reduce: fence due to result being value, not view",
1742 return_value);
1743}
1744
1745template <class FunctorType, class ReturnType>
1746inline std::enable_if_t<Kokkos::is_view<ReturnType>::value ||
1747 Kokkos::is_reducer<ReturnType>::value ||
1748 std::is_pointer<ReturnType>::value>
1749parallel_reduce(const std::string& label, const size_t& policy,
1750 const FunctorType& functor, const ReturnType& return_value) {
1751 using policy_type =
1752 typename Impl::ParallelReducePolicyType<void, size_t,
1753 FunctorType>::policy_type;
1754 ReturnType return_value_impl = return_value;
1755 Impl::ParallelReduceAdaptor<policy_type, FunctorType, ReturnType>::execute(
1756 label, policy_type(0, policy), functor, return_value_impl);
1757 Impl::ParallelReduceFence<typename policy_type::execution_space, ReturnType>::
1758 fence(
1759 typename policy_type::execution_space(),
1760 "Kokkos::parallel_reduce: fence due to result being value, not view",
1761 return_value);
1762}
1763
1764// No Return Argument
1765
1766template <class PolicyType, class FunctorType>
1767inline void parallel_reduce(
1768 const std::string& label, const PolicyType& policy,
1769 const FunctorType& functor,
1770 std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value>* =
1771 nullptr) {
1772 using FunctorAnalysis =
1773 Impl::FunctorAnalysis<Impl::FunctorPatternInterface::REDUCE, PolicyType,
1774 FunctorType>;
1775 using value_type = std::conditional_t<(FunctorAnalysis::StaticValueSize != 0),
1776 typename FunctorAnalysis::value_type,
1777 typename FunctorAnalysis::pointer_type>;
1778
1779 static_assert(
1780 FunctorAnalysis::has_final_member_function,
1781 "Calling parallel_reduce without either return value or final function.");
1782
1783 using result_view_type =
1785 result_view_type result_view;
1786
1787 Impl::ParallelReduceAdaptor<PolicyType, FunctorType,
1788 result_view_type>::execute(label, policy, functor,
1789 result_view);
1790}
1791
1792template <class PolicyType, class FunctorType>
1793inline void parallel_reduce(
1794 const PolicyType& policy, const FunctorType& functor,
1795 std::enable_if_t<Kokkos::is_execution_policy<PolicyType>::value>* =
1796 nullptr) {
1797 using FunctorAnalysis =
1798 Impl::FunctorAnalysis<Impl::FunctorPatternInterface::REDUCE, PolicyType,
1799 FunctorType>;
1800 using value_type = std::conditional_t<(FunctorAnalysis::StaticValueSize != 0),
1801 typename FunctorAnalysis::value_type,
1802 typename FunctorAnalysis::pointer_type>;
1803
1804 static_assert(
1805 FunctorAnalysis::has_final_member_function,
1806 "Calling parallel_reduce without either return value or final function.");
1807
1808 using result_view_type =
1810 result_view_type result_view;
1811
1812 Impl::ParallelReduceAdaptor<PolicyType, FunctorType,
1813 result_view_type>::execute("", policy, functor,
1814 result_view);
1815}
1816
1817template <class FunctorType>
1818inline void parallel_reduce(const size_t& policy, const FunctorType& functor) {
1819 using policy_type =
1820 typename Impl::ParallelReducePolicyType<void, size_t,
1821 FunctorType>::policy_type;
1822 using FunctorAnalysis =
1823 Impl::FunctorAnalysis<Impl::FunctorPatternInterface::REDUCE, policy_type,
1824 FunctorType>;
1825 using value_type = std::conditional_t<(FunctorAnalysis::StaticValueSize != 0),
1826 typename FunctorAnalysis::value_type,
1827 typename FunctorAnalysis::pointer_type>;
1828
1829 static_assert(
1830 FunctorAnalysis::has_final_member_function,
1831 "Calling parallel_reduce without either return value or final function.");
1832
1833 using result_view_type =
1835 result_view_type result_view;
1836
1837 Impl::ParallelReduceAdaptor<policy_type, FunctorType,
1838 result_view_type>::execute("",
1839 policy_type(0, policy),
1840 functor, result_view);
1841}
1842
1843template <class FunctorType>
1844inline void parallel_reduce(const std::string& label, const size_t& policy,
1845 const FunctorType& functor) {
1846 using policy_type =
1847 typename Impl::ParallelReducePolicyType<void, size_t,
1848 FunctorType>::policy_type;
1849 using FunctorAnalysis =
1850 Impl::FunctorAnalysis<Impl::FunctorPatternInterface::REDUCE, policy_type,
1851 FunctorType>;
1852 using value_type = std::conditional_t<(FunctorAnalysis::StaticValueSize != 0),
1853 typename FunctorAnalysis::value_type,
1854 typename FunctorAnalysis::pointer_type>;
1855
1856 static_assert(
1857 FunctorAnalysis::has_final_member_function,
1858 "Calling parallel_reduce without either return value or final function.");
1859
1860 using result_view_type =
1862 result_view_type result_view;
1863
1864 Impl::ParallelReduceAdaptor<policy_type, FunctorType,
1865 result_view_type>::execute(label,
1866 policy_type(0, policy),
1867 functor, result_view);
1868}
1869
1870} // namespace Kokkos
1871
1872#endif // KOKKOS_PARALLEL_REDUCE_HPP
View
Memory management for host memory.
Execution policy for work over a range of an integral type.
View to an array of data.
ReturnType