forked from root-project/root
-
Notifications
You must be signed in to change notification settings - Fork 0
/
TFormulaVecTests.h
118 lines (91 loc) · 3.22 KB
/
TFormulaVecTests.h
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
#include <TF1.h>
#include <TF2.h>
#include <iostream>
#include <TString.h>
bool verbose = true;
//#define VEC_IN_CTOR
bool CheckValues(const TString & testName, double r, double r0) {
bool ret = TMath::AreEqualAbs(r,r0,1.E-12);
if (!ret) {
std::cout << "Test failed for " << testName << " : " << r << " " << r0 << std::endl;
}
else
if (verbose) std::cout << testName << "\t ok\n";
return ret;
}
typedef double ( * FreeFunc1D) (double );
bool testVec1D(TF1 * f1, const TString & formula, FreeFunc1D func, double x ) {
auto y = f1->Eval(x);
auto y0 = func(x);
bool ret = CheckValues(formula, y, y0);
// check passing double_v interface
#ifdef R__HAS_VECCORE
ROOT::Double_v vx = x;
ROOT::Double_v vy = f1->EvalPar(&vx, nullptr);
double y2 = vecCore::Get(vy,0);
ret &= CheckValues(formula+TString("_v"), y2, y0);
#endif
return ret;
}
bool testVec1D(const TString & formula, FreeFunc1D func, double x = 1.) {
// test first by vectorizing in the constructor
auto f1 = new TF1("f",formula,0,1,"VEC");
bool ret = testVec1D(f1,formula,func,x);
// test by vectorizing afterwards calling SetVectorized
auto f2 = new TF1("f2",formula);
f2->SetVectorized(true);
//std::cout << "set vectprization" << std::endl;
//f2->Print("V");
ret &= testVec1D(f2,formula,func,x);
// test by removing vectorization
f2->SetVectorized(false);
ret &= testVec1D(f1,formula,func,x);
return ret;
}
typedef double ( * FreeFunc2D) (double, double );
bool testVec2D(TF2 * f1, const TString & formula, FreeFunc2D func, double x, double y ) {
auto r = f1->Eval(x,y);
auto r0 = func(x,y);
bool ret = CheckValues(formula,r,r0);
// check passing double_v interface
#ifdef R__HAS_VECCORE
ROOT::Double_v vx[2] = { x, y};
ROOT::Double_v vy = f1->EvalPar(vx, nullptr);
double r2 = vecCore::Get(vy,0);
ret &= CheckValues(formula+TString("_v"), r2, r0);
#endif
return ret;
}
bool testVec2D(const TString & formula, FreeFunc2D func, double x = 1., double y = 1.) {
auto f1 = new TF2("f",formula,0,1,0,1,"VEC");
bool ret = testVec2D(f1, formula, func, x, y);
auto f2 = new TF2("f",formula);
f2->SetVectorized(true);
ret &= testVec2D(f2,formula,func,x,y);
return ret;
}
double constant_function(double ) { return 3; }
bool testVecFormula() {
bool ok = true;
ok &= testVec1D("3+[0]",constant_function, 3.333);
ok &= testVec1D("3",constant_function, 3.333);
ok &= testVec1D("sin(x)",std::sin,1);
ok &= testVec1D("cos(x)",std::cos,1);
ok &= testVec1D("exp(x)",std::exp,1);
ok &= testVec1D("log(x)",std::log,2);
ok &= testVec1D("log10(x)",TMath::Log10,2);
ok &= testVec1D("tan(x)",std::tan,1);
//ok &= testVec1D("sinh(x)",std::sinh,1);
//ok &= testVec1D("cosh(x)",std::cosh,1);
//ok &= testVec1D("tanh(x)",std::tanh,1);
ok &= testVec1D("asin(x)",std::asin,.1);
ok &= testVec1D("acos(x)",std::acos,.1);
ok &= testVec1D("atan(x)",std::atan,.1);
ok &= testVec1D("sqrt(x)",std::sqrt,2);
ok &= testVec1D("abs(x)",std::abs,-1);
ok &= testVec2D("pow(x,y)",std::pow,2,3);
ok &= testVec2D("min(x,y)",TMath::Min,2,3);
ok &= testVec2D("max(x,y)",TMath::Max,2,3);
ok &= testVec2D("atan2(x,y)",TMath::ATan2,2,3);
return ok;
}