InnerProduct.hpp 3.14 KB
 Praetorius, Simon committed Dec 21, 2020 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 #pragma once #include #include #include #include #include #include #include namespace AMDiS { namespace Recursive { template struct InnerProduct;  Praetorius, Simon committed Dec 21, 2020 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 /** * \brief Recursive inner product of two containers [in1] and [in2] following the * signature of std::inner_product but applied to the leaf elements of the possibly * hierarchic containers. * * This implements the reduction operation of the form *  * T result = init; * for (i...) * result = op1(result, op2(in1[i], in2[i]); * return result; *  **/ template T innerProduct (In1 const& in1, In2 const& in2, T init, BinaryOp1 op1, BinaryOp2 op2)  Praetorius, Simon committed Dec 21, 2020 33 34 35 36 { return InnerProduct::impl(in1, in2, std::move(init), op1, op2); }  Praetorius, Simon committed Dec 21, 2020 37 38 39 40 41 42 43 44 45 46 47 /// Specialization of \ref innerProduct if no functors are given, use plus and multiplies /// by default. /** * This implements the standard Euclidean inner product *  * result = init + sum_i in1[i] * in2[i] *  * * Note, no conjugation of the first or second argument is used in this default. **/  Praetorius, Simon committed Dec 21, 2020 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 template T innerProduct (In1 const& in1, In2 const& in2, T init) { return InnerProduct::impl(in1, in2, std::move(init), std::plus<>{}, std::multiplies<>{}); } /// General implementation of recursive inner-product template struct InnerProduct { private: // dynamic ranges template ())), class = decltype(std::end(std::declval())), class = decltype(std::begin(std::declval()))>  Praetorius, Simon committed Dec 21, 2020 65 66  static T impl2 (Dune::PriorityTag<2>, In1 const& in1, In2 const& in2, T init, BinOp1 op1, BinOp2 op2)  Praetorius, Simon committed Dec 21, 2020 67 68 69 70 71 72 73 74 75 76 77  { auto first1 = std::begin(in1); auto first2 = std::begin(in2); for (; first1 != std::end(in1); ++first1, ++first2) init = Recursive::innerProduct(*first1, *first2, std::move(init), op1, op2); return init; } // ranges with static index access template ()[std::integral_constant{}])>  Praetorius, Simon committed Dec 21, 2020 78 79  static T impl2 (Dune::PriorityTag<1>, In1 const& in1, In2 const& in2, T init, BinOp1 op1, BinOp2 op2)  Praetorius, Simon committed Dec 21, 2020 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104  { static_assert(static_size_v == static_size_v); Ranges::forIndices>([&](auto ii) { init = Recursive::innerProduct(in1[ii], in2[ii], std::move(init), op1, op2); }); return init; } // no range template static T impl2 (Dune::PriorityTag<0>, In1 const& in1, In2 const& in2, T init, BinOp1 op1, BinOp2 op2) { return op1(std::move(init), op2(in1, in2)); } public: template static T impl (In1 const& in1, In2 const& in2, T init, BinOp1 op1, BinOp2 op2) { return impl2(Dune::PriorityTag<5>{}, in1, in2, init, op1, op2); } }; }} // end namespace AMDiS::Recursive