ROOT  6.06/09
Reference Guide
stressVdt.cxx
Go to the documentation of this file.
1 /// Simple program to benchmark vdt accuracy and cpu performance.
2 #include <vector>
3 #include <iostream>
4 #include <cmath> //for log2
5 #include <assert.h>
6 #include <limits>
7 #include <iostream>
8 #include <iomanip>
9 #include <map>
10 #include <string>
11 
12 #include "vdt/vdtMath.h"
13 
14 #include "TStopwatch.h"
15 #include "TRandom3.h"
16 #include "TH1F.h"
17 #include "TCanvas.h"
18 #include "TError.h"
19 
20 const double RANGE=3000.;
21 const uint32_t SIZE= 16777216;
22 
23 //------------------------------------------------------------------------------
24 // Not good for floating point, but just to get the bit difference
25 template <class T>
26 uint64_t fp2uint (T /*x*/)
27 {
28  T::not_implemented; // "Static assert" in C++03
29  return 0;
30 }
31 
32 template <>
33 uint64_t fp2uint<double> (double x)
34 {
35  return vdt::details::dp2uint64(x);
36 }
37 
38 template <>
39 uint64_t fp2uint<float> (float x)
40 {
41  return vdt::details::sp2uint32(x);
42 }
43 
44 //------------------------------------------------------------------------------
45 /// Returns most significative different bit
46 template <class T>
47 inline uint32_t diffbit(const T a,const T b )
48 {
49  uint64_t ia = fp2uint<T>(a);
50  uint64_t ib = fp2uint<T>(b);
51  uint64_t c = ia>ib? ia-ib : ib-ia;
52  return log2(c)+1;
53 }
54 
55 //------------------------------------------------------------------------------
56 // This allows to vectorise on very modern compilers (>=gcc 4.7)
57 // Note how the templating mechanism allows the compiler to *inline* the
58 // function. It is much more efficient than a void ptr or std::function!
59 template <typename T, typename F>
60 inline void calculateValues(F mathFunc,
61  const std::vector<T>& inputVector,
62  std::vector<T>& outputVector)
63 {
64 
65  const uint32_t size = inputVector.size();
66 
67  for (unsigned int i=0;i<size;++i){
68  outputVector[i]=mathFunc(inputVector[i]);
69  }
70 
71 }
72 
73 ////////////////////////////////////////////////////////////////////////////////
74 
75 template <typename T>
76 void compareOutputs(const std::vector<T>& inputVector1,
77  const std::vector<T>& inputVector2,
78  std::vector<uint32_t>& outputVector)
79 {
80  assert(inputVector1.size()==inputVector2.size() &&
81  inputVector1.size()==outputVector.size());
82 
83  const uint32_t size = inputVector1.size();
84 
85  for (unsigned int i=0;i<size;++i)
86  outputVector[i]=diffbit(inputVector1[i],inputVector2[i]);
87 }
88 
89 //------------------------------------------------------------------------------
90 
96 
97 template <typename T>
98 void fillRandom(std::vector<T>& randomV,
99  const rangeType theRangeType)
100 {
101  // Yeah, well, maybe it can be done better. But this is not a tutorial about
102  // random generation!
103  const uint32_t size=randomV.size();
104  static TRandom3 rndmGenerator(123);
105  T* arr = &(randomV[0]);
106  rndmGenerator.RndmArray(size,arr);
107  if (kReal == theRangeType ) for (uint32_t i=0;i<size;++i) randomV[i]=(randomV[i]-0.5)*2*RANGE;
108  if (kExp == theRangeType ) for (uint32_t i=0;i<size;++i) randomV[i]=(randomV[i]-0.5)*2*705.;
109  if (kExpf == theRangeType ) for (uint32_t i=0;i<size;++i) randomV[i]=(randomV[i]-0.5)*2*85.;
110  if (kRealPlus == theRangeType ) for (uint32_t i=0;i<size;++i) randomV[i]=randomV[i]*RANGE+0.000001;
111  if (km1p1 == theRangeType ) for (uint32_t i=0;i<size;++i) randomV[i]=(randomV[i]-0.5)*2;
112 
113 }
114 
115 //------------------------------------------------------------------------------
116 
117 template<typename T>
119  const std::vector<T>& VDTVals,
120  const std::vector<T>& SystemVals)
121 {
122  const uint32_t size = VDTVals.size();
123  std::vector<uint32_t> diff(size);
124  compareOutputs(VDTVals,SystemVals,diff);
125  uint32_t theDiff=0;
126  for (uint32_t i =0;i<size;++i){
127  theDiff=diff[i];
128  histo.Fill(theDiff);
129  }
130 }
131 
132 ////////////////////////////////////////////////////////////////////////////////
133 
134 template <typename T, typename F>
135 inline double measureTiming(F mathFunc,
136  const std::vector<T>& inputVector,
137  std::vector<T>& outputVector)
138 {
140  timer.Start();
141  calculateValues<T>(mathFunc,inputVector,outputVector);
142  timer.Stop();
143  return timer.RealTime();
144 }
145 
146 //------------------------------------------------------------------------------
147 
148  /*
149  Reference values on a Intel(R) Core(TM) i7 CPU 950 @ 3.07GHz
150  gcc 4.8.2 -Ofast
151 
152  Expf 0.110623 0.821863 7.4294 3.35636 6
153  Sinf 0.151893 9.83798 64.7692 0.235055 9
154  Cosf 0.11277 9.87508 87.5683 0.234271 8
155  Tanf 0.167273 9.84792 58.8733 0.519784 9
156  Atanf 0.0855272 0.529288 6.18854 0.366963 2
157  Logf 0.114023 0.465541 4.08287 0.270999 2
158  Isqrtf 0.0328619 0.275043 8.36965 4.35237 7
159  Asinf 0.0958891 0.415733 4.33556 0.60167 3
160  Acosf 0.099067 0.470179 4.74607 0.480427 10
161  Exp 0.204585 0.64904 3.17247 0.137142 2
162  Sin 0.327537 1.5099 4.60986 0.253579 2
163  Cos 0.299601 1.5038 5.01933 0.253664 2
164  Tan 0.276369 2.13009 7.7074 0.351065 5
165  Atan 0.355532 0.902413 2.5382 0.326134 2
166  Log 0.244172 1.25513 5.14034 0.385112 2
167  Isqrt 0.167836 0.283752 1.69065 0.453692 2
168  Asin 0.40379 0.869315 2.15289 0.318644 2
169  Acos 0.392566 0.864706 2.2027 0.391922 11
170  */
171 
172  // The reference values: speedup and accuracy. Some contingency is given
173 struct staticInitHelper{
174  std::map<std::string, std::pair<float,uint32_t> > referenceValues;
175 
176  staticInitHelper()
177  {
178  referenceValues["Expf"] = std::make_pair(1.f,8);
179  referenceValues["Sinf"] = std::make_pair(1.f,11);
180  referenceValues["Cosf"] = std::make_pair(1.f,10);
181  referenceValues["Tanf"] = std::make_pair(1.f,11);
182  referenceValues["Atanf"] = std::make_pair(1.f,4);
183  referenceValues["Logf"] = std::make_pair(1.f,4);
184  referenceValues["Isqrtf"]= std::make_pair(1.f,9);
185  referenceValues["Asinf"] = std::make_pair(1.f,5);
186  referenceValues["Acosf"] = std::make_pair(1.f,12);
187  referenceValues["Exp"] = std::make_pair(1.f,4);
188  referenceValues["Sin"] = std::make_pair(1.f,4);
189  referenceValues["Cos"] = std::make_pair(1.f,4);
190  referenceValues["Tan"] = std::make_pair(1.f,7);
191  referenceValues["Atan"] = std::make_pair(1.f,4);
192  referenceValues["Log"] = std::make_pair(1.f,4);
193  referenceValues["Isqrt"] = std::make_pair(.4f,4); // Fix fluctuation on x86_64-slc5-gcc47
194  referenceValues["Asin"] = std::make_pair(1.f,4);
195  referenceValues["Acos"] = std::make_pair(1.f,13);
196  }
197 } gbl;
198 
199 template <typename T, typename F1, typename F2>
200 inline void compareFunctions(const std::string& label,
201  F1 vdtFunc,
202  F2 systemFunc,
203  const std::vector<T>& inputVector,
204  std::vector<T>& outputVectorVDT,
205  std::vector<T>& outputVectorSystem,
206  float& speedup,
207  uint32_t& maxdiffBit,
208  TH1F& histo)
209 {
210  double timeVdt = measureTiming<T>(vdtFunc,inputVector,outputVectorVDT);
211  double timeSystem = measureTiming<T>(systemFunc,inputVector,outputVectorSystem);
212  std::string name(label);
213  std::string title(label);
214  title+=";Diff. Bit;#";
215  histo.Reset();
216  histo.SetName(label.c_str());
217  histo.SetTitle(title.c_str());
218  treatBinDiffHisto(histo,outputVectorVDT,outputVectorSystem);
219  double meandiffBit = histo.GetMean();
220  maxdiffBit = 0;
221  const uint32_t xmax=histo.GetXaxis()->GetXmax();
222 
223  for (uint32_t i=1;i<=xmax;i++){
224  if ( histo.GetBinContent(i) > 0.f )
225  maxdiffBit=i-1;
226  }
227 
228  speedup = timeSystem/timeVdt ;
229 
230  std::cout << std::setw(8)
231  << label << std::setw(10)
232  << timeVdt << std::setw(10)
233  << timeSystem << std::setw(15)
234  << speedup << std::setw(15)
235  << meandiffBit << std::setw(15)
236  << maxdiffBit << std::endl;
237 
238  // Draw it
239  TCanvas c;
240  c.cd();
241  histo.SetLineColor(kBlue);
242  histo.SetLineWidth(2);
243  histo.Draw("hist");
244  c.SetLogy();
245  name+=".png";
246  Int_t oldErrorIgnoreLevel = gErrorIgnoreLevel; // we know we are printing..
247  gErrorIgnoreLevel=1001;
248  c.Print(name.c_str());
249  gErrorIgnoreLevel=oldErrorIgnoreLevel;
250 }
251 
252 //
253 void checkFunction(const std::string& label,float /*speedup*/, uint32_t maxdiffBit)
254 {
255 // Remove check on the speed as routinely this program is ran on virtual build nodes
256 // and several factors may cause fluctuations in the result.
257 // if (gbl.referenceValues[label].first > speedup)
258 // std::cerr << "Note " << label << " was slower than the system library.\n";
259  if (gbl.referenceValues[label].second < maxdiffBit)
260  std::cerr << "WARNING " << label << " is too inaccurate!\n";
261 }
262 
263 //------------------------------------------------------------------------------
264 // Some helpers
265 inline float isqrtf(float x) {return 1.f/sqrt(x);};
266 inline double isqrt(double x) {return 1./sqrt(x);};
267 
268 //------------------------------------------------------------------------------
269 // Artificially chunck the work to help the compiler manage everything.
270 // If all the content is moved to a single function, the vectorization breaks.
271 // Same holds for the check of speed and accuracy. Technically it could be part
272 // of compareFunctions.
273 void spStep1()
274 {
275  std::vector<float> VDTVals(SIZE);
276  std::vector<float> SystemVals(SIZE);
277  std::vector<float> realNumbers(SIZE);
278  TH1F histo("bitDiffHisto","willbechanged",32,0,32);
279 
280  float speedup;
281  uint32_t maxdiffBit;
282 
283  fillRandom(realNumbers,kExpf);
284  compareFunctions<float>("Expf", vdt::fast_expf, expf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
285  checkFunction("Expf",speedup, maxdiffBit);
286 
287  fillRandom(realNumbers,kReal);
288  compareFunctions<float>("Sinf", vdt::fast_sinf, sinf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
289  checkFunction("Sinf",speedup, maxdiffBit);
290  compareFunctions<float>("Cosf", vdt::fast_cosf, cosf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
291  checkFunction("Cosf",speedup, maxdiffBit);
292  compareFunctions<float>("Tanf", vdt::fast_tanf, tanf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
293  checkFunction("Tanf",speedup, maxdiffBit);
294  compareFunctions<float>("Atanf", vdt::fast_atanf, atanf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
295  checkFunction("Atanf",speedup, maxdiffBit);
296 }
297 
298 void spStep2()
299 {
300  std::vector<float> VDTVals(SIZE);
301  std::vector<float> SystemVals(SIZE);
302  std::vector<float> realNumbers(SIZE);
303  TH1F histo("bitDiffHisto","willbechanged",32,0,32);
304 
305  float speedup;
306  uint32_t maxdiffBit;
307 
308  fillRandom(realNumbers,kRealPlus);
309  compareFunctions<float>("Logf", vdt::fast_logf, logf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
310  checkFunction("Logf",speedup, maxdiffBit);
311  compareFunctions<float>("Isqrtf", vdt::fast_isqrtf, isqrtf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
312  checkFunction("Isqrtf",speedup, maxdiffBit);
313 }
314 
315 void spStep3()
316 {
317  std::vector<float> VDTVals(SIZE);
318  std::vector<float> SystemVals(SIZE);
319  std::vector<float> realNumbers(SIZE);
320  TH1F histo("bitDiffHisto","willbechanged",32,0,32);
321 
322  float speedup;
323  uint32_t maxdiffBit;
324 
325  fillRandom(realNumbers,km1p1);
326  compareFunctions<float>("Asinf", vdt::fast_asinf, asinf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
327  checkFunction("Asinf",speedup, maxdiffBit);
328  compareFunctions<float>("Acosf", vdt::fast_acosf, acosf, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
329  checkFunction("Acosf",speedup, maxdiffBit);
330 }
331 
332 void dpStep1()
333 {
334  std::vector<double> VDTVals(SIZE);
335  std::vector<double> SystemVals(SIZE);
336  std::vector<double> realNumbers(SIZE);
337  TH1F histo("bitDiffHisto","willbechanged",64,0,64);
338 
339  float speedup;
340  uint32_t maxdiffBit;
341 
342  fillRandom(realNumbers,kExp);
343  compareFunctions<double>("Exp", vdt::fast_exp, exp, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
344  checkFunction("Exp",speedup, maxdiffBit);
345  fillRandom(realNumbers,kReal);
346  compareFunctions<double>("Sin", vdt::fast_sin, sin, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
347  checkFunction("Sin",speedup, maxdiffBit);
348  compareFunctions<double>("Cos", vdt::fast_cos, cos, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
349  checkFunction("Cos",speedup, maxdiffBit);
350  compareFunctions<double>("Tan", vdt::fast_tan, tan, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
351  checkFunction("Tan",speedup, maxdiffBit);
352  compareFunctions<double>("Atan", vdt::fast_atan, atan, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
353  checkFunction("Atan",speedup, maxdiffBit);
354 }
355 
356 void dpStep2()
357 {
358  std::vector<double> VDTVals(SIZE);
359  std::vector<double> SystemVals(SIZE);
360  std::vector<double> realNumbers(SIZE);
361  TH1F histo("bitDiffHisto","willbechanged",64,0,64);
362 
363  float speedup;
364  uint32_t maxdiffBit;
365 
366  fillRandom(realNumbers,kRealPlus);
367  compareFunctions<double>("Log", vdt::fast_log, log, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
368  checkFunction("Log",speedup, maxdiffBit);
369  compareFunctions<double>("Isqrt", vdt::fast_isqrt, isqrt, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
370  checkFunction("Isqrt",speedup, maxdiffBit);
371 }
372 
373 void dpStep3()
374 {
375  std::vector<double> VDTVals(SIZE);
376  std::vector<double> SystemVals(SIZE);
377  std::vector<double> realNumbers(SIZE);
378  TH1F histo("bitDiffHisto","willbechanged",64,0,64);
379 
380  float speedup;
381  uint32_t maxdiffBit;
382 
383  fillRandom(realNumbers,km1p1);
384  compareFunctions<double>("Asin", vdt::fast_asin, asin, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
385  checkFunction("Asin",speedup, maxdiffBit);
386  compareFunctions<double>("Acos", vdt::fast_acos, acos, realNumbers, VDTVals, SystemVals, speedup, maxdiffBit, histo);
387  checkFunction("Acos",speedup, maxdiffBit);
388 }
389 //------------------------------------------------------------------------------
390 
391 int main(){
392 
393  std::cout << "Test performed on " << SIZE << " random numbers\n"
394  << std::setw(8)
395  << "Name" << std::setw(10)
396  << "VDT (s)" << std::setw(10)
397  << "Sys (s)" << std::setw(15)
398  << "Speedup" << std::setw(15)
399  << "<diff Bit>" << std::setw(15)
400  << "max(diff Bit)" << std::endl;
401 
402  // Single precision ----
403  spStep1();
404  spStep2();
405  spStep3();
406 
407  // Double precision ----
408  dpStep1();
409  dpStep2();
410  dpStep3();
411 
412  return 0;
413 
414 }
415 
416 
417 
418 
419 
420 
421 
422 
423 
424 
425 
426 
427 
428 
429 
430 
431 
432 
433 
434 
435 
436 
437 
438 
439 
440 
441 
442 
443 
444 
445 
446 
447 
448 
449 
450 
451 
452 
453 
454 
455 
456 
457 
458 
459 
460 
461 
462 
463 
464 
void spStep3()
Definition: stressVdt.cxx:315
virtual void SetLineWidth(Width_t lwidth)
Definition: TAttLine.h:57
virtual Int_t Fill(Double_t x)
Increment bin with abscissa X by 1.
Definition: TH1.cxx:3159
float fast_expf(float initial_x)
Exponential Function single precision.
Definition: exp.h:118
float fast_asinf(float x)
Single Precision asin.
Definition: asin.h:157
float fast_tanf(float x)
Definition: tan.h:130
Random number generator class based on M.
Definition: TRandom3.h:29
Double_t RealTime()
Stop the stopwatch (if it is running) and return the realtime (in seconds) passed between the start a...
Definition: TStopwatch.cxx:108
virtual Double_t GetBinContent(Int_t bin) const
Return content of bin number bin.
Definition: TH1.cxx:4629
float fast_sinf(float x)
Definition: sin.h:41
void Start(Bool_t reset=kTRUE)
Start the stopwatch.
Definition: TStopwatch.cxx:56
R__EXTERN Int_t gErrorIgnoreLevel
Definition: TError.h:107
double fast_exp(double initial_x)
Exponential Function double precision.
Definition: exp.h:70
double T(double x)
Definition: ChebyshevPol.h:34
#define assert(cond)
Definition: unittest.h:542
TVirtualPad * cd(Int_t subpadnumber=0)
Set current canvas & pad.
Definition: TCanvas.cxx:659
uint32_t sp2uint32(float x)
Converts a float to an int.
const double RANGE
Simple program to benchmark vdt accuracy and cpu performance.
Definition: stressVdt.cxx:20
1-D histogram with a float per channel (see TH1 documentation)}
Definition: TH1.h:570
double measureTiming(F mathFunc, const std::vector< T > &inputVector, std::vector< T > &outputVector)
Definition: stressVdt.cxx:135
struct staticInitHelper gbl
int Int_t
Definition: RtypesCore.h:41
TArc * a
Definition: textangle.C:12
double fast_asin(double x)
Double Precision asin.
Definition: asin.h:120
virtual void Print(const char *filename="") const
Save Pad contents in a file in one of various formats.
Definition: TPad.cxx:4138
uint32_t diffbit(const T a, const T b)
Returns most significative different bit.
Definition: stressVdt.cxx:47
float isqrtf(float x)
Definition: stressVdt.cxx:265
double cos(double)
void treatBinDiffHisto(TH1F &histo, const std::vector< T > &VDTVals, const std::vector< T > &SystemVals)
Definition: stressVdt.cxx:118
virtual void Reset(Option_t *option="")
Reset.
Definition: TH1.cxx:9501
const uint32_t SIZE
Definition: stressVdt.cxx:21
double sqrt(double)
double fast_acos(double x)
Definition: asin.h:208
TStopwatch timer
Definition: pirndm.C:37
void Stop()
Stop the stopwatch.
Definition: TStopwatch.cxx:75
Double_t x[n]
Definition: legend1.C:17
double acos(double)
void dpStep2()
Definition: stressVdt.cxx:356
double fast_atan(double x)
Fast Atan implementation double precision.
Definition: atan.h:86
double fast_cos(double x)
Double precision cosine: just call sincos.
Definition: cos.h:22
void compareFunctions(const std::string &label, F1 vdtFunc, F2 systemFunc, const std::vector< T > &inputVector, std::vector< T > &outputVectorVDT, std::vector< T > &outputVectorSystem, float &speedup, uint32_t &maxdiffBit, TH1F &histo)
Definition: stressVdt.cxx:200
#define F1(x, y, z)
Definition: TMD5.cxx:265
void dpStep3()
Definition: stressVdt.cxx:373
double sin(double)
void checkFunction(const std::string &label, float, uint32_t maxdiffBit)
Definition: stressVdt.cxx:253
uint64_t fp2uint< double >(double x)
Definition: stressVdt.cxx:33
double fast_log(double x)
Definition: log.h:89
virtual void SetLineColor(Color_t lcolor)
Definition: TAttLine.h:54
#define F(x, y, z)
virtual Double_t GetMean(Int_t axis=1) const
For axis = 1,2 or 3 returns the mean value of the histogram along X,Y or Z axis.
Definition: TH1.cxx:7014
void spStep2()
Definition: stressVdt.cxx:298
double fast_sin(double x)
Double precision sine: just call sincos.
Definition: sin.h:37
virtual void Draw(Option_t *option="")
Draw this histogram with options.
Definition: TH1.cxx:2878
float fast_isqrtf(float x)
Two iterations.
Definition: sqrt.h:85
uint64_t dp2uint64(double x)
Converts a double to an unsigned long long.
#define F2(x, y, z)
Definition: TMD5.cxx:266
void spStep1()
Definition: stressVdt.cxx:273
float xmax
Definition: THbookFile.cxx:93
double asin(double)
void dpStep1()
Definition: stressVdt.cxx:332
virtual void RndmArray(Int_t n, Float_t *array)
Return an array of n random numbers uniformly distributed in ]0,1].
Definition: TRandom3.cxx:137
virtual void SetName(const char *name)
Change the name of this histogram.
Definition: TH1.cxx:8288
The Canvas class.
Definition: TCanvas.h:48
void fillRandom(std::vector< T > &randomV, const rangeType theRangeType)
Definition: stressVdt.cxx:98
rangeType
Definition: stressVdt.cxx:91
double f(double x)
uint64_t fp2uint(T)
Definition: stressVdt.cxx:26
double isqrt(double x)
Definition: stressVdt.cxx:266
uint64_t fp2uint< float >(float x)
Definition: stressVdt.cxx:39
Double_t GetXmax() const
Definition: TAxis.h:138
double atan(double)
int main()
Definition: stressVdt.cxx:391
#define name(a, b)
Definition: linkTestLib0.cpp:5
double tan(double)
float fast_logf(float x)
Definition: log.h:172
Definition: Rtypes.h:61
virtual void SetTitle(const char *title)
Change (i.e.
Definition: TH1.cxx:6268
float fast_cosf(float x)
Definition: cos.h:26
void compareOutputs(const std::vector< T > &inputVector1, const std::vector< T > &inputVector2, std::vector< uint32_t > &outputVector)
Definition: stressVdt.cxx:76
double exp(double)
double fast_tan(double x)
Double precision tangent implementation.
Definition: tan.h:84
void calculateValues(F mathFunc, const std::vector< T > &inputVector, std::vector< T > &outputVector)
Definition: stressVdt.cxx:60
double log(double)
TAxis * GetXaxis()
Definition: TH1.h:319
float fast_atanf(float xx)
Fast Atan implementation single precision.
Definition: atan.h:129
double fast_isqrt(double x)
Four iterations.
Definition: sqrt.h:55
float fast_acosf(float x)
Definition: asin.h:212
Stopwatch class.
Definition: TStopwatch.h:30
virtual void SetLogy(Int_t value=1)
Set Lin/Log scale for Y.
Definition: TPad.cxx:5318