Visual Servoing Platform  version 3.6.1 under development (2024-02-13)
perfColVectorOperations.cpp
1 /****************************************************************************
2  *
3  * ViSP, open source Visual Servoing Platform software.
4  * Copyright (C) 2005 - 2023 by Inria. All rights reserved.
5  *
6  * This software is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  * See the file LICENSE.txt at the root directory of this source
11  * distribution for additional information about the GNU GPL.
12  *
13  * For using ViSP with software that can not be combined with the GNU
14  * GPL, please contact Inria about acquiring a ViSP Professional
15  * Edition License.
16  *
17  * See https://visp.inria.fr for more information.
18  *
19  * This software was developed at:
20  * Inria Rennes - Bretagne Atlantique
21  * Campus Universitaire de Beaulieu
22  * 35042 Rennes Cedex
23  * France
24  *
25  * If you have questions regarding the use of this file, please contact
26  * Inria at visp@inria.fr
27  *
28  * This file is provided AS IS with NO WARRANTY OF ANY KIND, INCLUDING THE
29  * WARRANTY OF DESIGN, MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE.
30  *
31  * Description:
32  * Benchmark column vector operations.
33  *
34 *****************************************************************************/
35 
36 #include <visp3/core/vpConfig.h>
37 
38 #ifdef VISP_HAVE_CATCH2
39 #define CATCH_CONFIG_ENABLE_BENCHMARKING
40 #define CATCH_CONFIG_RUNNER
41 #include <catch.hpp>
42 
43 #include <visp3/core/vpColVector.h>
44 
45 namespace
46 {
47 static bool g_runBenchmark = false;
48 static const std::vector<unsigned int> g_sizes = { 23, 127, 233, 419, 1153, 2749 };
49 
50 double getRandomValues(double min, double max) { return (max - min) * ((double)rand() / (double)RAND_MAX) + min; }
51 
52 vpColVector generateRandomVector(unsigned int rows, double min = -100, double max = 100)
53 {
54  vpColVector v(rows);
55 
56  for (unsigned int i = 0; i < v.getRows(); i++) {
57  v[i] = getRandomValues(min, max);
58  }
59 
60  return v;
61 }
62 
63 double stddev(const std::vector<double> &vec)
64 {
65  double sum = std::accumulate(vec.begin(), vec.end(), 0.0);
66  double mean = sum / vec.size();
67 
68  std::vector<double> diff(vec.size());
69  std::transform(vec.begin(), vec.end(), diff.begin(), [mean](double x) {
70  return x - mean; });
71  double sq_sum = std::inner_product(diff.begin(), diff.end(), diff.begin(), 0.0);
72  return std::sqrt(sq_sum / vec.size());
73 }
74 
75 double computeRegularSum(const vpColVector &v)
76 {
77  double sum = 0.0;
78 
79  for (unsigned int i = 0; i < v.getRows(); i++) {
80  sum += v[i];
81  }
82 
83  return sum;
84 }
85 
86 double computeRegularSumSquare(const vpColVector &v)
87 {
88  double sum_square = 0.0;
89 
90  for (unsigned int i = 0; i < v.getRows(); i++) {
91  sum_square += v[i] * v[i];
92  }
93 
94  return sum_square;
95 }
96 
97 double computeRegularStdev(const vpColVector &v)
98 {
99  double mean_value = computeRegularSum(v) / v.getRows();
100  double sum_squared_diff = 0.0;
101 
102  for (unsigned int i = 0; i < v.size(); i++) {
103  sum_squared_diff += (v[i] - mean_value) * (v[i] - mean_value);
104  }
105 
106  double divisor = (double)v.size();
107 
108  return std::sqrt(sum_squared_diff / divisor);
109 }
110 
111 std::vector<double> computeHadamard(const std::vector<double> &v1, const std::vector<double> &v2)
112 {
113  std::vector<double> result;
114  std::transform(v1.begin(), v1.end(), v2.begin(), std::back_inserter(result), std::multiplies<double>());
115  return result;
116 }
117 
118 bool almostEqual(const vpColVector &v1, const vpColVector &v2, double tol = 1e-9)
119 {
120  if (v1.getRows() != v2.getRows()) {
121  return false;
122  }
123 
124  for (unsigned int i = 0; i < v1.getRows(); i++) {
125  if (!vpMath::equal(v1[i], v2[i], tol)) {
126  return false;
127  }
128  }
129 
130  return true;
131 }
132 } // namespace
133 
134 TEST_CASE("Benchmark vpColVector::sum()", "[benchmark]")
135 {
136  // Sanity checks
137  {
138  const double val = 11;
139  vpColVector v(1, val);
140  CHECK(v.sum() == Approx(val).epsilon(std::numeric_limits<double>::epsilon()));
141  }
142  {
143  const unsigned int size = 11;
144  std::vector<double> vec(size);
145  vpColVector v(size);
146  for (size_t i = 0; i < 11; i++) {
147  vec[i] = 2. * i;
148  v[static_cast<unsigned int>(i)] = vec[i];
149  }
150  CHECK(v.sum() ==
151  Approx(std::accumulate(vec.begin(), vec.end(), 0.0)).epsilon(std::numeric_limits<double>::epsilon()));
152  }
153 
154  if (g_runBenchmark) {
155  for (auto size : g_sizes) {
156  vpColVector v = generateRandomVector(size);
157  std::vector<double> vec = v.toStdVector();
158 
159  std::ostringstream oss;
160  oss << "Benchmark vpColVector::sum() with size: " << size << " (ViSP)";
161  double vp_sum = 0;
162  BENCHMARK(oss.str().c_str())
163  {
164  vp_sum = v.sum();
165  return vp_sum;
166  };
167 
168  oss.str("");
169  oss << "Benchmark std::accumulate() with size: " << size << " (C++)";
170  double std_sum = 0;
171  BENCHMARK(oss.str().c_str())
172  {
173  std_sum = std::accumulate(vec.begin(), vec.end(), 0.0);
174  return std_sum;
175  };
176  CHECK(vp_sum == Approx(std_sum));
177 
178  oss.str("");
179  oss << "Benchmark naive sum() with size: " << size;
180  double naive_sum = 0;
181  BENCHMARK(oss.str().c_str())
182  {
183  naive_sum = computeRegularSum(v);
184  return naive_sum;
185  };
186  CHECK(naive_sum == Approx(std_sum));
187  }
188  }
189 }
190 
191 TEST_CASE("Benchmark vpColVector::sumSquare()", "[benchmark]")
192 {
193  // Sanity checks
194  {
195  const double val = 11;
196  vpColVector v(1, val);
197  CHECK(v.sumSquare() == Approx(val * val).epsilon(std::numeric_limits<double>::epsilon()));
198  }
199  {
200  const unsigned int size = 11;
201  std::vector<double> vec(size);
202  vpColVector v(size);
203  for (size_t i = 0; i < 11; i++) {
204  vec[i] = 2. * i;
205  v[static_cast<unsigned int>(i)] = vec[i];
206  }
207  CHECK(v.sumSquare() == Approx(std::inner_product(vec.begin(), vec.end(), vec.begin(), 0.0))
208  .epsilon(std::numeric_limits<double>::epsilon()));
209  }
210 
211  if (g_runBenchmark) {
212  for (auto size : g_sizes) {
213  vpColVector v = generateRandomVector(size);
214  std::vector<double> vec = v.toStdVector();
215 
216  std::ostringstream oss;
217  oss << "Benchmark vpColVector::sumSquare() with size: " << size << " (ViSP)";
218  double vp_sq_sum = 0;
219  BENCHMARK(oss.str().c_str())
220  {
221  vp_sq_sum = v.sumSquare();
222  return vp_sq_sum;
223  };
224 
225  oss.str("");
226  oss << "Benchmark std::inner_product with size: " << size << " (C++)";
227  double std_sq_sum = 0;
228  BENCHMARK(oss.str().c_str())
229  {
230  std_sq_sum = std::inner_product(vec.begin(), vec.end(), vec.begin(), 0.0);
231  return std_sq_sum;
232  };
233  CHECK(vp_sq_sum == Approx(std_sq_sum));
234 
235  oss.str("");
236  oss << "Benchmark naive sumSquare() with size: " << size;
237  double naive_sq_sum = 0;
238  BENCHMARK(oss.str().c_str())
239  {
240  naive_sq_sum = computeRegularSumSquare(v);
241  return naive_sq_sum;
242  };
243  CHECK(naive_sq_sum == Approx(std_sq_sum));
244  }
245  }
246 }
247 
248 TEST_CASE("Benchmark vpColVector::stdev()", "[benchmark]")
249 {
250  // Sanity checks
251  {
252  vpColVector v(2);
253  v[0] = 11;
254  v[1] = 16;
255  std::vector<double> vec = v.toStdVector();
256  CHECK(vpColVector::stdev(v) == Approx(stddev(vec)).epsilon(std::numeric_limits<double>::epsilon()));
257  }
258  {
259  const unsigned int size = 11;
260  std::vector<double> vec(size);
261  vpColVector v(size);
262  for (size_t i = 0; i < 11; i++) {
263  vec[i] = 2. * i;
264  v[static_cast<unsigned int>(i)] = vec[i];
265  }
266  CHECK(vpColVector::stdev(v) == Approx(stddev(vec)).epsilon(std::numeric_limits<double>::epsilon()));
267  }
268 
269  if (g_runBenchmark) {
270  for (auto size : g_sizes) {
271  vpColVector v = generateRandomVector(size);
272  std::vector<double> vec = v.toStdVector();
273 
274  std::ostringstream oss;
275  oss << "Benchmark vpColVector::stdev() with size: " << size << " (ViSP)";
276  double vp_stddev = 0;
277  BENCHMARK(oss.str().c_str())
278  {
279  vp_stddev = vpColVector::stdev(v);
280  return vp_stddev;
281  };
282 
283  oss.str("");
284  oss << "Benchmark C++ stddev impl with size: " << size << " (C++)";
285  double std_stddev = 0;
286  BENCHMARK(oss.str().c_str())
287  {
288  std_stddev = stddev(vec);
289  return std_stddev;
290  };
291  CHECK(vp_stddev == Approx(std_stddev));
292 
293  oss.str("");
294  oss << "Benchmark naive stdev() with size: " << size;
295  double naive_stddev = 0;
296  BENCHMARK(oss.str().c_str())
297  {
298  naive_stddev = computeRegularStdev(v);
299  return naive_stddev;
300  };
301  CHECK(naive_stddev == Approx(std_stddev));
302  }
303  }
304 }
305 
306 TEST_CASE("Benchmark vpColVector::hadamard()", "[benchmark]")
307 {
308  // Sanity checks
309  {
310  vpColVector v1(2), v2(2);
311  v1[0] = 11;
312  v1[1] = 16;
313  v2[0] = 8;
314  v2[1] = 23;
315  vpColVector res1 = v1.hadamard(v2);
316  vpColVector res2(computeHadamard(v1.toStdVector(), v2.toStdVector()));
317  CHECK(almostEqual(res1, res2));
318  }
319  {
320  const unsigned int size = 11;
321  std::vector<double> vec1(size), vec2(size);
322  for (size_t i = 0; i < 11; i++) {
323  vec1[i] = 2. * i;
324  vec2[i] = 3. * i + 5.;
325  }
326  vpColVector v1(vec1), v2(vec2);
327  vpColVector res1 = v1.hadamard(v2);
328  vpColVector res2(computeHadamard(v1.toStdVector(), v2.toStdVector()));
329  CHECK(almostEqual(res1, res2));
330  }
331 
332  if (g_runBenchmark) {
333  for (auto size : g_sizes) {
334  vpColVector v1 = generateRandomVector(size);
335  vpColVector v2 = generateRandomVector(size);
336  std::vector<double> vec1 = v1.toStdVector();
337  std::vector<double> vec2 = v2.toStdVector();
338 
339  std::ostringstream oss;
340  oss << "Benchmark vpColVector::hadamard() with size: " << size << " (ViSP)";
341  vpColVector vp_res;
342  BENCHMARK(oss.str().c_str())
343  {
344  vp_res = v1.hadamard(v2);
345  return vp_res;
346  };
347 
348  oss.str("");
349  oss << "Benchmark C++ element wise multiplication with size: " << size << " (C++)";
350  std::vector<double> std_res;
351  BENCHMARK(oss.str().c_str())
352  {
353  std_res = computeHadamard(vec1, vec2);
354  return std_res;
355  };
356  CHECK(almostEqual(vp_res, vpColVector(std_res)));
357  }
358  }
359 }
360 
361 int main(int argc, char *argv[])
362 {
363  // Set random seed explicitly to avoid confusion
364  // See: https://en.cppreference.com/w/cpp/numeric/random/srand
365  // If rand() is used before any calls to srand(), rand() behaves as if it was seeded with srand(1).
366  srand(1);
367 
368  Catch::Session session; // There must be exactly one instance
369 
370  // Build a new parser on top of Catch's
371  using namespace Catch::clara;
372  auto cli = session.cli() // Get Catch's composite command line parser
373  | Opt(g_runBenchmark) // bind variable to a new option, with a hint string
374  ["--benchmark"] // the option names it will respond to
375  ("run benchmark?"); // description string for the help output
376 
377 // Now pass the new composite back to Catch so it uses that
378  session.cli(cli);
379 
380  // Let Catch (using Clara) parse the command line
381  session.applyCommandLine(argc, argv);
382 
383  int numFailed = session.run();
384 
385  // numFailed is clamped to 255 as some unices only use the lower 8 bits.
386  // This clamping has already been applied, so just return it here
387  // You can also do any post run clean-up here
388  return numFailed;
389 }
390 #else
391 #include <iostream>
392 
393 int main() { return EXIT_SUCCESS; }
394 #endif
unsigned int size() const
Return the number of elements of the 2D array.
Definition: vpArray2D.h:286
unsigned int getRows() const
Definition: vpArray2D.h:284
Implementation of column vector and the associated operations.
Definition: vpColVector.h:163
vpColVector hadamard(const vpColVector &v) const
double sumSquare() const
std::vector< double > toStdVector() const
static double stdev(const vpColVector &v, bool useBesselCorrection=false)
double sum() const
static bool equal(double x, double y, double threshold=0.001)
Definition: vpMath.h:449