Eigen  5.0.1
 
Loading...
Searching...
No Matches
Complex.h
1// This file is part of Eigen, a lightweight C++ template library
2// for linear algebra.
3//
4// Copyright (C) 2018 Wave Computing, Inc.
5// Written by:
6// Chris Larsen
7// Alexey Frunze (afrunze@wavecomp.com)
8//
9// This Source Code Form is subject to the terms of the Mozilla
10// Public License v. 2.0. If a copy of the MPL was not distributed
11// with this file, You can obtain one at http://mozilla.org/MPL/2.0/.
12// SPDX-License-Identifier: MPL-2.0
13
14#ifndef EIGEN_COMPLEX_MSA_H
15#define EIGEN_COMPLEX_MSA_H
16
17#include <iostream>
18
19// IWYU pragma: private
20#include "../../InternalHeaderCheck.h"
21
22namespace Eigen {
23
24namespace internal {
25
26EIGEN_GCC_FAST_MATH_COMPLEX_VECTORIZE_WORKAROUND_PUSH
27
28//---------- float ----------
29struct Packet2cf {
30 EIGEN_STRONG_INLINE Packet2cf() {}
31 EIGEN_STRONG_INLINE explicit Packet2cf(const std::complex<float>& a, const std::complex<float>& b) {
32 Packet4f t = {std::real(a), std::imag(a), std::real(b), std::imag(b)};
33 v = t;
34 }
35 EIGEN_STRONG_INLINE explicit Packet2cf(const Packet4f& a) : v(a) {}
36 EIGEN_STRONG_INLINE Packet2cf(const Packet2cf& a) : v(a.v) {}
37 EIGEN_STRONG_INLINE Packet2cf& operator=(const Packet2cf& b) {
38 v = b.v;
39 return *this;
40 }
41 EIGEN_STRONG_INLINE Packet2cf conjugate(void) const {
42 return Packet2cf((Packet4f)__builtin_msa_bnegi_d((v2u64)v, 63));
43 }
44 EIGEN_STRONG_INLINE Packet2cf& operator*=(const Packet2cf& b) {
45 Packet4f v1, v2;
46
47 // Get the real values of a | a1_re | a1_re | a2_re | a2_re |
48 v1 = (Packet4f)__builtin_msa_ilvev_w((v4i32)v, (v4i32)v);
49 // Get the imag values of a | a1_im | a1_im | a2_im | a2_im |
50 v2 = (Packet4f)__builtin_msa_ilvod_w((v4i32)v, (v4i32)v);
51 // Multiply the real a with b
52 v1 = pmul(v1, b.v);
53 // Multiply the imag a with b
54 v2 = pmul(v2, b.v);
55 // Conjugate v2
56 v2 = Packet2cf(v2).conjugate().v;
57 // Swap real/imag elements in v2.
58 v2 = (Packet4f)__builtin_msa_shf_w((v4i32)v2, EIGEN_MSA_SHF_I8(1, 0, 3, 2));
59 // Add and return the result
60 v = padd(v1, v2);
61 return *this;
62 }
63 EIGEN_STRONG_INLINE Packet2cf operator*(const Packet2cf& b) const { return Packet2cf(*this) *= b; }
64 EIGEN_STRONG_INLINE Packet2cf& operator+=(const Packet2cf& b) {
65 v = padd(v, b.v);
66 return *this;
67 }
68 EIGEN_STRONG_INLINE Packet2cf operator+(const Packet2cf& b) const { return Packet2cf(*this) += b; }
69 EIGEN_STRONG_INLINE Packet2cf& operator-=(const Packet2cf& b) {
70 v = psub(v, b.v);
71 return *this;
72 }
73 EIGEN_STRONG_INLINE Packet2cf operator-(const Packet2cf& b) const { return Packet2cf(*this) -= b; }
74 EIGEN_STRONG_INLINE Packet2cf operator/(const Packet2cf& b) const { return pdiv_complex(Packet2cf(*this), b); }
75 EIGEN_STRONG_INLINE Packet2cf& operator/=(const Packet2cf& b) {
76 *this = Packet2cf(*this) / b;
77 return *this;
78 }
79 EIGEN_STRONG_INLINE Packet2cf operator-(void) const { return Packet2cf(pnegate(v)); }
80
81 Packet4f v;
82};
83
84inline std::ostream& operator<<(std::ostream& os, const Packet2cf& value) {
85 os << "[ (" << value.v[0] << ", " << value.v[1]
86 << "i),"
87 " ("
88 << value.v[2] << ", " << value.v[3] << "i) ]";
89 return os;
90}
91
92template <>
93struct packet_traits<std::complex<float> > : default_packet_traits {
94 typedef Packet2cf type;
95 typedef Packet2cf half;
96 enum {
97 Vectorizable = 1,
98 AlignedOnScalar = 1,
99 size = 2,
100
101 HasAdd = 1,
102 HasSub = 1,
103 HasMul = 1,
104 HasDiv = 1,
105 HasNegate = 1,
106 HasSqrt = 1,
107 HasLog = 1,
108 HasExp = 1,
109 HasAbs = 0,
110 HasAbs2 = 0,
111 HasMin = 0,
112 HasMax = 0,
113 HasSetLinear = 0,
114 };
115};
116
117template <>
118struct unpacket_traits<Packet2cf> {
119 typedef std::complex<float> type;
120 typedef Packet4f as_real;
121 enum {
122 size = 2,
123 alignment = Aligned16,
124 vectorizable = true,
125 masked_load_available = false,
126 masked_store_available = false
127 };
128 typedef Packet2cf half;
129};
130
131template <>
132EIGEN_STRONG_INLINE Packet2cf pset1<Packet2cf>(const std::complex<float>& from) {
133 EIGEN_MSA_DEBUG;
134
135 float f0 = from.real(), f1 = from.imag();
136 Packet4f v0 = {f0, f0, f0, f0};
137 Packet4f v1 = {f1, f1, f1, f1};
138 return Packet2cf((Packet4f)__builtin_msa_ilvr_w((Packet4i)v1, (Packet4i)v0));
139}
140
141template <>
142EIGEN_STRONG_INLINE Packet2cf padd<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
143 EIGEN_MSA_DEBUG;
144
145 return a + b;
146}
147
148template <>
149EIGEN_STRONG_INLINE Packet2cf psub<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
150 EIGEN_MSA_DEBUG;
151
152 return a - b;
153}
154
155template <>
156EIGEN_STRONG_INLINE Packet2cf pnegate(const Packet2cf& a) {
157 EIGEN_MSA_DEBUG;
158
159 return -a;
160}
161
162template <>
163EIGEN_STRONG_INLINE Packet2cf pconj(const Packet2cf& a) {
164 EIGEN_MSA_DEBUG;
165
166 return a.conjugate();
167}
168
169template <>
170EIGEN_STRONG_INLINE Packet2cf pmul<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
171 EIGEN_MSA_DEBUG;
172
173 return a * b;
174}
175
176template <>
177EIGEN_STRONG_INLINE Packet2cf pand<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
178 EIGEN_MSA_DEBUG;
179
180 return Packet2cf(pand(a.v, b.v));
181}
182
183template <>
184EIGEN_STRONG_INLINE Packet2cf por<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
185 EIGEN_MSA_DEBUG;
186
187 return Packet2cf(por(a.v, b.v));
188}
189
190template <>
191EIGEN_STRONG_INLINE Packet2cf pxor<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
192 EIGEN_MSA_DEBUG;
193
194 return Packet2cf(pxor(a.v, b.v));
195}
196
197template <>
198EIGEN_STRONG_INLINE Packet2cf pandnot<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
199 EIGEN_MSA_DEBUG;
200
201 return Packet2cf(pandnot(a.v, b.v));
202}
203
204template <>
205EIGEN_STRONG_INLINE Packet2cf pload<Packet2cf>(const std::complex<float>* from) {
206 EIGEN_MSA_DEBUG;
207
208 EIGEN_DEBUG_ALIGNED_LOAD return Packet2cf(pload<Packet4f>((const float*)from));
209}
210
211template <>
212EIGEN_STRONG_INLINE Packet2cf ploadu<Packet2cf>(const std::complex<float>* from) {
213 EIGEN_MSA_DEBUG;
214
215 EIGEN_DEBUG_UNALIGNED_LOAD return Packet2cf(ploadu<Packet4f>((const float*)from));
216}
217
218template <>
219EIGEN_STRONG_INLINE Packet2cf ploaddup<Packet2cf>(const std::complex<float>* from) {
220 EIGEN_MSA_DEBUG;
221
222 return pset1<Packet2cf>(*from);
223}
224
225template <>
226EIGEN_STRONG_INLINE void pstore<std::complex<float> >(std::complex<float>* to, const Packet2cf& from) {
227 EIGEN_MSA_DEBUG;
228
229 EIGEN_DEBUG_ALIGNED_STORE pstore<float>((float*)to, from.v);
230}
231
232template <>
233EIGEN_STRONG_INLINE void pstoreu<std::complex<float> >(std::complex<float>* to, const Packet2cf& from) {
234 EIGEN_MSA_DEBUG;
235
236 EIGEN_DEBUG_UNALIGNED_STORE pstoreu<float>((float*)to, from.v);
237}
238
239template <>
240EIGEN_DEVICE_FUNC inline Packet2cf pgather<std::complex<float>, Packet2cf>(const std::complex<float>* from,
241 Index stride) {
242 EIGEN_MSA_DEBUG;
243
244 return Packet2cf(from[0 * stride], from[1 * stride]);
245}
246
247template <>
248EIGEN_DEVICE_FUNC inline void pscatter<std::complex<float>, Packet2cf>(std::complex<float>* to, const Packet2cf& from,
249 Index stride) {
250 EIGEN_MSA_DEBUG;
251
252 *to = std::complex<float>(from.v[0], from.v[1]);
253 to += stride;
254 *to = std::complex<float>(from.v[2], from.v[3]);
255}
256
257template <>
258EIGEN_STRONG_INLINE void prefetch<std::complex<float> >(const std::complex<float>* addr) {
259 EIGEN_MSA_DEBUG;
260
261 prefetch(reinterpret_cast<const float*>(addr));
262}
263
264template <>
265EIGEN_STRONG_INLINE std::complex<float> pfirst<Packet2cf>(const Packet2cf& a) {
266 EIGEN_MSA_DEBUG;
267
268 return std::complex<float>(a.v[0], a.v[1]);
269}
270
271template <>
272EIGEN_STRONG_INLINE Packet2cf preverse(const Packet2cf& a) {
273 EIGEN_MSA_DEBUG;
274
275 return Packet2cf((Packet4f)__builtin_msa_shf_w((v4i32)a.v, EIGEN_MSA_SHF_I8(2, 3, 0, 1)));
276}
277
278template <>
279EIGEN_STRONG_INLINE Packet2cf pcplxflip<Packet2cf>(const Packet2cf& a) {
280 EIGEN_MSA_DEBUG;
281
282 return Packet2cf((Packet4f)__builtin_msa_shf_w((v4i32)a.v, EIGEN_MSA_SHF_I8(1, 0, 3, 2)));
283}
284
285template <>
286EIGEN_STRONG_INLINE std::complex<float> predux<Packet2cf>(const Packet2cf& a) {
287 EIGEN_MSA_DEBUG;
288
289 Packet4f value = (Packet4f)preverse((Packet2d)a.v);
290 value += a.v;
291 return std::complex<float>(value[0], value[1]);
292}
293
294template <>
295EIGEN_STRONG_INLINE std::complex<float> predux_mul<Packet2cf>(const Packet2cf& a) {
296 EIGEN_MSA_DEBUG;
297
298 return std::complex<float>((a.v[0] * a.v[2]) - (a.v[1] * a.v[3]), (a.v[0] * a.v[3]) + (a.v[1] * a.v[2]));
299}
300
301EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet2cf, Packet4f)
302
303template <>
304EIGEN_STRONG_INLINE Packet2cf pdiv<Packet2cf>(const Packet2cf& a, const Packet2cf& b) {
305 EIGEN_MSA_DEBUG;
306
307 return a / b;
308}
309
310inline std::ostream& operator<<(std::ostream& os, const PacketBlock<Packet2cf, 2>& value) {
311 os << "[ " << value.packet[0] << ", " << std::endl << " " << value.packet[1] << " ]";
312 return os;
313}
314
315EIGEN_DEVICE_FUNC inline void ptranspose(PacketBlock<Packet2cf, 2>& kernel) {
316 EIGEN_MSA_DEBUG;
317
318 Packet4f tmp = (Packet4f)__builtin_msa_ilvl_d((v2i64)kernel.packet[1].v, (v2i64)kernel.packet[0].v);
319 kernel.packet[0].v = (Packet4f)__builtin_msa_ilvr_d((v2i64)kernel.packet[1].v, (v2i64)kernel.packet[0].v);
320 kernel.packet[1].v = tmp;
321}
322
323//---------- double ----------
324
325struct Packet1cd {
326 EIGEN_STRONG_INLINE Packet1cd() {}
327 EIGEN_STRONG_INLINE explicit Packet1cd(const std::complex<double>& a) {
328 v[0] = std::real(a);
329 v[1] = std::imag(a);
330 }
331 EIGEN_STRONG_INLINE explicit Packet1cd(const Packet2d& a) : v(a) {}
332 EIGEN_STRONG_INLINE Packet1cd(const Packet1cd& a) : v(a.v) {}
333 EIGEN_STRONG_INLINE Packet1cd& operator=(const Packet1cd& b) {
334 v = b.v;
335 return *this;
336 }
337 EIGEN_STRONG_INLINE Packet1cd conjugate(void) const {
338 static const v2u64 p2ul_CONJ_XOR = {0x0, 0x8000000000000000};
339 return (Packet1cd)pxor(v, (Packet2d)p2ul_CONJ_XOR);
340 }
341 EIGEN_STRONG_INLINE Packet1cd& operator*=(const Packet1cd& b) {
342 Packet2d v1, v2;
343
344 // Get the real values of a | a1_re | a1_re
345 v1 = (Packet2d)__builtin_msa_ilvev_d((v2i64)v, (v2i64)v);
346 // Get the imag values of a | a1_im | a1_im
347 v2 = (Packet2d)__builtin_msa_ilvod_d((v2i64)v, (v2i64)v);
348 // Multiply the real a with b
349 v1 = pmul(v1, b.v);
350 // Multiply the imag a with b
351 v2 = pmul(v2, b.v);
352 // Conjugate v2
353 v2 = Packet1cd(v2).conjugate().v;
354 // Swap real/imag elements in v2.
355 v2 = (Packet2d)__builtin_msa_shf_w((v4i32)v2, EIGEN_MSA_SHF_I8(2, 3, 0, 1));
356 // Add and return the result
357 v = padd(v1, v2);
358 return *this;
359 }
360 EIGEN_STRONG_INLINE Packet1cd operator*(const Packet1cd& b) const { return Packet1cd(*this) *= b; }
361 EIGEN_STRONG_INLINE Packet1cd& operator+=(const Packet1cd& b) {
362 v = padd(v, b.v);
363 return *this;
364 }
365 EIGEN_STRONG_INLINE Packet1cd operator+(const Packet1cd& b) const { return Packet1cd(*this) += b; }
366 EIGEN_STRONG_INLINE Packet1cd& operator-=(const Packet1cd& b) {
367 v = psub(v, b.v);
368 return *this;
369 }
370 EIGEN_STRONG_INLINE Packet1cd operator-(const Packet1cd& b) const { return Packet1cd(*this) -= b; }
371 EIGEN_STRONG_INLINE Packet1cd& operator/=(const Packet1cd& b) {
372 *this *= b.conjugate();
373 Packet2d s = pmul<Packet2d>(b.v, b.v);
374 s = padd(s, preverse<Packet2d>(s));
375 v = pdiv(v, s);
376 return *this;
377 }
378 EIGEN_STRONG_INLINE Packet1cd operator/(const Packet1cd& b) const { return Packet1cd(*this) /= b; }
379 EIGEN_STRONG_INLINE Packet1cd operator-(void) const { return Packet1cd(pnegate(v)); }
380
381 Packet2d v;
382};
383
384inline std::ostream& operator<<(std::ostream& os, const Packet1cd& value) {
385 os << "[ (" << value.v[0] << ", " << value.v[1] << "i) ]";
386 return os;
387}
388
389template <>
390struct packet_traits<std::complex<double> > : default_packet_traits {
391 typedef Packet1cd type;
392 typedef Packet1cd half;
393 enum {
394 Vectorizable = 1,
395 AlignedOnScalar = 0,
396 size = 1,
397
398 HasAdd = 1,
399 HasSub = 1,
400 HasMul = 1,
401 HasDiv = 1,
402 HasNegate = 1,
403 HasSqrt = 1,
404 HasLog = 1,
405 HasAbs = 0,
406 HasAbs2 = 0,
407 HasMin = 0,
408 HasMax = 0,
409 HasSetLinear = 0
410 };
411};
412
413template <>
414struct unpacket_traits<Packet1cd> {
415 typedef std::complex<double> type;
416 typedef Packet2d as_real;
417 enum {
418 size = 1,
419 alignment = Aligned16,
420 vectorizable = true,
421 masked_load_available = false,
422 masked_store_available = false
423 };
424 typedef Packet1cd half;
425};
426
427template <>
428EIGEN_STRONG_INLINE Packet1cd pload<Packet1cd>(const std::complex<double>* from) {
429 EIGEN_MSA_DEBUG;
430
431 EIGEN_DEBUG_ALIGNED_LOAD return Packet1cd(pload<Packet2d>((const double*)from));
432}
433
434template <>
435EIGEN_STRONG_INLINE Packet1cd ploadu<Packet1cd>(const std::complex<double>* from) {
436 EIGEN_MSA_DEBUG;
437
438 EIGEN_DEBUG_UNALIGNED_LOAD return Packet1cd(ploadu<Packet2d>((const double*)from));
439}
440
441template <>
442EIGEN_STRONG_INLINE Packet1cd pset1<Packet1cd>(const std::complex<double>& from) {
443 EIGEN_MSA_DEBUG;
444
445 return Packet1cd(from);
446}
447
448template <>
449EIGEN_STRONG_INLINE Packet1cd padd<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
450 EIGEN_MSA_DEBUG;
451
452 return a + b;
453}
454
455template <>
456EIGEN_STRONG_INLINE Packet1cd psub<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
457 EIGEN_MSA_DEBUG;
458
459 return a - b;
460}
461
462template <>
463EIGEN_STRONG_INLINE Packet1cd pnegate(const Packet1cd& a) {
464 EIGEN_MSA_DEBUG;
465
466 return -a;
467}
468
469template <>
470EIGEN_STRONG_INLINE Packet1cd pconj(const Packet1cd& a) {
471 EIGEN_MSA_DEBUG;
472
473 return a.conjugate();
474}
475
476template <>
477EIGEN_STRONG_INLINE Packet1cd pmul<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
478 EIGEN_MSA_DEBUG;
479
480 return a * b;
481}
482
483template <>
484EIGEN_STRONG_INLINE Packet1cd pand<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
485 EIGEN_MSA_DEBUG;
486
487 return Packet1cd(pand(a.v, b.v));
488}
489
490template <>
491EIGEN_STRONG_INLINE Packet1cd por<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
492 EIGEN_MSA_DEBUG;
493
494 return Packet1cd(por(a.v, b.v));
495}
496
497template <>
498EIGEN_STRONG_INLINE Packet1cd pxor<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
499 EIGEN_MSA_DEBUG;
500
501 return Packet1cd(pxor(a.v, b.v));
502}
503
504template <>
505EIGEN_STRONG_INLINE Packet1cd pandnot<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
506 EIGEN_MSA_DEBUG;
507
508 return Packet1cd(pandnot(a.v, b.v));
509}
510
511template <>
512EIGEN_STRONG_INLINE Packet1cd ploaddup<Packet1cd>(const std::complex<double>* from) {
513 EIGEN_MSA_DEBUG;
514
515 return pset1<Packet1cd>(*from);
516}
517
518template <>
519EIGEN_STRONG_INLINE void pstore<std::complex<double> >(std::complex<double>* to, const Packet1cd& from) {
520 EIGEN_MSA_DEBUG;
521
522 EIGEN_DEBUG_ALIGNED_STORE pstore<double>((double*)to, from.v);
523}
524
525template <>
526EIGEN_STRONG_INLINE void pstoreu<std::complex<double> >(std::complex<double>* to, const Packet1cd& from) {
527 EIGEN_MSA_DEBUG;
528
529 EIGEN_DEBUG_UNALIGNED_STORE pstoreu<double>((double*)to, from.v);
530}
531
532template <>
533EIGEN_STRONG_INLINE void prefetch<std::complex<double> >(const std::complex<double>* addr) {
534 EIGEN_MSA_DEBUG;
535
536 prefetch(reinterpret_cast<const double*>(addr));
537}
538
539template <>
540EIGEN_DEVICE_FUNC inline Packet1cd pgather<std::complex<double>, Packet1cd>(const std::complex<double>* from,
541 Index stride __attribute__((unused))) {
542 EIGEN_MSA_DEBUG;
543
544 Packet1cd res;
545 res.v[0] = std::real(from[0]);
546 res.v[1] = std::imag(from[0]);
547 return res;
548}
549
550template <>
551EIGEN_DEVICE_FUNC inline void pscatter<std::complex<double>, Packet1cd>(std::complex<double>* to, const Packet1cd& from,
552 Index stride __attribute__((unused))) {
553 EIGEN_MSA_DEBUG;
554
555 pstore(to, from);
556}
557
558template <>
559EIGEN_STRONG_INLINE std::complex<double> pfirst<Packet1cd>(const Packet1cd& a) {
560 EIGEN_MSA_DEBUG;
561
562 return std::complex<double>(a.v[0], a.v[1]);
563}
564
565EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet1cd, Packet2d)
566
567template <>
568EIGEN_STRONG_INLINE Packet1cd pdiv<Packet1cd>(const Packet1cd& a, const Packet1cd& b) {
569 EIGEN_MSA_DEBUG;
570
571 return a / b;
572}
573
574EIGEN_STRONG_INLINE Packet1cd pcplxflip /*<Packet1cd>*/ (const Packet1cd& x) {
575 EIGEN_MSA_DEBUG;
576
577 return Packet1cd(preverse(Packet2d(x.v)));
578}
579
580inline std::ostream& operator<<(std::ostream& os, const PacketBlock<Packet1cd, 2>& value) {
581 os << "[ " << value.packet[0] << ", " << std::endl << " " << value.packet[1] << " ]";
582 return os;
583}
584
585EIGEN_STRONG_INLINE void ptranspose(PacketBlock<Packet1cd, 2>& kernel) {
586 EIGEN_MSA_DEBUG;
587
588 Packet2d v1, v2;
589
590 v1 = (Packet2d)__builtin_msa_ilvev_d((v2i64)kernel.packet[0].v, (v2i64)kernel.packet[1].v);
591 // Get the imag values of a
592 v2 = (Packet2d)__builtin_msa_ilvod_d((v2i64)kernel.packet[0].v, (v2i64)kernel.packet[1].v);
593
594 kernel.packet[0].v = v1;
595 kernel.packet[1].v = v2;
596}
597
598EIGEN_INSTANTIATE_COMPLEX_MATH_FUNCS(Packet2cf)
599EIGEN_INSTANTIATE_COMPLEX_MATH_FUNCS_NO_EXP(Packet1cd)
600
601EIGEN_GCC_FAST_MATH_COMPLEX_VECTORIZE_WORKAROUND_POP
602
603} // end namespace internal
604
605} // end namespace Eigen
606
607#endif // EIGEN_COMPLEX_MSA_H
Eigen::Aligned16
@ Aligned16
Definition Constants.h:238