/home/runner/work/DiFfRG_current/DiFfRG_current/DiFfRG/include/DiFfRG/common/complex_math.hh Source File#

DiFfRG: /home/runner/work/DiFfRG_current/DiFfRG_current/DiFfRG/include/DiFfRG/common/complex_math.hh Source File
DiFfRG
complex_math.hh
Go to the documentation of this file.
1#pragma once
2
3// DiFfRG
5
6// external libraries
7#include <autodiff/common/numbertraits.hpp>
8
9namespace DiFfRG
10{
11 using ::Kokkos::complex;
12
13 using ::Kokkos::imag;
14 using ::Kokkos::real;
15} // namespace DiFfRG
16
17// Specialize isArithmetic for complex to make it compatible with real
19{
20 using ::Kokkos::complex;
21 template <typename T> struct ArithmeticTraits<::Kokkos::complex<T>> : ArithmeticTraits<T> {
22 };
23} // namespace autodiff::detail
24
25namespace autodiff::detail
26{
27 // operators for multiplication of float and complex
28 template <typename T1, typename T2>
29 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
30 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const T1 x, const complex<T2> y)
31 {
32 return complex<decltype(T1(1.) * T2(1.))>(x * y.real(), x * y.imag());
33 }
34 template <typename T1, typename T2>
35 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
36 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const complex<T1> &x, const T2 &y)
37 {
38 return complex<decltype(T1(1.) * T2(1.))>(y * x.real(), y * x.imag());
39 }
40 // operators for addition of real and complex
41 template <typename T1, typename T2>
42 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
43 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const T1 &x, const complex<T2> &y)
44 {
45 return complex<decltype(T1(1.) + T2(1.))>(x + y.real(), y.imag());
46 }
47 template <typename T1, typename T2>
48 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
49 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const complex<T1> &x, const T2 &y)
50 {
51 return complex<decltype(T1(1.) + T2(1.))>(x.real() + y, x.imag());
52 }
53 // operators for subtraction of real and complex
54 template <typename T1, typename T2>
55 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
56 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const T1 &x, const complex<T2> &y)
57 {
58 return complex<decltype(T1(1.) - T2(1.))>(x - y.real(), -y.imag());
59 }
60 template <typename T1, typename T2>
61 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
62 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const complex<T1> &x, const T2 &y)
63 {
64 return complex<decltype(T1(1.) - T2(1.))>(x.real() - y, x.imag());
65 }
66 // operators for division of real and complex
67 template <typename T1, typename T2>
68 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
69 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const T1 &x, const complex<T2> &y)
70 {
71 return complex<decltype(T1(1.) / T2(1.))>(x * y.real(), -x * y.imag()) / (powr<2>(y.real()) + powr<2>(y.imag()));
72 }
73 template <typename T1, typename T2>
74 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
75 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const complex<T1> x, const T2 &y)
76 {
77 return complex<decltype(T1(1.) / T2(1.))>(x.real() / y, x.imag() / y);
78 }
79} // namespace autodiff::detail
80
81// external libraries
82#include <autodiff/forward/real.hpp>
83
84namespace DiFfRG
85{
86 template <size_t N, typename T> using cxReal = autodiff::Real<N, complex<T>>;
87 using cxreal = autodiff::Real<1, complex<double>>;
88
89 template <typename T> struct is_complex : public std::false_type {
90 };
91 template <typename T> struct is_complex<complex<T>> : public std::true_type {
92 };
93 template <size_t N, typename T> struct is_complex<cxReal<N, T>> : public std::true_type {
94 };
95
96 template <size_t N, typename T> KOKKOS_FORCEINLINE_FUNCTION auto real(const autodiff::Real<N, T> &a) { return a; }
97 template <size_t N, typename T> constexpr KOKKOS_FORCEINLINE_FUNCTION auto imag(const autodiff::Real<N, T> &)
98 {
99 return 0.;
100 }
101
102 template <size_t N, typename T> KOKKOS_FORCEINLINE_FUNCTION auto real(const cxReal<N, T> &x)
103 {
104 autodiff::Real<N, T> res;
105 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = real(x[i]); });
106 return res;
107 }
108 template <size_t N, typename T> KOKKOS_FORCEINLINE_FUNCTION auto imag(const cxReal<N, T> &x)
109 {
110 autodiff::Real<N, T> res;
111 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = imag(x[i]); });
112 return res;
113 }
114
115 // operators for multiplication of real and complex
116 template <size_t N, typename T>
117 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const autodiff::Real<N, T> &x, const complex<double> &y)
118 {
119 cxReal<N, T> res;
120 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = x[i] * y; });
121 return res;
122 }
123 template <size_t N, typename T>
124 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const complex<double> &x, const autodiff::Real<N, T> &y)
125 {
126 cxReal<N, T> res;
127 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = x * y[i]; });
128 return res;
129 }
130 // operators for addition of real and complex
131 template <size_t N, typename T>
132 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const autodiff::Real<N, T> &x, const complex<double> &y)
133 {
134 cxReal<N, T> res;
135 res[0] = x[0] + y;
136 autodiff::detail::For<1, N + 1>([&](auto i) constexpr { res[i] = x[i]; });
137 return res;
138 }
139 template <size_t N, typename T>
140 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const complex<double> &x, const autodiff::Real<N, T> &y)
141 {
142 cxReal<N, T> res;
143 res[0] = x + y[0];
144 autodiff::detail::For<1, N + 1>([&](auto i) constexpr { res[i] = y[i]; });
145 return res;
146 }
147 // operators for subtraction of real and complex
148 template <size_t N, typename T>
149 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const autodiff::Real<N, T> &x, const complex<double> &y)
150 {
151 cxReal<N, T> res;
152 res[0] = x[0] - y;
153 autodiff::detail::For<1, N + 1>([&](auto i) constexpr { res[i] = x[i]; });
154 return res;
155 }
156 template <size_t N, typename T>
157 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const complex<double> &x, const autodiff::Real<N, T> &y)
158 {
159 cxReal<N, T> res;
160 res[0] = x - y[0];
161 autodiff::detail::For<1, N + 1>([&](auto i) constexpr { res[i] = -y[i]; });
162 return res;
163 }
164 // operators for division of real and complex
165 template <size_t N, typename T>
166 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const autodiff::Real<N, T> &x, const complex<double> &y)
167 {
168 cxReal<N, T> res;
169 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = x[i] / y; });
170 return res;
171 }
172 template <size_t N, typename T>
173 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const complex<double> &x, const autodiff::Real<N, T> &y)
174 {
175 cxReal<N, T> ux(x);
176 cxReal<N, T> uy(y);
177 return ux /= uy;
178 }
179
180 // operators for multiplication of real and cxReal
181 template <size_t N, typename T>
182 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const autodiff::Real<N, T> &x, const cxReal<N, T> &y)
183 {
184 cxReal<N, T> res(x);
185 return res *= y;
186 }
187 template <size_t N, typename T>
188 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const cxReal<N, T> &x, const autodiff::Real<N, T> &y)
189 {
190 cxReal<N, T> res(y);
191 return res *= x;
192 }
193 // operators for addition of real and cxReal
194 template <size_t N, typename T>
195 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const autodiff::Real<N, T> &x, const cxReal<N, T> &y)
196 {
197 cxReal<N, T> res(x);
198 return res += y;
199 }
200 template <size_t N, typename T>
201 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const cxReal<N, T> &x, const autodiff::Real<N, T> &y)
202 {
203 cxReal<N, T> res(y);
204 return res += x;
205 }
206 // operators for subtraction of real and cxReal
207 template <size_t N, typename T>
208 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const autodiff::Real<N, T> &x, const cxReal<N, T> &y)
209 {
210 cxReal<N, T> res(x);
211 return res -= y;
212 }
213 template <size_t N, typename T>
214 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const cxReal<N, T> &x, const autodiff::Real<N, T> &y)
215 {
216 cxReal<N, T> ux(x);
217 cxReal<N, T> uy(y);
218 return ux -= uy;
219 }
220 // operators for division of real and cxReal
221 template <size_t N, typename T>
222 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const autodiff::Real<N, T> &x, const cxReal<N, T> &y)
223 {
224 cxReal<N, T> ux(x);
225 return ux /= y;
226 }
227 template <size_t N, typename T>
228 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const cxReal<N, T> &x, const autodiff::Real<N, T> &y)
229 {
230 cxReal<N, T> uy(y);
231 return x / uy;
232 }
233
234 // operators for multiplication of complex and cxReal
235 template <size_t N, typename T>
236 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const complex<double> &x, const cxReal<N, T> &y)
237 {
238 cxReal<N, T> res;
239 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = y[i] * x; });
240 return res;
241 }
242 template <size_t N, typename T>
243 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const cxReal<N, T> &x, const complex<double> &y)
244 {
245 cxReal<N, T> res;
246 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = x[i] * y; });
247 return res;
248 }
249 // operators for addition of complex and cxReal
250 template <size_t N, typename T>
251 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const complex<double> &x, const cxReal<N, T> &y)
252 {
253 cxReal<N, T> res(y);
254 res[0] += x;
255 return res;
256 }
257 template <size_t N, typename T>
258 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const cxReal<N, T> &x, const complex<double> &y)
259 {
260 cxReal<N, T> res(x);
261 res[0] += y;
262 return res;
263 }
264 // operators for subtraction of complex and cxReal
265 template <size_t N, typename T>
266 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const complex<double> &x, const cxReal<N, T> &y)
267 {
268 cxReal<N, T> ux(x);
269 return ux -= y;
270 }
271 template <size_t N, typename T>
272 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const cxReal<N, T> &x, const complex<double> &y)
273 {
274 cxReal<N, T> res(x);
275 res[0] -= y;
276 return res;
277 }
278 // operators for division of complex and cxReal
279 template <size_t N, typename T>
280 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const cxReal<N, T> &x, const complex<double> &y)
281 {
282 cxReal<N, T> res;
283 autodiff::detail::For<0, N + 1>([&](auto i) constexpr { res[i] = x[i] / y; });
284 return res;
285 }
286 template <size_t N, typename T>
287 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const complex<double> &x, const cxReal<N, T> &y)
288 {
289 cxReal<N, T> ux(x);
290 cxReal<N, T> uy(y);
291 return ux /= uy;
292 }
293
294 // operators for division of double and cxReal
295 template <size_t N, typename T> KOKKOS_FORCEINLINE_FUNCTION auto operator/(const double x, const cxReal<N, T> &y)
296 {
297 cxReal<N, T> res;
298 res[0] = x;
299 autodiff::detail::For<N + 1>([&](auto i) constexpr {
300 res[i] -= autodiff::detail::Sum<0, i>([&](auto j) constexpr {
301 constexpr auto c = autodiff::detail::BinomialCoefficient<i.index, j.index>;
302 return c * res[j] * y[i - j];
303 });
304 res[i] /= y[0];
305 });
306 return res;
307 }
308
309 // ------------------------------------------------------------------
310 // operators for complex<T1> and T2
311 // ------------------------------------------------------------------
312
313 // operators for multiplication of float and complex
314 /*
315 template <typename T1, typename T2>
316 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
317 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const T1 x, const complex<T2> y)
318 {
319 return complex<decltype(T1(1.) * T2(1.))>(x * y.real(), x * y.imag());
320 }
321 template <typename T1, typename T2>
322 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
323 KOKKOS_FORCEINLINE_FUNCTION auto operator*(const complex<T1> &x, const T2 &y)
324 {
325 return complex<decltype(T1(1.) * T2(1.))>(y * x.real(), y * x.imag());
326 }
327 // operators for addition of real and complex
328 template <typename T1, typename T2>
329 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
330 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const T1 &x, const complex<T2> &y)
331 {
332 return complex<decltype(T1(1.) + T2(1.))>(x + y.real(), y.imag());
333 }
334 template <typename T1, typename T2>
335 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
336 KOKKOS_FORCEINLINE_FUNCTION auto operator+(const complex<T1> &x, const T2 &y)
337 {
338 return complex<decltype(T1(1.) + T2(1.))>(x.real() + y, x.imag());
339 }
340 // operators for subtraction of real and complex
341 template <typename T1, typename T2>
342 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
343 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const T1 &x, const complex<T2> &y)
344 {
345 return complex<decltype(T1(1.) + T2(1.))>(x - y.real(), -y.imag());
346 }
347 template <typename T1, typename T2>
348 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
349 KOKKOS_FORCEINLINE_FUNCTION auto operator-(const complex<T1> &x, const T2 &y)
350 {
351 return complex<decltype(T1(1.) + T2(1.))>(x.real() - y, x.imag());
352 }
353 // operators for division of real and complex
354 template <typename T1, typename T2>
355 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
356 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const T1 &x, const complex<T2> &y)
357 {
358 return complex<decltype(T1(1.) + T2(1.))>(x * y.real(), -x * y.imag()) / (powr<2>(y.real()) + powr<2>(y.imag()));
359 }
360 template <typename T1, typename T2>
361 requires(!std::is_same_v<T1, T2>) && std::is_arithmetic_v<T1> && std::is_arithmetic_v<T2>
362 KOKKOS_FORCEINLINE_FUNCTION auto operator/(const complex<T1> x, const T2 &y)
363 {
364 return complex<decltype(T1(1.) + T2(1.))>(x.real() / y, x.imag() / y);
365 }
366*/
367} // namespace DiFfRG
Definition complex_math.hh:10
KOKKOS_FORCEINLINE_FUNCTION auto operator/(const autodiff::Real< N, T > &x, const complex< double > &y)
Definition complex_math.hh:166
KOKKOS_FORCEINLINE_FUNCTION auto operator+(const autodiff::Real< N, T > &x, const complex< double > &y)
Definition complex_math.hh:132
autodiff::Real< 1, complex< double > > cxreal
Definition complex_math.hh:87
constexpr KOKKOS_FORCEINLINE_FUNCTION auto imag(const autodiff::Real< N, T > &)
Definition complex_math.hh:97
KOKKOS_FORCEINLINE_FUNCTION auto operator-(const autodiff::Real< N, T > &x, const complex< double > &y)
Definition complex_math.hh:149
KOKKOS_FORCEINLINE_FUNCTION auto operator*(const autodiff::Real< N, T > &x, const complex< double > &y)
Definition complex_math.hh:117
KOKKOS_FORCEINLINE_FUNCTION auto real(const autodiff::Real< N, T > &a)
Definition complex_math.hh:96
autodiff::Real< N, complex< T > > cxReal
Definition complex_math.hh:86
Definition complex_math.hh:19
Definition complex_math.hh:89