OpenFPM_pdata  4.1.0
Project that contain the implementation of distributed structures
 
Loading...
Searching...
No Matches
vector_performance_test.hpp
1/*
2 * vector_performance_test.hpp
3 *
4 * Created on: Jul 15, 2019
5 * Author: i-bird
6 */
7
8#ifndef VECTOR_PERFORMANCE_TEST_HPP_
9#define VECTOR_PERFORMANCE_TEST_HPP_
10
11/*
12 * vector_performance_test.hpp
13 *
14 * Created on: Jan 11, 2016
15 * Author: i-bird
16 */
17
18#ifndef OPENFPM_DATA_SRC_VECTOR_VECTOR_PERFORMANCE_TEST_HPP_
19#define OPENFPM_DATA_SRC_VECTOR_VECTOR_PERFORMANCE_TEST_HPP_
20
21#define NADD 128*128*128
22
23// Property tree
25{
26 boost::property_tree::ptree graphs;
27};
28
29report_vector_func_tests report_vector_funcs;
30
31BOOST_AUTO_TEST_SUITE( vector_performance )
32
33BOOST_AUTO_TEST_CASE(vector_performance)
34{
35 report_vector_funcs.graphs.put("performance.vector(0).funcs.nele",NADD);
36 report_vector_funcs.graphs.put("performance.vector(0).funcs.name","add");
37
38 report_vector_funcs.graphs.put("performance.vector(1).funcs.nele",NADD);
39 report_vector_funcs.graphs.put("performance.vector(1).funcs.name","get");
40
41 std::vector<double> times(N_STAT + 1);
42 std::vector<double> times_g(N_STAT + 1);
43
44 // get test
45 double tot_accu = 0.0;
46
47 for (size_t i = 0 ; i < N_STAT+1 ; i++)
48 {
49 timer t;
50 t.start();
51
52 // create a vector
54
55 // Point
57 p.setx(1.0);
58 p.sety(2.0);
59 p.setz(3.0);
60 p.sets(4.0);
61
62 p.get<P::v>()[0] = 1.0;
63 p.get<P::v>()[1] = 2.0;
64 p.get<P::v>()[2] = 7.0;
65
66 p.get<P::t>()[0][0] = 10.0;
67 p.get<P::t>()[0][1] = 13.0;
68 p.get<P::t>()[0][2] = 8.0;
69 p.get<P::t>()[1][0] = 19.0;
70 p.get<P::t>()[1][1] = 23.0;
71 p.get<P::t>()[1][2] = 5.0;
72 p.get<P::t>()[2][0] = 4.0;
73 p.get<P::t>()[2][1] = 3.0;
74 p.get<P::t>()[2][2] = 11.0;
75
76 // Add test
77
78 for (size_t i = 0 ; i < NADD ; i++)
79 {
80 v1.add(p);
81 }
82
83 t.stop();
84 times[i] = t.getwct();
85
86 timer tg;
87 tg.start();
88
89 for (size_t i = 0 ; i < NADD ; i++)
90 {
91 double accu1 = v1.template get<P::x>(i);
92 double accu2 = v1.template get<P::y>(i);
93 double accu3 = v1.template get<P::z>(i);
94 double accu4 = v1.template get<P::s>(i);
95
96 double accu5 = v1.template get<P::v>(i)[0];
97 double accu6 = v1.template get<P::v>(i)[1];
98 double accu7 = v1.template get<P::v>(i)[2];
99
100 double accu8 = v1.template get<P::t>(i)[0][0];
101 double accu9 = v1.template get<P::t>(i)[0][1];
102 double accu10 = v1.template get<P::t>(i)[0][2];
103 double accu11 = v1.template get<P::t>(i)[1][0];
104 double accu12 = v1.template get<P::t>(i)[1][1];
105 double accu13 = v1.template get<P::t>(i)[1][2];
106 double accu14 = v1.template get<P::t>(i)[2][0];
107 double accu15 = v1.template get<P::t>(i)[2][1];
108 double accu16 = v1.template get<P::t>(i)[2][2];
109
110 tot_accu += accu1 + accu2 + accu3 + accu4 + accu5 + accu6 + accu7 + accu8 + accu9 + accu10 + accu11 + accu12 +
111 accu13 + accu14 + accu15 + accu16;
112 }
113
114 tg.stop();
115
116 times_g[i] = tg.getwct();
117 }
118
119 double mean;
120 double dev;
121 standard_deviation(times,mean,dev);
122
123 report_vector_funcs.graphs.put("performance.vector(0).y.data.mean",mean);
124 report_vector_funcs.graphs.put("performance.vector(0).y.data.dev",dev);
125
126 standard_deviation(times_g,mean,dev);
127
128 report_vector_funcs.graphs.put("performance.vector(1).y.data.mean",mean);
129 report_vector_funcs.graphs.put("performance.vector(1).y.data.dev",dev);
130}
131
132template<typename vector_prop_type, typename vector_pos_type>
133__device__ __host__ void read_write(vector_prop_type & vd_prop, vector_pos_type & vd_pos, unsigned int p)
134{
135 vd_prop.template get<0>(p) = vd_pos.template get<0>(p)[0] + vd_pos.template get<0>(p)[1];
136
137 vd_prop.template get<1>(p)[0] = vd_pos.template get<0>(p)[0];
138 vd_prop.template get<1>(p)[1] = vd_pos.template get<0>(p)[1];
139
140 vd_prop.template get<2>(p)[0][0] = vd_pos.template get<0>(p)[0];
141 vd_prop.template get<2>(p)[0][1] = vd_pos.template get<0>(p)[1];
142 vd_prop.template get<2>(p)[1][0] = vd_pos.template get<0>(p)[0] +
143 vd_pos.template get<0>(p)[1];
144 vd_prop.template get<2>(p)[1][1] = vd_pos.template get<0>(p)[1] -
145 vd_pos.template get<0>(p)[0];
146
147 vd_pos.template get<0>(p)[0] += 0.01f;
148 vd_pos.template get<0>(p)[1] += 0.01f;
149}
150
151
152struct ele
153{
154 double s;
155 double v[2];
156 double t[2][2];
157};
158
159__device__ __host__ void read_write_lin(double * pos, ele * prp, unsigned int p)
160{
161 prp[p].s = pos[2*p] + pos[2*p+1];
162
163 prp[p].v[0] = pos[2*p];
164 prp[p].v[1] = pos[2*p+1];
165
166 prp[p].t[0][0] = pos[2*p];
167 prp[p].t[0][1] = pos[2*p+1];
168 prp[p].t[1][0] = pos[2*p] +
169 pos[2*p+1];
170 prp[p].t[1][1] = pos[2*p+1] -
171 pos[2*p];
172
173 pos[2*p] += 0.01f;
174 pos[2*p+1] += 0.01f;
175}
176
177__device__ __host__ void read_write_inte(double * pos, double * prp0, double * prp1, double * prp2, unsigned int p, unsigned int n_pos)
178{
179 prp0[0*n_pos + p] = pos[0*n_pos + p] + pos[1*n_pos+p];
180
181 prp1[0*n_pos + p] = pos[0*n_pos + p];
182 prp1[1*n_pos + p] = pos[1*n_pos + p];
183
184 prp2[0*n_pos*2+0*n_pos + p] = pos[0*n_pos + p];
185 prp2[0*n_pos*2+1*n_pos + p] = pos[1*n_pos + p];
186 prp2[1*n_pos*2+0*n_pos + p] = pos[0*n_pos + p] +
187 pos[1*n_pos + p];
188 prp2[1*n_pos*2+1*n_pos + p] = pos[1*n_pos + p] -
189 pos[0*n_pos + p];
190
191 pos[0*n_pos + p] += 0.01f;
192 pos[1*n_pos + p] += 0.01f;
193}
194
195BOOST_AUTO_TEST_CASE(vector_performance_layout_vs_plain_array)
196{
197 std::vector<double> times(N_STAT + 1);
198 std::vector<double> times_g(N_STAT + 1);
199
200 std::vector<double> times2(N_STAT + 1);
201 std::vector<double> times2_g(N_STAT + 1);
202
203 // get test
204 double tot_accu = 0.0;
205
206 report_vector_funcs.graphs.put("performance.vector_layout(0).funcs.nele",NADD);
207 report_vector_funcs.graphs.put("performance.vector_layout(0).funcs.name","read_write_lin");
208
209 for (size_t i = 0 ; i < N_STAT+1 ; i++)
210 {
211 // create a vector
214
215 // Point
217 p.get<0>()[0] = 1.0;
218 p.get<0>()[1] = 2.0;
219
221 pa.get<0>() = 1.0;
222
223 pa.get<1>()[0] = 1.0;
224 pa.get<1>()[1] = 1.0;
225
226 pa.get<2>()[0][0] = 1.0;
227 pa.get<2>()[0][1] = 1.0;
228 pa.get<2>()[1][0] = 1.0;
229 pa.get<2>()[1][1] = 1.0;
230
231 // Add test
232
233 for (size_t i = 0 ; i < NADD ; i++)
234 {
235 v1.add(pa);
236 v2.add(p);
237 }
238
239 timer tg;
240 tg.start();
241
242 for (size_t i = 0 ; i < NADD ; i++)
243 {
244 read_write(v1,v2,i);
245 }
246
247 tg.stop();
248
249 times_g[i] = tg.getwct();
250
251 timer tga;
252 tga.start();
253
254 double * prp = (double *)v1.getPointer<0>();
255 double * pos = (double *)v2.getPointer<0>();
256
257 for (size_t i = 0 ; i < NADD ; i++)
258 {
259 read_write_lin(pos,(struct ele *)prp,i);
260 }
261
262 tga.stop();
263
264 times[i] = tga.getwct();
265 }
266
267 double mean;
268 double dev;
269 standard_deviation(times_g,mean,dev);
270
271 double mean_;
272 double dev_;
273 standard_deviation(times,mean_,dev_);
274
275 report_vector_funcs.graphs.put("performance.vector_layout(0).y.data.mean",mean_/mean);
276
277 // Deviation od x/y = x/y^2 dy + 1/y dx
278
279 report_vector_funcs.graphs.put("performance.vector_layout(0).y.data.dev",mean_/(mean*mean)*dev + dev_ / mean );
280
281 report_vector_funcs.graphs.put("performance.vector_layout(1).funcs.nele",NADD);
282 report_vector_funcs.graphs.put("performance.vector_layout(1).funcs.name","read_write_inte");
283
284 for (size_t i = 0 ; i < N_STAT+1 ; i++)
285 {
286 // create a vector
289
290 // Point
292 p.get<0>()[0] = 1.0;
293 p.get<0>()[1] = 2.0;
294
296 pa.get<0>() = 1.0;
297
298 pa.get<1>()[0] = 1.0;
299 pa.get<1>()[1] = 1.0;
300
301 pa.get<2>()[0][0] = 1.0;
302 pa.get<2>()[0][1] = 1.0;
303 pa.get<2>()[1][0] = 1.0;
304 pa.get<2>()[1][1] = 1.0;
305
306 // Add test
307
308 for (size_t i = 0 ; i < NADD ; i++)
309 {
310 v1.add(pa);
311 v2.add(p);
312 }
313
314 timer tg;
315 tg.start();
316
317 for (size_t i = 0 ; i < NADD ; i++)
318 {
319 read_write(v1,v2,i);
320 }
321
322 tg.stop();
323
324 times2_g[i] = tg.getwct();
325 int sz = v1.size();
326
327 timer tga;
328 tga.start();
329
330 double * prp0 = (double *)v1.getPointer<0>();
331 double * prp1 = (double *)v1.getPointer<1>();
332 double * prp2 = (double *)v1.getPointer<2>();
333
334 double * pos = (double *)v2.getPointer<0>();
335
336 for (size_t i = 0 ; i < NADD ; i++)
337 {
338 read_write_inte(pos,prp0,prp1,prp2,i,sz);
339 }
340
341 tga.stop();
342
343 times2[i] = tga.getwct();
344 }
345
346 double mean2;
347 double dev2;
348 standard_deviation(times2_g,mean2,dev2);
349
350 double mean2_;
351 double dev2_;
352 standard_deviation(times2,mean2_,dev2_);
353
354 report_vector_funcs.graphs.put("performance.vector_layout(1).y.data.mean",mean2_/mean2);
355
356 // Deviation od x/y = x/y^2 dy + 1/y dx
357
358 report_vector_funcs.graphs.put("performance.vector_layout(1).y.data.dev",mean2_/(mean2*mean2)*dev2 + dev2_ / mean2 );
359}
360
362
363BOOST_AUTO_TEST_CASE(vector_performance_write_report)
364{
365 // Create a graphs
366
367 report_vector_funcs.graphs.put("graphs.graph(0).type","line");
368 report_vector_funcs.graphs.add("graphs.graph(0).title","Vector add and get");
369 report_vector_funcs.graphs.add("graphs.graph(0).x.title","Tests");
370 report_vector_funcs.graphs.add("graphs.graph(0).y.title","Time seconds");
371 report_vector_funcs.graphs.add("graphs.graph(0).y.data(0).source","performance.vector(#).y.data.mean");
372 report_vector_funcs.graphs.add("graphs.graph(0).x.data(0).source","performance.vector(#).funcs.name");
373 report_vector_funcs.graphs.add("graphs.graph(0).y.data(0).title","Actual");
374 report_vector_funcs.graphs.add("graphs.graph(0).interpolation","lines");
375
376 report_vector_funcs.graphs.put("graphs.graph(1).type","line");
377 report_vector_funcs.graphs.add("graphs.graph(1).title","Vector read write");
378 report_vector_funcs.graphs.add("graphs.graph(1).x.title","Layout");
379 report_vector_funcs.graphs.add("graphs.graph(1).y.title","Time seconds");
380 report_vector_funcs.graphs.add("graphs.graph(1).y.data(0).source","performance.vector_layout(#).y.data.mean");
381 report_vector_funcs.graphs.add("graphs.graph(1).x.data(0).source","performance.vector_layout(#).funcs.name");
382 report_vector_funcs.graphs.add("graphs.graph(1).y.data(0).title","Actual");
383 report_vector_funcs.graphs.add("graphs.graph(1).interpolation","lines");
384
385 boost::property_tree::xml_writer_settings<std::string> settings(' ', 4);
386 boost::property_tree::write_xml("vector_performance_funcs.xml", report_vector_funcs.graphs,std::locale(),settings);
387
388 GoogleChart cg;
389
390 std::string file_xml_ref(test_dir);
391 file_xml_ref += std::string("/openfpm_data/vector_performance_funcs_ref.xml");
392
393 StandardXMLPerformanceGraph("vector_performance_funcs.xml",file_xml_ref,cg);
394
395 addUpdtateTime(cg,1,"data","vector_performance_funcs");
396
397 cg.write("vector_performance_funcs.html");
398}
399
400BOOST_AUTO_TEST_SUITE_END()
401
402#endif /* OPENFPM_DATA_SRC_VECTOR_VECTOR_PERFORMANCE_TEST_HPP_ */
403
404
405
406#endif /* VECTOR_PERFORMANCE_TEST_HPP_ */
Small class to produce graph with Google chart in HTML.
void write(std::string file)
It write the graphs on file in html format using Google charts.
This class allocate, and destroy CPU memory.
Test structure used for several test.
void sety(T y_)
set the y property
void setz(T z_)
set the z property
auto get() -> decltype(boost::fusion::at_c< i >(data))
getter method for a general property i
static const unsigned int v
v property is at position 4 in the boost::fusion::vector
void sets(T s_)
set the s property
void setx(T x_)
set the x property
static const unsigned int t
t property is at position 5 in the boost::fusion::vector
Implementation of 1-D std::vector like structure.
size_t size()
Stub size.
Class for cpu time benchmarking.
Definition timer.hpp:28
void stop()
Stop the timer.
Definition timer.hpp:119
void start()
Start the timer.
Definition timer.hpp:90
double getwct()
Return the elapsed real time.
Definition timer.hpp:130
aggregate of properties, from a list of object if create a struct that follow the OPENFPM native stru...
__device__ __host__ boost::mpl::at< type, boost::mpl::int_< i > >::type & get()
get the properties i
Transform the boost::fusion::vector into memory specification (memory_traits)