1 // This file is part of Eigen, a lightweight C++ template library
2 // for linear algebra.
3 //
4 // Copyright (C) 2020 Everton Constantino <[email protected]>
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>
compare(const Packet & a,const Packet & b)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
setPacketBlockPacketBlockSet47 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
comparePacketBlockPacketBlockSet55 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>
run_bdmp_spec_1()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>
run_test()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
EIGEN_DECLARE_TEST(blasutil)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 }
211