|
|
|
Functions | |
| double | cosint (double x) |
| Private version of the cosine and sine integral. | |
| template<typename F > | |
| int | dzero (double a, double b, double &x0, double &rv, double eps, int mxf, F func) |
| double | expint (double x) |
| Private version of the sine integral. | |
| double | f1 (double x, double const *h_) |
| Private version of the exponential integral. | |
| double | f2 (double x, double const *h_) |
| void | sincosint (double x, double &sint, double &cint) |
| double | sinint (double x) |
| Private version of the cosine integral. | |
| double VVIObjDetails::cosint | ( | double | x | ) |
Private version of the cosine and sine integral.
Definition at line 194 of file VVIObj.cc.
References funct::cos(), funct::log(), AlCaHLTBitMon_ParallelJobs::p, createTree::pp, lumiQueryAPI::q, funct::sin(), x, detailsBasic3DVector::y, and zero.
{
// Initialized data
const double zero = 0.;
const double one = 1.;
const double two = 2.;
const double eight = 8.;
const double ce = .57721566490153;
const double c__[14] = { 1.9405491464836,.9413409132865,
-.579845034293,.3091572011159,-.0916101792208,.0164437407515,
-.0019713091952,1.692538851e-4,-1.09393296e-5,5.522386e-7,
-2.23995e-8,7.465e-10,-2.08e-11,5e-13 };
const double p[23] = { .96074783975204,-.0371138962124,
.00194143988899,-1.7165988425e-4,2.112637753e-5,-3.27163257e-6,
6.0069212e-7,-1.2586794e-7,2.932563e-8,-7.45696e-9,2.04105e-9,
-5.9502e-10,1.8323e-10,-5.921e-11,1.997e-11,-7e-12,2.54e-12,
-9.5e-13,3.7e-13,-1.4e-13,6e-14,-2e-14,1e-14 };
const double q[20] = { .98604065696238,-.0134717382083,
4.5329284117e-4,-3.067288652e-5,3.13199198e-6,-4.2110196e-7,
6.907245e-8,-1.318321e-8,2.83697e-9,-6.7329e-10,1.734e-10,
-4.787e-11,1.403e-11,-4.33e-12,1.4e-12,-4.7e-13,1.7e-13,-6e-14,
2e-14,-1e-14 };
// System generated locals
double d__1;
// Local variables
double h__;
int i__;
double r__, y, b0, b1, b2, pp, qq, alfa;
// If x==0, return same
if (x == zero) {
return zero;
}
if (fabs(x) <= eight) {
y = x / eight;
// Computing 2nd power
d__1 = y;
h__ = two * (d__1 * d__1) - one;
alfa = -two * h__;
b1 = zero;
b2 = zero;
for (i__ = 13; i__ >= 0; --i__) {
b0 = c__[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
b1 = ce + log((fabs(x))) - b0 + h__ * b2;
} else {
r__ = one / x;
y = eight * r__;
// Computing 2nd power
d__1 = y;
h__ = two * (d__1 * d__1) - one;
alfa = -two * h__;
b1 = zero;
b2 = zero;
for (i__ = 22; i__ >= 0; --i__) {
b0 = p[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
pp = b0 - h__ * b2;
b1 = zero;
b2 = zero;
for (i__ = 19; i__ >= 0; --i__) {
b0 = q[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
qq = b0 - h__ * b2;
b1 = r__ * (qq * sin(x) - r__ * pp * cos(x));
}
return b1;
} // cosint
| int VVIObjDetails::dzero | ( | double | a, |
| double | b, | ||
| double & | x0, | ||
| double & | rv, | ||
| double | eps, | ||
| int | mxf, | ||
| F | func | ||
| ) |
Definition at line 571 of file VVIObj.cc.
References f, f1(), f2(), connectstrParser::f3, benchmark_cfg::fa, benchmark_cfg::fb, createTree::ff, max(), and min.
Referenced by VVIObj::VVIObj().
{
/* System generated locals */
double d__1, d__2, d__3, d__4;
// Local variables
double f1, f2, f3, u1, u2, x1, x2, u3, u4, x3, ca, cb, cc, fa, fb, ee, ff;
int mc;
double xa, xb, fx, xx, su4;
xa = std::min(a,b);
xb = std::max(a,b);
fa = func(xa);
fb = func(xb);
if (fa * fb > 0.) {
rv = (xb - xa) * -2;
x0 = 0.;
return 1;
}
mc = 0;
L1:
x0 = (xa + xb) * .5;
rv = x0 - xa;
ee = eps * (fabs(x0) + 1);
if (rv <= ee) {
rv = ee;
ff = func(x0);
return 0;
}
f1 = fa;
x1 = xa;
f2 = fb;
x2 = xb;
L2:
fx = func(x0);
++mc;
if (mc > mxf) {
rv = (d__1 = xb - xa, fabs(d__1)) * -.5;
x0 = 0.;
return 0;
}
if (fx * fa > 0.) {
xa = x0;
fa = fx;
} else {
xb = x0;
fb = fx;
}
L3:
u1 = f1 - f2;
u2 = x1 - x2;
u3 = f2 - fx;
u4 = x2 - x0;
if (u2 == 0. || u4 == 0.) {goto L1;}
f3 = fx;
x3 = x0;
u1 /= u2;
u2 = u3 / u4;
ca = u1 - u2;
cb = (x1 + x2) * u2 - (x2 + x0) * u1;
cc = (x1 - x0) * f1 - x1 * (ca * x1 + cb);
if (ca == 0.) {
if (cb == 0.) {goto L1;}
x0 = -cc / cb;
} else {
u3 = cb / (ca * 2);
u4 = u3 * u3 - cc / ca;
if (u4 < 0.) {goto L1;}
su4 = fabs(u4);
if (x0 + u3 < 0.f) {su4 = -su4;}
x0 = -u3 + su4;
}
if (x0 < xa || x0 > xb) {goto L1;}
// Computing MIN
d__3 = (d__1 = x0 - x3, fabs(d__1)), d__4 = (d__2 = x0 - x2, fabs(d__2));
rv = std::min(d__3,d__4);
ee = eps * (fabs(x0) + 1);
if (rv > ee) {
f1 = f2;
x1 = x2;
f2 = f3;
x2 = x3;
goto L2;
}
fx = func(x0);
if (fx == 0.) {
rv = ee;
ff = func(x0);
return 0;
}
if (fx * fa < 0.) {
xx = x0 - ee;
if (xx <= xa) {
rv = ee;
ff = func(x0);
return 0;
}
ff = func(xx);
fb = ff;
xb = xx;
} else {
xx = x0 + ee;
if (xx >= xb) {
rv = ee;
ff = func(x0);
return 0;
}
ff = func(xx);
fa = ff;
xa = xx;
}
if (fx * ff > 0.) {
mc += 2;
if (mc > mxf) {
rv = (d__1 = xb - xa, fabs(d__1)) * -.5;
x0 = 0.;
return 0;
}
f1 = f3;
x1 = x3;
f2 = fx;
x2 = x0;
x0 = xx;
fx = ff;
goto L3;
}
/* L4: */
rv = ee;
ff = func(x0);
return 0;
} // dzero
| double VVIObjDetails::expint | ( | double | x | ) |
Private version of the sine integral.
Definition at line 449 of file VVIObj.cc.
References funct::exp(), funct::log(), p1, p2, p3, p4, q1, q2, v, x, detailsBasic3DVector::y, and zero.
Referenced by f2(), and VVIObj::VVIObj().
{
// Initialized data
const double zero = 0.;
const double q2[7] = { .10340013040487,3.319092135933,
20.449478501379,41.280784189142,32.426421069514,10.041164382905,
1. };
const double p3[6] = { -2.3909964453136,-147.98219500504,
-254.3763397689,-119.55761038372,-19.630408535939,-.9999999999036
};
const double q3[6] = { 177.60070940351,530.68509610812,
462.23027156148,156.81843364539,21.630408494238,1. };
const double p4[8] = { -8.6693733995107,-549.14226552109,
-4210.0161535707,-249301.39345865,-119623.66934925,
-22174462.775885,3892804.213112,-391546073.8091 };
const double q4[8] = { 34.171875,-1607.0892658722,35730.029805851,
-483547.43616216,4285596.2461175,-24903337.574054,89192576.757561,
-165254299.72521 };
const double a1[8] = { -2.1808638152072,-21.901023385488,
9.3081638566217,25.076281129356,-33.184253199722,60.121799083008,
-43.253113287813,1.0044310922808 };
const double b1[8] = { 0.,3.9370770185272,300.89264837292,
-6.2504116167188,1003.6743951673,14.325673812194,2736.2411988933,
.52746885196291 };
const double a2[8] = { -3.4833465360285,-18.65454548834,
-8.2856199414064,-32.34673303054,17.960168876925,1.7565631546961,
-1.9502232128966,.99999429607471 };
const double b2[8] = { 0.,69.500065588743,57.283719383732,
25.777638423844,760.76114800773,28.951672792514,-3.4394226689987,
1.0008386740264 };
const double a3[6] = { -27.780928934438,-10.10479081576,
-9.1483008216736,-5.0223317461851,-3.0000077799358,
1.0000000000704 };
const double one = 1.;
const double b3[6] = { 0.,122.39993926823,2.7276100778779,
-7.1897518395045,-2.9990118065262,1.999999942826 };
const double two = 2.;
const double three = 3.;
const double x0 = .37250741078137;
const double xl[6] = { -24.,-12.,-6.,0.,1.,4. };
const double p1[5] = { 4.293125234321,39.894153870321,
292.52518866921,425.69682638592,-434.98143832952 };
const double q1[5] = { 1.,18.899288395003,150.95038744251,
568.05252718987,753.58564359843 };
const double p2[7] = { .43096783946939,6.9052252278444,
23.019255939133,24.378408879132,9.0416155694633,.99997957705159,
4.656271079751e-7 };
/* Local variables */
double v, y, ap, bp, aq, dp, bq, dq;
if (x <= xl[0]) {
ap = a3[0] - x;
for ( int i__ = 2; i__ <= 5; ++i__) {
/* L1: */
ap = a3[i__ - 1] - x + b3[i__ - 1] / ap;
}
y = exp(-x) / x * (one - (a3[5] + b3[5] / ap) / x);
} else if (x <= xl[1]) {
ap = a2[0] - x;
for ( int i__ = 2; i__ <= 7; ++i__) {
ap = a2[i__ - 1] - x + b2[i__ - 1] / ap;
}
y = exp(-x) / x * (a2[7] + b2[7] / ap);
} else if (x <= xl[2]) {
ap = a1[0] - x;
for ( int i__ = 2; i__ <= 7; ++i__) {
ap = a1[i__ - 1] - x + b1[i__ - 1] / ap;
}
y = exp(-x) / x * (a1[7] + b1[7] / ap);
} else if (x < xl[3]) {
v = -two * (x / three + one);
bp = zero;
dp = p4[0];
for ( int i__ = 2; i__ <= 8; ++i__) {
ap = bp;
bp = dp;
dp = p4[i__ - 1] - ap + v * bp;
}
bq = zero;
dq = q4[0];
for ( int i__ = 2; i__ <= 8; ++i__) {
aq = bq;
bq = dq;
dq = q4[i__ - 1] - aq + v * bq;
}
y = -log(-x / x0) + (x + x0) * (dp - ap) / (dq - aq);
} else if (x == xl[3]) {
return zero;
} else if (x < xl[4]) {
ap = p1[0];
aq = q1[0];
for ( int i__ = 2; i__ <= 5; ++i__) {
ap = p1[i__ - 1] + x * ap;
aq = q1[i__ - 1] + x * aq;
}
y = -log(x) + ap / aq;
} else if (x <= xl[5]) {
y = one / x;
ap = p2[0];
aq = q2[0];
for ( int i__ = 2; i__ <= 7; ++i__) {
ap = p2[i__ - 1] + y * ap;
aq = q2[i__ - 1] + y * aq;
}
y = exp(-x) * ap / aq;
} else {
y = one / x;
ap = p3[0];
aq = q3[0];
for ( int i__ = 2; i__ <= 6; ++i__) {
ap = p3[i__ - 1] + y * ap;
aq = q3[i__ - 1] + y * aq;
}
y = exp(-x) * y * (one + y * ap / aq);
}
return y;
} // expint
| double VVIObjDetails::f1 | ( | double | x, |
| double const * | h_ | ||
| ) | [inline] |
Private version of the exponential integral.
Definition at line 32 of file VVIObj.cc.
References funct::log(), and x.
Referenced by dzero(), and VVIObj::VVIObj().
| double VVIObjDetails::f2 | ( | double | x, |
| double const * | h_ | ||
| ) | [inline] |
| void VVIObjDetails::sincosint | ( | double | x, |
| double & | sint, | ||
| double & | cint | ||
| ) |
Definition at line 346 of file VVIObj.cc.
References funct::cos(), funct::log(), AlCaHLTBitMon_ParallelJobs::p, createTree::pp, lumiQueryAPI::q, funct::sin(), x, detailsBasic3DVector::y, and zero.
Referenced by VVIObj::VVIObj().
{
// Initialized data
const double zero = 0.;
const double one = 1.;
const double two = 2.;
const double eight = 8.;
const double ce = .57721566490153;
const double pih = 1.5707963267949;
const double s__[14] = { 1.9522209759531,-.6884042321257,
.4551855132256,-.1804571236838,.0410422133759,-.0059586169556,
6.001427414e-4,-4.44708329e-5,2.5300782e-6,-1.141308e-7,4.1858e-9,
-1.273e-10,3.3e-12,-1e-13 };
const double c__[14] = { 1.9405491464836,.9413409132865,
-.579845034293,.3091572011159,-.0916101792208,.0164437407515,
-.0019713091952,1.692538851e-4,-1.09393296e-5,5.522386e-7,
-2.23995e-8,7.465e-10,-2.08e-11,5e-13 };
const double p[23] = { .96074783975204,-.0371138962124,
.00194143988899,-1.7165988425e-4,2.112637753e-5,-3.27163257e-6,
6.0069212e-7,-1.2586794e-7,2.932563e-8,-7.45696e-9,2.04105e-9,
-5.9502e-10,1.8323e-10,-5.921e-11,1.997e-11,-7e-12,2.54e-12,
-9.5e-13,3.7e-13,-1.4e-13,6e-14,-2e-14,1e-14 };
const double q[20] = { .98604065696238,-.0134717382083,
4.5329284117e-4,-3.067288652e-5,3.13199198e-6,-4.2110196e-7,
6.907245e-8,-1.318321e-8,2.83697e-9,-6.7329e-10,1.734e-10,
-4.787e-11,1.403e-11,-4.33e-12,1.4e-12,-4.7e-13,1.7e-13,-6e-14,
2e-14,-1e-14 };
// System generated locals
double d__1;
// Local variables
double h__;
int i__;
double r__, y, b0, b1, b2, pp, qq, alfa;
sint=0;
cint=0;
if (fabs(x) <= eight) {
y = x / eight;
// Computing 2nd power
d__1 = y;
h__ = two * (d__1 * d__1) - one;
alfa = -two * h__;
// cos
if (x!=0) {
b1 = zero;
b2 = zero;
for (i__ = 13; i__ >= 0; --i__) {
b0 = c__[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
cint = ce + log((fabs(x))) - b0 + h__ * b2;
}
// sin
b1 = zero;
b2 = zero;
for (i__ = 13; i__ >= 0; --i__) {
b0 = s__[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
sint = y * (b0 - b2);
} else {
r__ = one / x;
y = eight * r__;
// Computing 2nd power
d__1 = y;
h__ = two * (d__1 * d__1) - one;
alfa = -two * h__;
b1 = zero;
b2 = zero;
for (i__ = 22; i__ >= 0; --i__) {
b0 = p[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
pp = b0 - h__ * b2;
b1 = zero;
b2 = zero;
for (i__ = 19; i__ >= 0; --i__) {
b0 = q[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
qq = b0 - h__ * b2;
// cos
cint = r__ * (qq * sin(x) - r__ * pp * cos(x));
// sin
d__1 = pih;
if(x < 0.) d__1 = -d__1;
sint = d__1 - r__ * (r__ * pp * sin(x) + qq * cos(x));
}
}
| double VVIObjDetails::sinint | ( | double | x | ) |
Private version of the cosine integral.
Definition at line 272 of file VVIObj.cc.
References funct::cos(), AlCaHLTBitMon_ParallelJobs::p, createTree::pp, lumiQueryAPI::q, asciidump::s, funct::sin(), x, detailsBasic3DVector::y, and zero.
{
// Initialized data
const double zero = 0.;
const double one = 1.;
const double two = 2.;
const double eight = 8.;
const double pih = 1.5707963267949;
const double s[14] = { 1.9522209759531,-.6884042321257,
.4551855132256,-.1804571236838,.0410422133759,-.0059586169556,
6.001427414e-4,-4.44708329e-5,2.5300782e-6,-1.141308e-7,4.1858e-9,
-1.273e-10,3.3e-12,-1e-13 };
const double p[23] = { .96074783975204,-.0371138962124,
.00194143988899,-1.7165988425e-4,2.112637753e-5,-3.27163257e-6,
6.0069212e-7,-1.2586794e-7,2.932563e-8,-7.45696e-9,2.04105e-9,
-5.9502e-10,1.8323e-10,-5.921e-11,1.997e-11,-7e-12,2.54e-12,
-9.5e-13,3.7e-13,-1.4e-13,6e-14,-2e-14,1e-14 };
const double q[20] = { .98604065696238,-.0134717382083,
4.5329284117e-4,-3.067288652e-5,3.13199198e-6,-4.2110196e-7,
6.907245e-8,-1.318321e-8,2.83697e-9,-6.7329e-10,1.734e-10,
-4.787e-11,1.403e-11,-4.33e-12,1.4e-12,-4.7e-13,1.7e-13,-6e-14,
2e-14,-1e-14 };
// System generated locals
double d__1;
// Local variables
double h__;
int i__;
double r__, y, b0, b1, b2, pp, qq, alfa;
if (fabs(x) <= eight) {
y = x / eight;
d__1 = y;
h__ = two * (d__1 * d__1) - one;
alfa = -two * h__;
b1 = zero;
b2 = zero;
for (i__ = 13; i__ >= 0; --i__) {
b0 = s[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
b1 = y * (b0 - b2);
} else {
r__ = one / x;
y = eight * r__;
d__1 = y;
h__ = two * (d__1 * d__1) - one;
alfa = -two * h__;
b1 = zero;
b2 = zero;
for (i__ = 22; i__ >= 0; --i__) {
b0 = p[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
pp = b0 - h__ * b2;
b1 = zero;
b2 = zero;
for (i__ = 19; i__ >= 0; --i__) {
b0 = q[i__] - alfa * b1 - b2;
b2 = b1;
b1 = b0;
}
qq = b0 - h__ * b2;
d__1 = fabs(pih);
if(x < 0.) d__1 = -d__1;
b1 = d__1 - r__ * (r__ * pp * sin(x) + qq * cos(x));
}
return b1;
} // sinint
1.7.3