1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2010 Gael Guennebaud <[email protected]> 5 // Copyright (C) 2016 Konstantinos Margaritis <[email protected]> 6 // 7 // This Source Code Form is subject to the terms of the Mozilla 8 // Public License v. 2.0. If a copy of the MPL was not distributed 9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10 11 #ifndef EIGEN_COMPLEX32_ALTIVEC_H 12 #define EIGEN_COMPLEX32_ALTIVEC_H 13 14 namespace Eigen { 15 16 namespace internal { 17 18 #if !defined(__ARCH__) || (defined(__ARCH__) && __ARCH__ >= 12) 19 static Packet4ui p4ui_CONJ_XOR = { 0x00000000, 0x80000000, 0x00000000, 0x80000000 }; //vec_mergeh((Packet4ui)p4i_ZERO, (Packet4ui)p4f_MZERO); 20 #endif 21 22 static Packet2ul p2ul_CONJ_XOR1 = (Packet2ul) vec_sld((Packet4ui) p2d_ZERO_, (Packet4ui) p2l_ZERO, 8);//{ 0x8000000000000000, 0x0000000000000000 }; 23 static Packet2ul p2ul_CONJ_XOR2 = (Packet2ul) vec_sld((Packet4ui) p2l_ZERO, (Packet4ui) p2d_ZERO_, 8);//{ 0x8000000000000000, 0x0000000000000000 }; 24 25 struct Packet1cd 26 { Packet1cdPacket1cd27 EIGEN_STRONG_INLINE Packet1cd() {} Packet1cdPacket1cd28 EIGEN_STRONG_INLINE explicit Packet1cd(const Packet2d& a) : v(a) {} 29 Packet2d v; 30 }; 31 32 struct Packet2cf 33 { Packet2cfPacket2cf34 EIGEN_STRONG_INLINE Packet2cf() {} Packet2cfPacket2cf35 EIGEN_STRONG_INLINE explicit Packet2cf(const Packet4f& a) : v(a) {} 36 #if !defined(__ARCH__) || (defined(__ARCH__) && __ARCH__ < 12) 37 union { 38 Packet4f v; 39 Packet1cd cd[2]; 40 }; 41 #else 42 Packet4f v; 43 #endif 44 }; 45 46 template<> struct packet_traits<std::complex<float> > : default_packet_traits 47 { 48 typedef Packet2cf type; 49 typedef Packet2cf half; 50 enum { 51 Vectorizable = 1, 52 AlignedOnScalar = 1, 53 size = 2, 54 HasHalfPacket = 0, 55 56 HasAdd = 1, 57 HasSub = 1, 58 HasMul = 1, 59 HasDiv = 1, 60 HasNegate = 1, 61 HasAbs = 0, 62 HasAbs2 = 0, 63 HasMin = 0, 64 HasMax = 0, 65 HasBlend = 1, 66 HasSetLinear = 0 67 }; 68 }; 69 70 71 template<> struct packet_traits<std::complex<double> > : default_packet_traits 72 { 73 typedef Packet1cd type; 74 typedef Packet1cd half; 75 enum { 76 Vectorizable = 1, 77 AlignedOnScalar = 1, 78 size = 1, 79 HasHalfPacket = 0, 80 81 HasAdd = 1, 82 HasSub = 1, 83 HasMul = 1, 84 HasDiv = 1, 85 HasNegate = 1, 86 HasAbs = 0, 87 HasAbs2 = 0, 88 HasMin = 0, 89 HasMax = 0, 90 HasSetLinear = 0 91 }; 92 }; 93 94 template<> struct unpacket_traits<Packet2cf> { typedef std::complex<float> type; enum {size=2, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet2cf half; }; 95 template<> struct unpacket_traits<Packet1cd> { typedef std::complex<double> type; enum {size=1, alignment=Aligned16, vectorizable=true, masked_load_available=false, masked_store_available=false}; typedef Packet1cd half; }; 96 97 /* Forward declaration */ 98 EIGEN_STRONG_INLINE void ptranspose(PacketBlock<Packet2cf,2>& kernel); 99 100 /* complex<double> first */ 101 template<> EIGEN_STRONG_INLINE Packet1cd pload <Packet1cd>(const std::complex<double>* from) { EIGEN_DEBUG_ALIGNED_LOAD return Packet1cd(pload<Packet2d>((const double*)from)); } 102 template<> EIGEN_STRONG_INLINE Packet1cd ploadu<Packet1cd>(const std::complex<double>* from) { EIGEN_DEBUG_UNALIGNED_LOAD return Packet1cd(ploadu<Packet2d>((const double*)from)); } 103 template<> EIGEN_STRONG_INLINE void pstore <std::complex<double> >(std::complex<double> * to, const Packet1cd& from) { EIGEN_DEBUG_ALIGNED_STORE pstore((double*)to, from.v); } 104 template<> EIGEN_STRONG_INLINE void pstoreu<std::complex<double> >(std::complex<double> * to, const Packet1cd& from) { EIGEN_DEBUG_UNALIGNED_STORE pstoreu((double*)to, from.v); } 105 106 template<> EIGEN_STRONG_INLINE Packet1cd pset1<Packet1cd>(const std::complex<double>& from) 107 { /* here we really have to use unaligned loads :( */ return ploadu<Packet1cd>(&from); } 108 109 template<> EIGEN_DEVICE_FUNC inline Packet1cd pgather<std::complex<double>, Packet1cd>(const std::complex<double>* from, Index stride EIGEN_UNUSED) 110 { 111 return pload<Packet1cd>(from); 112 } 113 template<> EIGEN_DEVICE_FUNC inline void pscatter<std::complex<double>, Packet1cd>(std::complex<double>* to, const Packet1cd& from, Index stride EIGEN_UNUSED) 114 { 115 pstore<std::complex<double> >(to, from); 116 } 117 template<> EIGEN_STRONG_INLINE Packet1cd padd<Packet1cd>(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(a.v + b.v); } 118 template<> EIGEN_STRONG_INLINE Packet1cd psub<Packet1cd>(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(a.v - b.v); } 119 template<> EIGEN_STRONG_INLINE Packet1cd pnegate(const Packet1cd& a) { return Packet1cd(pnegate(Packet2d(a.v))); } 120 template<> EIGEN_STRONG_INLINE Packet1cd pconj(const Packet1cd& a) { return Packet1cd((Packet2d)vec_xor((Packet2d)a.v, (Packet2d)p2ul_CONJ_XOR2)); } 121 template<> EIGEN_STRONG_INLINE Packet1cd pmul<Packet1cd>(const Packet1cd& a, const Packet1cd& b) 122 { 123 Packet2d a_re, a_im, v1, v2; 124 125 // Permute and multiply the real parts of a and b 126 a_re = vec_perm(a.v, a.v, p16uc_PSET64_HI); 127 // Get the imaginary parts of a 128 a_im = vec_perm(a.v, a.v, p16uc_PSET64_LO); 129 // multiply a_re * b 130 v1 = vec_madd(a_re, b.v, p2d_ZERO); 131 // multiply a_im * b and get the conjugate result 132 v2 = vec_madd(a_im, b.v, p2d_ZERO); 133 v2 = (Packet2d) vec_sld((Packet4ui)v2, (Packet4ui)v2, 8); 134 v2 = (Packet2d) vec_xor((Packet2d)v2, (Packet2d) p2ul_CONJ_XOR1); 135 136 return Packet1cd(v1 + v2); 137 } 138 template<> EIGEN_STRONG_INLINE Packet1cd pand <Packet1cd>(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(vec_and(a.v,b.v)); } 139 template<> EIGEN_STRONG_INLINE Packet1cd por <Packet1cd>(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(vec_or(a.v,b.v)); } 140 template<> EIGEN_STRONG_INLINE Packet1cd pxor <Packet1cd>(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(vec_xor(a.v,b.v)); } 141 template<> EIGEN_STRONG_INLINE Packet1cd pandnot <Packet1cd>(const Packet1cd& a, const Packet1cd& b) { return Packet1cd(vec_and(a.v, vec_nor(b.v,b.v))); } 142 template<> EIGEN_STRONG_INLINE Packet1cd ploaddup<Packet1cd>(const std::complex<double>* from) { return pset1<Packet1cd>(*from); } 143 template<> EIGEN_STRONG_INLINE Packet1cd pcmp_eq(const Packet1cd& a, const Packet1cd& b) { 144 Packet2d eq = vec_cmpeq (a.v, b.v); 145 Packet2d tmp = { eq[1], eq[0] }; 146 return (Packet1cd)pand<Packet2d>(eq, tmp); 147 } 148 149 template<> EIGEN_STRONG_INLINE void prefetch<std::complex<double> >(const std::complex<double> * addr) { EIGEN_ZVECTOR_PREFETCH(addr); } 150 151 template<> EIGEN_STRONG_INLINE std::complex<double> pfirst<Packet1cd>(const Packet1cd& a) 152 { 153 std::complex<double> EIGEN_ALIGN16 res; 154 pstore<std::complex<double> >(&res, a); 155 156 return res; 157 } 158 159 template<> EIGEN_STRONG_INLINE Packet1cd preverse(const Packet1cd& a) { return a; } 160 template<> EIGEN_STRONG_INLINE std::complex<double> predux<Packet1cd>(const Packet1cd& a) 161 { 162 return pfirst(a); 163 } 164 template<> EIGEN_STRONG_INLINE std::complex<double> predux_mul<Packet1cd>(const Packet1cd& a) 165 { 166 return pfirst(a); 167 } 168 EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet1cd,Packet2d) 169 170 template<> EIGEN_STRONG_INLINE Packet1cd pdiv<Packet1cd>(const Packet1cd& a, const Packet1cd& b) 171 { 172 // TODO optimize it for AltiVec 173 Packet1cd res = pmul(a,pconj(b)); 174 Packet2d s = vec_madd(b.v, b.v, p2d_ZERO_); 175 return Packet1cd(pdiv(res.v, s + vec_perm(s, s, p16uc_REVERSE64))); 176 } 177 178 EIGEN_STRONG_INLINE Packet1cd pcplxflip/*<Packet1cd>*/(const Packet1cd& x) 179 { 180 return Packet1cd(preverse(Packet2d(x.v))); 181 } 182 183 EIGEN_STRONG_INLINE void ptranspose(PacketBlock<Packet1cd,2>& kernel) 184 { 185 Packet2d tmp = vec_perm(kernel.packet[0].v, kernel.packet[1].v, p16uc_TRANSPOSE64_HI); 186 kernel.packet[1].v = vec_perm(kernel.packet[0].v, kernel.packet[1].v, p16uc_TRANSPOSE64_LO); 187 kernel.packet[0].v = tmp; 188 } 189 190 /* complex<float> follows */ 191 template<> EIGEN_STRONG_INLINE Packet2cf pload <Packet2cf>(const std::complex<float>* from) { EIGEN_DEBUG_ALIGNED_LOAD return Packet2cf(pload<Packet4f>((const float*)from)); } 192 template<> EIGEN_STRONG_INLINE Packet2cf ploadu<Packet2cf>(const std::complex<float>* from) { EIGEN_DEBUG_UNALIGNED_LOAD return Packet2cf(ploadu<Packet4f>((const float*)from)); } 193 template<> EIGEN_STRONG_INLINE void pstore <std::complex<float> >(std::complex<float> * to, const Packet2cf& from) { EIGEN_DEBUG_ALIGNED_STORE pstore((float*)to, from.v); } 194 template<> EIGEN_STRONG_INLINE void pstoreu<std::complex<float> >(std::complex<float> * to, const Packet2cf& from) { EIGEN_DEBUG_UNALIGNED_STORE pstoreu((float*)to, from.v); } 195 196 template<> EIGEN_STRONG_INLINE std::complex<float> pfirst<Packet2cf>(const Packet2cf& a) 197 { 198 std::complex<float> EIGEN_ALIGN16 res[2]; 199 pstore<std::complex<float> >(res, a); 200 201 return res[0]; 202 } 203 204 205 #if !defined(__ARCH__) || (defined(__ARCH__) && __ARCH__ < 12) 206 template<> EIGEN_STRONG_INLINE Packet2cf pset1<Packet2cf>(const std::complex<float>& from) 207 { 208 Packet2cf res; 209 res.cd[0] = Packet1cd(vec_ld2f((const float *)&from)); 210 res.cd[1] = res.cd[0]; 211 return res; 212 } 213 #else 214 template<> EIGEN_STRONG_INLINE Packet2cf pset1<Packet2cf>(const std::complex<float>& from) 215 { 216 Packet2cf res; 217 if((std::ptrdiff_t(&from) % 16) == 0) 218 res.v = pload<Packet4f>((const float *)&from); 219 else 220 res.v = ploadu<Packet4f>((const float *)&from); 221 res.v = vec_perm(res.v, res.v, p16uc_PSET64_HI); 222 return res; 223 } 224 #endif 225 226 template<> EIGEN_DEVICE_FUNC inline Packet2cf pgather<std::complex<float>, Packet2cf>(const std::complex<float>* from, Index stride) 227 { 228 std::complex<float> EIGEN_ALIGN16 af[2]; 229 af[0] = from[0*stride]; 230 af[1] = from[1*stride]; 231 return pload<Packet2cf>(af); 232 } 233 template<> EIGEN_DEVICE_FUNC inline void pscatter<std::complex<float>, Packet2cf>(std::complex<float>* to, const Packet2cf& from, Index stride) 234 { 235 std::complex<float> EIGEN_ALIGN16 af[2]; 236 pstore<std::complex<float> >((std::complex<float> *) af, from); 237 to[0*stride] = af[0]; 238 to[1*stride] = af[1]; 239 } 240 241 template<> EIGEN_STRONG_INLINE Packet2cf padd<Packet2cf>(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(padd<Packet4f>(a.v, b.v)); } 242 template<> EIGEN_STRONG_INLINE Packet2cf psub<Packet2cf>(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(psub<Packet4f>(a.v, b.v)); } 243 template<> EIGEN_STRONG_INLINE Packet2cf pnegate(const Packet2cf& a) { return Packet2cf(pnegate(Packet4f(a.v))); } 244 245 template<> EIGEN_STRONG_INLINE Packet2cf pand <Packet2cf>(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(pand<Packet4f>(a.v,b.v)); } 246 template<> EIGEN_STRONG_INLINE Packet2cf por <Packet2cf>(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(por<Packet4f>(a.v,b.v)); } 247 template<> EIGEN_STRONG_INLINE Packet2cf pxor <Packet2cf>(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(pxor<Packet4f>(a.v,b.v)); } 248 template<> EIGEN_STRONG_INLINE Packet2cf pandnot<Packet2cf>(const Packet2cf& a, const Packet2cf& b) { return Packet2cf(pandnot<Packet4f>(a.v,b.v)); } 249 250 template<> EIGEN_STRONG_INLINE Packet2cf ploaddup<Packet2cf>(const std::complex<float>* from) { return pset1<Packet2cf>(*from); } 251 252 template<> EIGEN_STRONG_INLINE void prefetch<std::complex<float> >(const std::complex<float> * addr) { EIGEN_ZVECTOR_PREFETCH(addr); } 253 254 255 #if !defined(__ARCH__) || (defined(__ARCH__) && __ARCH__ < 12) 256 257 template<> EIGEN_STRONG_INLINE Packet2cf pcmp_eq(const Packet2cf& a, const Packet2cf& b) { 258 Packet4f eq = pcmp_eq<Packet4f> (a.v, b.v); 259 Packet2cf res; 260 Packet2d tmp1 = { eq.v4f[0][1], eq.v4f[0][0] }; 261 Packet2d tmp2 = { eq.v4f[1][1], eq.v4f[1][0] }; 262 res.v.v4f[0] = pand<Packet2d>(eq.v4f[0], tmp1); 263 res.v.v4f[1] = pand<Packet2d>(eq.v4f[1], tmp2); 264 return res; 265 } 266 267 template<> EIGEN_STRONG_INLINE Packet2cf pconj(const Packet2cf& a) 268 { 269 Packet2cf res; 270 res.v.v4f[0] = pconj(Packet1cd(reinterpret_cast<Packet2d>(a.v.v4f[0]))).v; 271 res.v.v4f[1] = pconj(Packet1cd(reinterpret_cast<Packet2d>(a.v.v4f[1]))).v; 272 return res; 273 } 274 275 template<> EIGEN_STRONG_INLINE Packet2cf pmul<Packet2cf>(const Packet2cf& a, const Packet2cf& b) 276 { 277 Packet2cf res; 278 res.v.v4f[0] = pmul(Packet1cd(reinterpret_cast<Packet2d>(a.v.v4f[0])), Packet1cd(reinterpret_cast<Packet2d>(b.v.v4f[0]))).v; 279 res.v.v4f[1] = pmul(Packet1cd(reinterpret_cast<Packet2d>(a.v.v4f[1])), Packet1cd(reinterpret_cast<Packet2d>(b.v.v4f[1]))).v; 280 return res; 281 } 282 283 template<> EIGEN_STRONG_INLINE Packet2cf preverse(const Packet2cf& a) 284 { 285 Packet2cf res; 286 res.cd[0] = a.cd[1]; 287 res.cd[1] = a.cd[0]; 288 return res; 289 } 290 291 template<> EIGEN_STRONG_INLINE std::complex<float> predux<Packet2cf>(const Packet2cf& a) 292 { 293 std::complex<float> res; 294 Packet1cd b = padd<Packet1cd>(a.cd[0], a.cd[1]); 295 vec_st2f(b.v, (float*)&res); 296 return res; 297 } 298 299 template<> EIGEN_STRONG_INLINE std::complex<float> predux_mul<Packet2cf>(const Packet2cf& a) 300 { 301 std::complex<float> res; 302 Packet1cd b = pmul<Packet1cd>(a.cd[0], a.cd[1]); 303 vec_st2f(b.v, (float*)&res); 304 return res; 305 } 306 307 EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet2cf,Packet4f) 308 309 template<> EIGEN_STRONG_INLINE Packet2cf pdiv<Packet2cf>(const Packet2cf& a, const Packet2cf& b) 310 { 311 // TODO optimize it for AltiVec 312 Packet2cf res; 313 res.cd[0] = pdiv<Packet1cd>(a.cd[0], b.cd[0]); 314 res.cd[1] = pdiv<Packet1cd>(a.cd[1], b.cd[1]); 315 return res; 316 } 317 318 EIGEN_STRONG_INLINE Packet2cf pcplxflip/*<Packet2cf>*/(const Packet2cf& x) 319 { 320 Packet2cf res; 321 res.cd[0] = pcplxflip(x.cd[0]); 322 res.cd[1] = pcplxflip(x.cd[1]); 323 return res; 324 } 325 326 EIGEN_STRONG_INLINE void ptranspose(PacketBlock<Packet2cf,2>& kernel) 327 { 328 Packet1cd tmp = kernel.packet[0].cd[1]; 329 kernel.packet[0].cd[1] = kernel.packet[1].cd[0]; 330 kernel.packet[1].cd[0] = tmp; 331 } 332 333 template<> EIGEN_STRONG_INLINE Packet2cf pblend(const Selector<2>& ifPacket, const Packet2cf& thenPacket, const Packet2cf& elsePacket) { 334 Packet2cf result; 335 const Selector<4> ifPacket4 = { ifPacket.select[0], ifPacket.select[0], ifPacket.select[1], ifPacket.select[1] }; 336 result.v = pblend<Packet4f>(ifPacket4, thenPacket.v, elsePacket.v); 337 return result; 338 } 339 #else 340 template<> EIGEN_STRONG_INLINE Packet2cf pcmp_eq(const Packet2cf& a, const Packet2cf& b) { 341 Packet4f eq = vec_cmpeq (a.v, b.v); 342 Packet4f tmp = { eq[1], eq[0], eq[3], eq[2] }; 343 return (Packet2cf)pand<Packet4f>(eq, tmp); 344 } 345 template<> EIGEN_STRONG_INLINE Packet2cf pconj(const Packet2cf& a) { return Packet2cf(pxor<Packet4f>(a.v, reinterpret_cast<Packet4f>(p4ui_CONJ_XOR))); } 346 template<> EIGEN_STRONG_INLINE Packet2cf pmul<Packet2cf>(const Packet2cf& a, const Packet2cf& b) 347 { 348 Packet4f a_re, a_im, prod, prod_im; 349 350 // Permute and multiply the real parts of a and b 351 a_re = vec_perm(a.v, a.v, p16uc_PSET32_WODD); 352 353 // Get the imaginary parts of a 354 a_im = vec_perm(a.v, a.v, p16uc_PSET32_WEVEN); 355 356 // multiply a_im * b and get the conjugate result 357 prod_im = a_im * b.v; 358 prod_im = pxor<Packet4f>(prod_im, reinterpret_cast<Packet4f>(p4ui_CONJ_XOR)); 359 // permute back to a proper order 360 prod_im = vec_perm(prod_im, prod_im, p16uc_COMPLEX32_REV); 361 362 // multiply a_re * b, add prod_im 363 prod = pmadd<Packet4f>(a_re, b.v, prod_im); 364 365 return Packet2cf(prod); 366 } 367 368 template<> EIGEN_STRONG_INLINE Packet2cf preverse(const Packet2cf& a) 369 { 370 Packet4f rev_a; 371 rev_a = vec_perm(a.v, a.v, p16uc_COMPLEX32_REV2); 372 return Packet2cf(rev_a); 373 } 374 375 template<> EIGEN_STRONG_INLINE std::complex<float> predux<Packet2cf>(const Packet2cf& a) 376 { 377 Packet4f b; 378 b = vec_sld(a.v, a.v, 8); 379 b = padd<Packet4f>(a.v, b); 380 return pfirst<Packet2cf>(Packet2cf(b)); 381 } 382 383 template<> EIGEN_STRONG_INLINE std::complex<float> predux_mul<Packet2cf>(const Packet2cf& a) 384 { 385 Packet4f b; 386 Packet2cf prod; 387 b = vec_sld(a.v, a.v, 8); 388 prod = pmul<Packet2cf>(a, Packet2cf(b)); 389 390 return pfirst<Packet2cf>(prod); 391 } 392 393 EIGEN_MAKE_CONJ_HELPER_CPLX_REAL(Packet2cf,Packet4f) 394 395 template<> EIGEN_STRONG_INLINE Packet2cf pdiv<Packet2cf>(const Packet2cf& a, const Packet2cf& b) 396 { 397 // TODO optimize it for AltiVec 398 Packet2cf res = pmul(a, pconj(b)); 399 Packet4f s = pmul<Packet4f>(b.v, b.v); 400 return Packet2cf(pdiv(res.v, padd<Packet4f>(s, vec_perm(s, s, p16uc_COMPLEX32_REV)))); 401 } 402 403 template<> EIGEN_STRONG_INLINE Packet2cf pcplxflip<Packet2cf>(const Packet2cf& x) 404 { 405 return Packet2cf(vec_perm(x.v, x.v, p16uc_COMPLEX32_REV)); 406 } 407 408 EIGEN_STRONG_INLINE void ptranspose(PacketBlock<Packet2cf,2>& kernel) 409 { 410 Packet4f tmp = vec_perm(kernel.packet[0].v, kernel.packet[1].v, p16uc_TRANSPOSE64_HI); 411 kernel.packet[1].v = vec_perm(kernel.packet[0].v, kernel.packet[1].v, p16uc_TRANSPOSE64_LO); 412 kernel.packet[0].v = tmp; 413 } 414 415 template<> EIGEN_STRONG_INLINE Packet2cf pblend(const Selector<2>& ifPacket, const Packet2cf& thenPacket, const Packet2cf& elsePacket) { 416 Packet2cf result; 417 result.v = reinterpret_cast<Packet4f>(pblend<Packet2d>(ifPacket, reinterpret_cast<Packet2d>(thenPacket.v), reinterpret_cast<Packet2d>(elsePacket.v))); 418 return result; 419 } 420 #endif 421 422 } // end namespace internal 423 424 } // end namespace Eigen 425 426 #endif // EIGEN_COMPLEX32_ALTIVEC_H 427