1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. Eigen itself is part of the KDE project. 3 // 4 // Copyright (C) 2006-2008 Benoit Jacob <jacob.benoit.1 (at) gmail.com> 5 // 6 // This Source Code Form is subject to the terms of the Mozilla 7 // Public License v. 2.0. If a copy of the MPL was not distributed 8 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 9 10 #include "main.h" 11 12 // using namespace Eigen; 13 14 template<typename Scalar> bool areApprox(const Scalar* a, const Scalar* b, int size) 15 { 16 for (int i=0; i<size; ++i) 17 if (!ei_isApprox(a[i],b[i])) return false; 18 return true; 19 } 20 21 #define CHECK_CWISE(REFOP, POP) { \ 22 for (int i=0; i<PacketSize; ++i) \ 23 ref[i] = REFOP(data1[i], data1[i+PacketSize]); \ 24 ei_pstore(data2, POP(ei_pload(data1), ei_pload(data1+PacketSize))); \ 25 VERIFY(areApprox(ref, data2, PacketSize) && #POP); \ 26 } 27 28 #define REF_ADD(a,b) ((a)+(b)) 29 #define REF_SUB(a,b) ((a)-(b)) 30 #define REF_MUL(a,b) ((a)*(b)) 31 #define REF_DIV(a,b) ((a)/(b)) 32 33 namespace std { 34 35 template<> const complex<float>& min(const complex<float>& a, const complex<float>& b) 36 { return a.real() < b.real() ? a : b; } 37 38 template<> const complex<float>& max(const complex<float>& a, const complex<float>& b) 39 { return a.real() < b.real() ? b : a; } 40 41 } 42 43 template<typename Scalar> void packetmath() 44 { 45 typedef typename ei_packet_traits<Scalar>::type Packet; 46 const int PacketSize = ei_packet_traits<Scalar>::size; 47 48 const int size = PacketSize*4; 49 EIGEN_ALIGN_128 Scalar data1[ei_packet_traits<Scalar>::size*4]; 50 EIGEN_ALIGN_128 Scalar data2[ei_packet_traits<Scalar>::size*4]; 51 EIGEN_ALIGN_128 Packet packets[PacketSize*2]; 52 EIGEN_ALIGN_128 Scalar ref[ei_packet_traits<Scalar>::size*4]; 53 for (int i=0; i<size; ++i) 54 { 55 data1[i] = ei_random<Scalar>(); 56 data2[i] = ei_random<Scalar>(); 57 } 58 59 ei_pstore(data2, ei_pload(data1)); 60 VERIFY(areApprox(data1, data2, PacketSize) && "aligned load/store"); 61 62 for (int offset=0; offset<PacketSize; ++offset) 63 { 64 ei_pstore(data2, ei_ploadu(data1+offset)); 65 VERIFY(areApprox(data1+offset, data2, PacketSize) && "ei_ploadu"); 66 } 67 68 for (int offset=0; offset<PacketSize; ++offset) 69 { 70 ei_pstoreu(data2+offset, ei_pload(data1)); 71 VERIFY(areApprox(data1, data2+offset, PacketSize) && "ei_pstoreu"); 72 } 73 74 for (int offset=0; offset<PacketSize; ++offset) 75 { 76 packets[0] = ei_pload(data1); 77 packets[1] = ei_pload(data1+PacketSize); 78 if (offset==0) ei_palign<0>(packets[0], packets[1]); 79 else if (offset==1) ei_palign<1>(packets[0], packets[1]); 80 else if (offset==2) ei_palign<2>(packets[0], packets[1]); 81 else if (offset==3) ei_palign<3>(packets[0], packets[1]); 82 ei_pstore(data2, packets[0]); 83 84 for (int i=0; i<PacketSize; ++i) 85 ref[i] = data1[i+offset]; 86 87 typedef Matrix<Scalar, PacketSize, 1> Vector; 88 VERIFY(areApprox(ref, data2, PacketSize) && "ei_palign"); 89 } 90 91 CHECK_CWISE(REF_ADD, ei_padd); 92 CHECK_CWISE(REF_SUB, ei_psub); 93 CHECK_CWISE(REF_MUL, ei_pmul); 94 #ifndef EIGEN_VECTORIZE_ALTIVEC 95 if (!ei_is_same_type<Scalar,int>::ret) 96 CHECK_CWISE(REF_DIV, ei_pdiv); 97 #endif 98 CHECK_CWISE(std::min, ei_pmin); 99 CHECK_CWISE(std::max, ei_pmax); 100 101 for (int i=0; i<PacketSize; ++i) 102 ref[i] = data1[0]; 103 ei_pstore(data2, ei_pset1(data1[0])); 104 VERIFY(areApprox(ref, data2, PacketSize) && "ei_pset1"); 105 106 VERIFY(ei_isApprox(data1[0], ei_pfirst(ei_pload(data1))) && "ei_pfirst"); 107 108 ref[0] = 0; 109 for (int i=0; i<PacketSize; ++i) 110 ref[0] += data1[i]; 111 VERIFY(ei_isApprox(ref[0], ei_predux(ei_pload(data1))) && "ei_predux"); 112 113 for (int j=0; j<PacketSize; ++j) 114 { 115 ref[j] = 0; 116 for (int i=0; i<PacketSize; ++i) 117 ref[j] += data1[i+j*PacketSize]; 118 packets[j] = ei_pload(data1+j*PacketSize); 119 } 120 ei_pstore(data2, ei_preduxp(packets)); 121 VERIFY(areApprox(ref, data2, PacketSize) && "ei_preduxp"); 122 } 123 124 void test_eigen2_packetmath() 125 { 126 for(int i = 0; i < g_repeat; i++) { 127 CALL_SUBTEST_1( packetmath<float>() ); 128 CALL_SUBTEST_2( packetmath<double>() ); 129 CALL_SUBTEST_3( packetmath<int>() ); 130 CALL_SUBTEST_4( packetmath<std::complex<float> >() ); 131 } 132 } 133
