cart-elc

Source code for CART-ELC
git clone git://git.laack.co/cart-elc.git
Log | Files | Refs | README | LICENSE

blasutil.cpp (6413B)


      1 // This file is part of Eigen, a lightweight C++ template library
      2 // for linear algebra.
      3 //
      4 // Copyright (C) 2020 Everton Constantino <everton.constantino@ibm.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 // Disable "ignoring attributes on template argument"
     13 // for packet_traits<Packet*>
     14 // => The only workaround would be to wrap _m128 and the likes
     15 //    within wrappers.
     16 #if EIGEN_GNUC_AT_LEAST(6,0)
     17     #pragma GCC diagnostic ignored "-Wignored-attributes"
     18 #endif
     19 
     20 #define GET(i,j) (StorageOrder == RowMajor ? (i)*stride + (j) : (i) + (j)*stride)
     21 #define SCATTER(i,j,k) (StorageOrder == RowMajor ? ((i)+(k))*stride + (j) : (i) + ((j)+(k))*stride)
     22 
     23 template<typename Scalar, typename Packet>
     24 void compare(const Packet& a, const Packet& b)
     25 {
     26     int pktsz = internal::packet_traits<Scalar>::size;
     27     Scalar *buffA = new Scalar[pktsz];
     28     Scalar *buffB = new Scalar[pktsz];
     29 
     30     internal::pstoreu<Scalar, Packet>(buffA, a);
     31     internal::pstoreu<Scalar, Packet>(buffB, b);
     32 
     33     for(int i = 0; i < pktsz; i++)
     34     {
     35         VERIFY_IS_EQUAL(buffA[i], buffB[i]);
     36     }
     37 
     38     delete[] buffA;
     39     delete[] buffB;
     40 }
     41 
     42 template<typename Scalar, int StorageOrder, int n>
     43 struct PacketBlockSet
     44 {
     45     typedef typename internal::packet_traits<Scalar>::type Packet;
     46 
     47     void setPacketBlock(internal::PacketBlock<Packet,n>& block, Scalar value)
     48     {
     49         for(int idx = 0; idx < n; idx++)
     50         {
     51             block.packet[idx] = internal::pset1<Packet>(value);
     52         }
     53     }
     54 
     55     void comparePacketBlock(Scalar *data, int i, int j, int stride, internal::PacketBlock<Packet, n>& block)
     56     {
     57         for(int idx = 0; idx < n; idx++)
     58         {
     59             Packet line = internal::ploadu<Packet>(data + SCATTER(i,j,idx));
     60             compare<Scalar, Packet>(block.packet[idx], line);
     61         }
     62     }
     63 };
     64 
     65 template<typename Scalar, int StorageOrder, int BlockSize>
     66 void run_bdmp_spec_1()
     67 {
     68     typedef internal::blas_data_mapper<Scalar, int, StorageOrder> BlasDataMapper;
     69     int packetSize = internal::packet_traits<Scalar>::size;
     70     int minSize = std::max<int>(packetSize, BlockSize);
     71     typedef typename internal::packet_traits<Scalar>::type Packet;
     72 
     73     int szm = internal::random<int>(minSize,500), szn = internal::random<int>(minSize,500);
     74     int stride = StorageOrder == RowMajor ? szn : szm;
     75     Scalar *d = new Scalar[szn*szm];
     76 
     77     // Initializing with random entries
     78     for(int i = 0; i < szm*szn; i++)
     79     {
     80         d[i] = internal::random<Scalar>(static_cast<Scalar>(3), static_cast<Scalar>(10));
     81     }
     82 
     83     BlasDataMapper bdm(d, stride);
     84 
     85     // Testing operator()
     86     for(int i = 0; i < szm; i++)
     87     {
     88         for(int j = 0; j < szn; j++)
     89         {
     90             VERIFY_IS_EQUAL(d[GET(i,j)], bdm(i,j));
     91         }
     92     }
     93 
     94     // Testing getSubMapper and getLinearMapper
     95     int i0 = internal::random<int>(0,szm-2);
     96     int j0 = internal::random<int>(0,szn-2);
     97     for(int i = i0; i < szm; i++)
     98     {
     99         for(int j = j0; j < szn; j++)
    100         {
    101             const BlasDataMapper& bdmSM = bdm.getSubMapper(i0,j0);
    102             const internal::BlasLinearMapper<Scalar, int, 0>& bdmLM = bdm.getLinearMapper(i0,j0);
    103 
    104             Scalar v = bdmSM(i - i0, j - j0);
    105             Scalar vd = d[GET(i,j)];
    106             VERIFY_IS_EQUAL(vd, v);
    107             VERIFY_IS_EQUAL(vd, bdmLM(GET(i-i0, j-j0)));
    108         }
    109     }
    110 
    111     // Testing loadPacket
    112     for(int i = 0; i < szm - minSize; i++)
    113     {
    114         for(int j = 0; j < szn - minSize; j++)
    115         {
    116             Packet pktBDM = bdm.template loadPacket<Packet>(i,j);
    117             Packet pktD = internal::ploadu<Packet>(d + GET(i,j));
    118 
    119             compare<Scalar, Packet>(pktBDM, pktD);
    120         }
    121     }
    122 
    123     // Testing gatherPacket
    124     Scalar *buff = new Scalar[packetSize];
    125     for(int i = 0; i < szm - minSize; i++)
    126     {
    127         for(int j = 0; j < szn - minSize; j++)
    128         {
    129             Packet p = bdm.template gatherPacket<Packet>(i,j);
    130             internal::pstoreu<Scalar, Packet>(buff, p);
    131 
    132             for(int k = 0; k < packetSize; k++)
    133             {
    134                 VERIFY_IS_EQUAL(d[SCATTER(i,j,k)], buff[k]);
    135             }
    136 
    137         }
    138     }
    139     delete[] buff;
    140 
    141     // Testing scatterPacket
    142     for(int i = 0; i < szm - minSize; i++)
    143     {
    144         for(int j = 0; j < szn - minSize; j++)
    145         {
    146             Packet p = internal::pset1<Packet>(static_cast<Scalar>(1));
    147             bdm.template scatterPacket<Packet>(i,j,p);
    148             for(int k = 0; k < packetSize; k++)
    149             {
    150                 VERIFY_IS_EQUAL(d[SCATTER(i,j,k)], static_cast<Scalar>(1));
    151             }
    152         }
    153     }
    154 
    155     //Testing storePacketBlock
    156     internal::PacketBlock<Packet, BlockSize> block;
    157 
    158     PacketBlockSet<Scalar, StorageOrder, BlockSize> pbs;
    159     pbs.setPacketBlock(block, static_cast<Scalar>(2));
    160 
    161     for(int i = 0; i < szm - minSize; i++)
    162     {
    163         for(int j = 0; j < szn - minSize; j++)
    164         {
    165             bdm.template storePacketBlock<Packet, BlockSize>(i, j, block);
    166 
    167             pbs.comparePacketBlock(d, i, j, stride, block);
    168         }
    169     }
    170 
    171     delete[] d;
    172 }
    173 
    174 template<typename Scalar>
    175 void run_test()
    176 {
    177     run_bdmp_spec_1<Scalar, RowMajor, 1>();
    178     run_bdmp_spec_1<Scalar, ColMajor, 1>();
    179     run_bdmp_spec_1<Scalar, RowMajor, 2>();
    180     run_bdmp_spec_1<Scalar, ColMajor, 2>();
    181     run_bdmp_spec_1<Scalar, RowMajor, 4>();
    182     run_bdmp_spec_1<Scalar, ColMajor, 4>();
    183     run_bdmp_spec_1<Scalar, RowMajor, 8>();
    184     run_bdmp_spec_1<Scalar, ColMajor, 8>();
    185     run_bdmp_spec_1<Scalar, RowMajor, 16>();
    186     run_bdmp_spec_1<Scalar, ColMajor, 16>();
    187 }
    188 
    189 EIGEN_DECLARE_TEST(blasutil)
    190 {
    191     for(int i = 0; i < g_repeat; i++)
    192     {
    193         CALL_SUBTEST_1(run_test<numext::int8_t>());
    194         CALL_SUBTEST_2(run_test<numext::int16_t>());
    195         CALL_SUBTEST_3(run_test<numext::int32_t>());
    196 
    197 // TODO: Replace this by a call to numext::int64_t as soon as we have a way to
    198 // detect the typedef for int64_t on all platforms
    199 #if EIGEN_HAS_CXX11
    200         CALL_SUBTEST_4(run_test<signed long long>());
    201 #else
    202         CALL_SUBTEST_4(run_test<signed long>());
    203 #endif
    204 
    205         CALL_SUBTEST_5(run_test<float_t>());
    206         CALL_SUBTEST_6(run_test<double_t>());
    207         CALL_SUBTEST_7(run_test<std::complex<float> >());
    208         CALL_SUBTEST_8(run_test<std::complex<double> >());
    209     }
    210 }