1// Copyright 2016 Ismael Jimenez Martinez. All rights reserved.
2//
3// Licensed under the Apache License, Version 2.0 (the "License");
4// you may not use this file except in compliance with the License.
5// You may obtain a copy of the License at
6//
7//     http://www.apache.org/licenses/LICENSE-2.0
8//
9// Unless required by applicable law or agreed to in writing, software
10// distributed under the License is distributed on an "AS IS" BASIS,
11// WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
12// See the License for the specific language governing permissions and
13// limitations under the License.
14
15// Source project : https://github.com/ismaelJimenez/cpp.leastsq
16// Adapted to be used with google benchmark
17
18#include "benchmark/benchmark.h"
19
20#include <algorithm>
21#include <cmath>
22#include "check.h"
23#include "complexity.h"
24
25namespace benchmark {
26
27// Internal function to calculate the different scalability forms
28BigOFunc* FittingCurve(BigO complexity) {
29  switch (complexity) {
30    case oN:
31      return [](int64_t n) -> double { return static_cast<double>(n); };
32    case oNSquared:
33      return [](int64_t n) -> double { return std::pow(n, 2); };
34    case oNCubed:
35      return [](int64_t n) -> double { return std::pow(n, 3); };
36    case oLogN:
37      return [](int64_t n) { return log2(n); };
38    case oNLogN:
39      return [](int64_t n) { return n * log2(n); };
40    case o1:
41    default:
42      return [](int64_t) { return 1.0; };
43  }
44}
45
46// Function to return an string for the calculated complexity
47std::string GetBigOString(BigO complexity) {
48  switch (complexity) {
49    case oN:
50      return "N";
51    case oNSquared:
52      return "N^2";
53    case oNCubed:
54      return "N^3";
55    case oLogN:
56      return "lgN";
57    case oNLogN:
58      return "NlgN";
59    case o1:
60      return "(1)";
61    default:
62      return "f(N)";
63  }
64}
65
66// Find the coefficient for the high-order term in the running time, by
67// minimizing the sum of squares of relative error, for the fitting curve
68// given by the lambda expression.
69//   - n             : Vector containing the size of the benchmark tests.
70//   - time          : Vector containing the times for the benchmark tests.
71//   - fitting_curve : lambda expression (e.g. [](int64_t n) {return n; };).
72
73// For a deeper explanation on the algorithm logic, look the README file at
74// http://github.com/ismaelJimenez/Minimal-Cpp-Least-Squared-Fit
75
76LeastSq MinimalLeastSq(const std::vector<int64_t>& n,
77                       const std::vector<double>& time,
78                       BigOFunc* fitting_curve) {
79  double sigma_gn = 0.0;
80  double sigma_gn_squared = 0.0;
81  double sigma_time = 0.0;
82  double sigma_time_gn = 0.0;
83
84  // Calculate least square fitting parameter
85  for (size_t i = 0; i < n.size(); ++i) {
86    double gn_i = fitting_curve(n[i]);
87    sigma_gn += gn_i;
88    sigma_gn_squared += gn_i * gn_i;
89    sigma_time += time[i];
90    sigma_time_gn += time[i] * gn_i;
91  }
92
93  LeastSq result;
94  result.complexity = oLambda;
95
96  // Calculate complexity.
97  result.coef = sigma_time_gn / sigma_gn_squared;
98
99  // Calculate RMS
100  double rms = 0.0;
101  for (size_t i = 0; i < n.size(); ++i) {
102    double fit = result.coef * fitting_curve(n[i]);
103    rms += pow((time[i] - fit), 2);
104  }
105
106  // Normalized RMS by the mean of the observed values
107  double mean = sigma_time / n.size();
108  result.rms = sqrt(rms / n.size()) / mean;
109
110  return result;
111}
112
113// Find the coefficient for the high-order term in the running time, by
114// minimizing the sum of squares of relative error.
115//   - n          : Vector containing the size of the benchmark tests.
116//   - time       : Vector containing the times for the benchmark tests.
117//   - complexity : If different than oAuto, the fitting curve will stick to
118//                  this one. If it is oAuto, it will be calculated the best
119//                  fitting curve.
120LeastSq MinimalLeastSq(const std::vector<int64_t>& n,
121                       const std::vector<double>& time, const BigO complexity) {
122  CHECK_EQ(n.size(), time.size());
123  CHECK_GE(n.size(), 2);  // Do not compute fitting curve is less than two
124                          // benchmark runs are given
125  CHECK_NE(complexity, oNone);
126
127  LeastSq best_fit;
128
129  if (complexity == oAuto) {
130    std::vector<BigO> fit_curves = {oLogN, oN, oNLogN, oNSquared, oNCubed};
131
132    // Take o1 as default best fitting curve
133    best_fit = MinimalLeastSq(n, time, FittingCurve(o1));
134    best_fit.complexity = o1;
135
136    // Compute all possible fitting curves and stick to the best one
137    for (const auto& fit : fit_curves) {
138      LeastSq current_fit = MinimalLeastSq(n, time, FittingCurve(fit));
139      if (current_fit.rms < best_fit.rms) {
140        best_fit = current_fit;
141        best_fit.complexity = fit;
142      }
143    }
144  } else {
145    best_fit = MinimalLeastSq(n, time, FittingCurve(complexity));
146    best_fit.complexity = complexity;
147  }
148
149  return best_fit;
150}
151
152std::vector<BenchmarkReporter::Run> ComputeBigO(
153    const std::vector<BenchmarkReporter::Run>& reports) {
154  typedef BenchmarkReporter::Run Run;
155  std::vector<Run> results;
156
157  if (reports.size() < 2) return results;
158
159  // Accumulators.
160  std::vector<int64_t> n;
161  std::vector<double> real_time;
162  std::vector<double> cpu_time;
163
164  // Populate the accumulators.
165  for (const Run& run : reports) {
166    CHECK_GT(run.complexity_n, 0) << "Did you forget to call SetComplexityN?";
167    n.push_back(run.complexity_n);
168    real_time.push_back(run.real_accumulated_time / run.iterations);
169    cpu_time.push_back(run.cpu_accumulated_time / run.iterations);
170  }
171
172  LeastSq result_cpu;
173  LeastSq result_real;
174
175  if (reports[0].complexity == oLambda) {
176    result_cpu = MinimalLeastSq(n, cpu_time, reports[0].complexity_lambda);
177    result_real = MinimalLeastSq(n, real_time, reports[0].complexity_lambda);
178  } else {
179    result_cpu = MinimalLeastSq(n, cpu_time, reports[0].complexity);
180    result_real = MinimalLeastSq(n, real_time, result_cpu.complexity);
181  }
182  std::string benchmark_name =
183      reports[0].benchmark_name.substr(0, reports[0].benchmark_name.find('/'));
184
185  // Get the data from the accumulator to BenchmarkReporter::Run's.
186  Run big_o;
187  big_o.benchmark_name = benchmark_name + "_BigO";
188  big_o.iterations = 0;
189  big_o.real_accumulated_time = result_real.coef;
190  big_o.cpu_accumulated_time = result_cpu.coef;
191  big_o.report_big_o = true;
192  big_o.complexity = result_cpu.complexity;
193
194  // All the time results are reported after being multiplied by the
195  // time unit multiplier. But since RMS is a relative quantity it
196  // should not be multiplied at all. So, here, we _divide_ it by the
197  // multiplier so that when it is multiplied later the result is the
198  // correct one.
199  double multiplier = GetTimeUnitMultiplier(reports[0].time_unit);
200
201  // Only add label to mean/stddev if it is same for all runs
202  Run rms;
203  big_o.report_label = reports[0].report_label;
204  rms.benchmark_name = benchmark_name + "_RMS";
205  rms.report_label = big_o.report_label;
206  rms.iterations = 0;
207  rms.real_accumulated_time = result_real.rms / multiplier;
208  rms.cpu_accumulated_time = result_cpu.rms / multiplier;
209  rms.report_rms = true;
210  rms.complexity = result_cpu.complexity;
211  // don't forget to keep the time unit, or we won't be able to
212  // recover the correct value.
213  rms.time_unit = reports[0].time_unit;
214
215  results.push_back(big_o);
216  results.push_back(rms);
217  return results;
218}
219
220}  // end namespace benchmark
221