Line data Source code
1 1 : // Distributed under the MIT License.
2 : // See LICENSE.txt for details.
3 :
4 : /// \file
5 : /// Defines functions euclidean dot_product and dot_product with a metric
6 :
7 : #pragma once
8 :
9 : #include <cstddef>
10 :
11 : #include "DataStructures/Tensor/Tensor.hpp"
12 : #include "Utilities/ContainerHelpers.hpp"
13 : #include "Utilities/Kokkos/KokkosMacros.hpp"
14 : #include "Utilities/MakeWithValue.hpp"
15 :
16 : /// @{
17 : /*!
18 : * \ingroup TensorGroup
19 : * \brief Compute the Euclidean dot product of two vectors or one forms
20 : *
21 : * \details
22 : * Returns \f$A^a B^b \delta_{ab}\f$ for input vectors \f$A^a\f$ and \f$B^b\f$
23 : * or \f$A_a B_b \delta^{ab}\f$ for input one forms \f$A_a\f$ and \f$B_b\f$.
24 : *
25 : * These functions can be called in Kokkos kernels with `double`s.
26 : */
27 : SPECTRE_KOKKOS_DISABLE_EXEC_CHECK
28 : template <typename DataTypeLhs, typename DataTypeRhs, typename Index,
29 : typename DataTypeResult = decltype(blaze::evaluate(DataTypeLhs() *
30 : DataTypeRhs()))>
31 1 : KOKKOS_FUNCTION void dot_product(
32 : const gsl::not_null<Scalar<DataTypeResult>*> dot_product,
33 : const Tensor<DataTypeLhs, Symmetry<1>, index_list<Index>>& vector_a,
34 : const Tensor<DataTypeRhs, Symmetry<1>, index_list<Index>>& vector_b) {
35 : get(*dot_product) = get<0>(vector_a) * get<0>(vector_b);
36 : for (size_t d = 1; d < Index::dim; ++d) {
37 : get(*dot_product) += vector_a.get(d) * vector_b.get(d);
38 : }
39 : }
40 :
41 : SPECTRE_KOKKOS_DISABLE_EXEC_CHECK
42 : template <typename DataTypeLhs, typename DataTypeRhs, typename Index,
43 : typename DataTypeResult = decltype(blaze::evaluate(DataTypeLhs() *
44 : DataTypeRhs()))>
45 1 : KOKKOS_FUNCTION Scalar<DataTypeResult> dot_product(
46 : const Tensor<DataTypeLhs, Symmetry<1>, index_list<Index>>& vector_a,
47 : const Tensor<DataTypeRhs, Symmetry<1>, index_list<Index>>& vector_b) {
48 : Scalar<DataTypeResult> dot_product{};
49 : ::dot_product(make_not_null(&dot_product), vector_a, vector_b);
50 : return dot_product;
51 : }
52 : /// @}
53 :
54 : /// @{
55 : /*!
56 : * \ingroup TensorGroup
57 : * \brief Compute the dot product of a vector and a one form
58 : *
59 : * \details
60 : * Returns \f$A^a B_b \delta_{a}^b\f$ for input vector \f$A^a\f$ and
61 : * input one form \f$B_b\f$
62 : * or \f$A_a B^b \delta^a_b\f$ for input one form \f$A_a\f$ and
63 : * input vector \f$B^b\f$.
64 : *
65 : * These functions can be called in Kokkos kernels with `double`s.
66 : */
67 : SPECTRE_KOKKOS_DISABLE_EXEC_CHECK
68 : template <typename DataTypeLhs, typename DataTypeRhs, typename Index,
69 : typename DataTypeResult = decltype(blaze::evaluate(DataTypeLhs() *
70 : DataTypeRhs()))>
71 1 : KOKKOS_FUNCTION void dot_product(
72 : // NOLINTNEXTLINE(readability-avoid-const-params-in-decls) false positive
73 : const gsl::not_null<Scalar<DataTypeResult>*> dot_product,
74 : const Tensor<DataTypeLhs, Symmetry<1>, index_list<Index>>& vector_a,
75 : const Tensor<DataTypeRhs, Symmetry<1>,
76 : index_list<change_index_up_lo<Index>>>& vector_b) {
77 : get(*dot_product) = get<0>(vector_a) * get<0>(vector_b);
78 : for (size_t d = 1; d < Index::dim; ++d) {
79 : get(*dot_product) += vector_a.get(d) * vector_b.get(d);
80 : }
81 : }
82 :
83 : SPECTRE_KOKKOS_DISABLE_EXEC_CHECK
84 : template <typename DataTypeLhs, typename DataTypeRhs, typename Index,
85 : typename DataTypeResult = decltype(blaze::evaluate(DataTypeLhs() *
86 : DataTypeRhs()))>
87 1 : KOKKOS_FUNCTION Scalar<DataTypeResult> dot_product(
88 : const Tensor<DataTypeLhs, Symmetry<1>, index_list<Index>>& vector_a,
89 : const Tensor<DataTypeRhs, Symmetry<1>,
90 : index_list<change_index_up_lo<Index>>>& vector_b) {
91 : Scalar<DataTypeResult> dot_product{};
92 : ::dot_product(make_not_null(&dot_product), vector_a, vector_b);
93 : return dot_product;
94 : }
95 : /// @}
96 :
97 : /// @{
98 : /*!
99 : * \ingroup TensorGroup
100 : * \brief Compute the dot_product of two vectors or one forms
101 : *
102 : * \details
103 : * Returns \f$g_{ab} A^a B^b\f$, where \f$g_{ab}\f$ is the metric,
104 : * \f$A^a\f$ is vector_a, and \f$B^b\f$ is vector_b.
105 : * Or, returns \f$g^{ab} A_a B_b\f$ when given one forms \f$A_a\f$
106 : * and \f$B_b\f$ with an inverse metric \f$g^{ab}\f$.
107 : *
108 : * These functions can be called in Kokkos kernels with `double`s.
109 : */
110 : SPECTRE_KOKKOS_DISABLE_EXEC_CHECK
111 : template <typename DataTypeLhs, typename DataTypeRhs, typename DataTypeMetric,
112 : typename Index,
113 : typename DataTypeResult = decltype(blaze::evaluate(
114 : DataTypeLhs() * DataTypeRhs() * DataTypeMetric()))>
115 1 : KOKKOS_FUNCTION void dot_product(
116 : const gsl::not_null<Scalar<DataTypeResult>*> dot_product,
117 : const Tensor<DataTypeLhs, Symmetry<1>, index_list<Index>>& vector_a,
118 : const Tensor<DataTypeRhs, Symmetry<1>, index_list<Index>>& vector_b,
119 : const Tensor<DataTypeMetric, Symmetry<1, 1>,
120 : index_list<change_index_up_lo<Index>,
121 : change_index_up_lo<Index>>>& metric) {
122 : if constexpr (Index::dim == 1) {
123 : get(*dot_product) = get<0>(vector_a) * get<0>(vector_b) * get<0, 0>(metric);
124 : } else {
125 : if (&vector_a == &vector_b) {
126 : get(*dot_product) =
127 : get<0>(vector_a) * get<1>(vector_a) * get<0, 1>(metric);
128 : for (size_t j = 2; j < Index::dim; ++j) {
129 : get(*dot_product) +=
130 : vector_a.get(0) * vector_a.get(j) * metric.get(0, j);
131 : }
132 : for (size_t i = 1; i < Index::dim; ++i) {
133 : for (size_t j = i + 1; j < Index::dim; ++j) {
134 : get(*dot_product) +=
135 : vector_a.get(i) * vector_a.get(j) * metric.get(i, j);
136 : }
137 : }
138 : get(*dot_product) *= 2.0;
139 :
140 : for (size_t i = 0; i < Index::dim; ++i) {
141 : get(*dot_product) += square(vector_a.get(i)) * metric.get(i, i);
142 : }
143 : } else {
144 : get(*dot_product) =
145 : get<0>(vector_a) * get<0>(vector_b) * get<0, 0>(metric);
146 : for (size_t b = 1; b < Index::dim; ++b) {
147 : get(*dot_product) +=
148 : get<0>(vector_a) * vector_b.get(b) * metric.get(0, b);
149 : }
150 :
151 : for (size_t a = 1; a < Index::dim; ++a) {
152 : for (size_t b = 0; b < Index::dim; ++b) {
153 : get(*dot_product) +=
154 : vector_a.get(a) * vector_b.get(b) * metric.get(a, b);
155 : }
156 : }
157 : }
158 : }
159 : }
160 :
161 : SPECTRE_KOKKOS_DISABLE_EXEC_CHECK
162 : template <typename DataTypeLhs, typename DataTypeRhs, typename DataTypeMetric,
163 : typename Index,
164 : typename DataTypeResult = decltype(blaze::evaluate(
165 : DataTypeLhs() * DataTypeRhs() * DataTypeMetric()))>
166 1 : KOKKOS_FUNCTION Scalar<DataTypeResult> dot_product(
167 : const Tensor<DataTypeLhs, Symmetry<1>, index_list<Index>>& vector_a,
168 : const Tensor<DataTypeRhs, Symmetry<1>, index_list<Index>>& vector_b,
169 : const Tensor<DataTypeMetric, Symmetry<1, 1>,
170 : index_list<change_index_up_lo<Index>,
171 : change_index_up_lo<Index>>>& metric) {
172 : Scalar<DataTypeResult> dot_product{};
173 : ::dot_product(make_not_null(&dot_product), vector_a, vector_b, metric);
174 : return dot_product;
175 : }
176 : /// @}
|
