Accuracy of Integral with epsilon

[lrdheat]

I was surprised that my TI-30X Pro MathPrint came up with a more accurate solution to the integral (2*e^x)/(cosh (x) + sinh (x)) from -10 to 10 using an epsilon of 10^-4 or greater. I was naively thinking that the stock epsilon of 10^-5 or smaller would yield a more accurate result. The actual answer is 40. My TI cx ii was accurate to 5 digits as is the TI-30X Pro MathPrint and CASIO fx-991X. The CASIO fx-CG50 is accurate to 6 digits, the HP Prime in home mode thought about it a bit, and was accurate to 7 digits. The WP 34S on my SwissMicros in double precision quickly produced an answer accurate to 27 digits!

[toml_12953]

(10-08-2021 09:57 PM)lrdheat Wrote:  I was surprised that my TI-30X Pro MathPrint came up with a more accurate solution to the integral (2*e^x)/(cosh (x) + sinh (x)) from -10 to 10 using an epsilon of 10^-4 or greater. I was naively thinking that the stock epsilon of 10^-5 or smaller would yield a more accurate result. The actual answer is 40. My TI cx ii was accurate to 5 digits as is the TI-30X Pro MathPrint and CASIO fx-991X. The CASIO fx-CG50 is accurate to 6 digits, the HP Prime in home mode thought about it a bit, and was accurate to 7 digits. The WP 34S on my SwissMicros in double precision quickly produced an answer accurate to 27 digits!

I don't understand what you mean about Prime being accurate to 7 digits. I get 40 exactly on Prime in home mode.

[Thomas Okken]

I don't even understand what's supposed to be challenging about that integral, since that function should be the constant 2...

[toml_12953]

(10-08-2021 10:18 PM)Thomas Okken Wrote:  I don't even understand what's supposed to be challenging about that integral, since that function should be the constant 2...

Maybe he was hoping we wouldn't notice.

[lrdheat]

Did not think the integral was remarkable, but was fascinated by how different calculators handled it as written. My Prime in home shows 40.0000000143 (!). Impressed with the WP 34S.

[toml_12953]

(10-08-2021 11:28 PM)lrdheat Wrote:  Did not think the integral was remarkable, but was fascinated by how different calculators handled it as written. My Prime in home shows 40.0000000143 (!). Impressed with the WP 34S.

What firmware version are you using? I use the 2021 10 01 version on a G2.

I wonder why the big difference between yours and mine.

[robve]

(10-08-2021 11:28 PM)lrdheat Wrote:  Did not think the integral was remarkable, but was fascinated by how different calculators handled it as written. My Prime in home shows 40.0000000143 (!). Impressed with the WP 34S.

This integral is efficiently integrated with the Tanh-Sinh method (HP forum post). This method is particularly suitable for your function.

I did a quick check with my implementation: for a specified accuracy of eps=10^-4 Tanh-Sinh only performs 27 integrand evaluations to converge to 40 exactly, or at least very close (depends on the floating point format, single/double IEEE 754 or BCD etc.)

The WP 34S uses Tanh-Sinh. That implementation can be improved though. Tanh-Sinh is not considered a general integration method as it may perform quite poorly in certain (known) cases.

- Rob

[Valentin Albillo]

(10-09-2021 01:41 AM)robve Wrote:  

(10-08-2021 11:28 PM)lrdheat Wrote:  Did not think the integral was remarkable, but was fascinated by how different calculators handled it as written. My Prime in home shows 40.0000000143 (!). Impressed with the WP 34S.

This integral is efficiently integrated with the Tanh-Sinh method (HP forum post). This method is particularly suitable for your function.

I did a quick check with my implementation: for a specified accuracy of eps=10^-4 Tanh-Sinh only performs 27 integrand evaluations to converge to 40 exactly, or at least very close (depends on the floating point format, single/double IEEE 754 or BCD etc.)

Someone here is not paying attention.

Are you saying that a complex state-of-the-art Tanh-Sinh method requires 27 evaluations of the integrand (the constant 2) to integrate the constant 2 from -10 to 10 ? Really !?

What a lousy implementation !

V.

[Paul Dale]

(10-09-2021 02:22 AM)Valentin Albillo Wrote:  Are you saying that a complex state-of-the-art Tanh-Sinh method requires 27 evaluations of the integrand (the constant 2) to integrate the constant 2 from -10 to 10 ? Really !?

It's neither complex nor state of the art.

How many function evaluations should be performed?

Clearly one evaluation is sufficient to get the correct answer here. Two might be prudent to verify that it is a constant function. Three perhaps?


Pauli

[lrdheat]

I’m on 2.1.14588

The Prime does not seem to directly know how to recognize that the integrand is a constant, even when simplify is applied to the integrand.

The WP 34S was quick/accurate, regardless of number of iterations. The MathPrint and 991X took ~15 seconds.

I messed with a few home settings, but still came up with the 40.0000000143 answer on the Prime G2

[lrdheat]

What surprises me most is my CASIO Classpad readily simplifies the integrand to a constant (2), yet, when I integrate the integral with the unsimplified integrand, it fails to recognize the constant integrand, and produces 40.00000008

[toml_12953]

(10-09-2021 03:13 AM)lrdheat Wrote:  I’m on 2.1.14588

The Prime does not seem to directly know how to recognize that the integrand is a constant, even when simplify is applied to the integrand.

The WP 34S was quick/accurate, regardless of number of iterations. The MathPrint and 991X took ~15 seconds.

I messed with a few home settings, but still came up with the 40.0000000143 answer on the Prime G2

You should get the latest firmware, then: 2.1.14596 (2021 10 01) It looks like that has been fixed.

[Thomas Okken]

(10-09-2021 02:27 AM)Paul Dale Wrote:  Clearly one evaluation is sufficient to get the correct answer here. Two might be prudent to verify that it is a constant function. Three perhaps?

It seems what's really being tested here is how the integrator handles the slight noisiness of the function at the low end of the range, when sinh(x) and cosh(x) become large but with opposite signs, so numerical cancellation becomes an issue.

[toml_12953]

(10-09-2021 03:13 AM)lrdheat Wrote:  I’m on 2.1.14588

The Prime does not seem to directly know how to recognize that the integrand is a constant, even when simplify is applied to the integrand.

The WP 34S was quick/accurate, regardless of number of iterations. The MathPrint and 991X took ~15 seconds.

I messed with a few home settings, but still came up with the 40.0000000143 answer on the Prime G2

What function are you using to calculate the integral? I use int.

[robve]

(10-09-2021 01:41 AM)robve Wrote:  

(10-08-2021 11:28 PM)lrdheat Wrote:  Did not think the integral was remarkable, but was fascinated by how different calculators handled it as written. My Prime in home shows 40.0000000143 (!). Impressed with the WP 34S.

I did a quick check with my implementation: for a specified accuracy of eps=10^-4 Tanh-Sinh only performs 27 integrand evaluations to converge to 40 exactly, or at least very close (depends on the floating point format, single/double IEEE 754 or BCD etc.)

The WP 34S uses Tanh-Sinh. That implementation can be improved though.

WP 34S results with a IEEE 754 based implementation matching the Python code given by the author:
eps=10^-4: abs error = 0.0051939; est. error = 0.0155524; evals = 7
eps=10^-5: abs error = 0.000130107; est. error = 0.0156835; evals = 9
eps=10^-6: abs error = 0.000121986; est. error = 6.30233e-06; evals = 17
eps=10^-7: abs error = 2.14894e-06; est. error = 3.30639e-06; evals = 19
eps=10^-8: abs error = 2.14894e-06; est. error = 3.30639e-06; evals = 19
eps=10^-9: abs error = 5.81537e-09; est. error = 3.35997e-06; evals = 21

Note that the given eps is not necessarily smaller than the abs error obtained ("under delivers"). Only the last result delivers a decent estimate. Whether or not the convergence check is satisfied with the estimated error is a detail that may differ in the implementation(s) of Tanh-Sinh, thus eps=10^4 may produce a better answer for some, but to get the better answer 19, 21 or more evaluations are needed to get close enough.

The qthsh implementation (discussed on the HP forum thread) of the same algorithm with some differences has 27 evaluations for eps=10^4 with small abs error 5.63299e-09. Qthsh approaches integration boundaries more consistently, evaluating additional boundary points by handling singularities instead of terminating early. Otherwise the algorithms are the same. This explains why qthsh has 27 versus WP 34S 21 evaluations. Like the WP 34S, qthsh was verified and compared to WP 34S and other algorithms, including a verification of accuracy with over 800 integrals.

- Rob

[lrdheat]

The Prime and TI Nspire CAS emulators both (with integrand not simplified) produce an exact result of 40 (not 40.). Surprisingly, the CASIO Classpad, while readily simplifying the integrand by itself to “2”, it’s CAS, when presented with the integral with the integrand not simplified will not give an exact result.

[John Keith]

(10-09-2021 08:47 AM)Thomas Okken Wrote:  It seems what's really being tested here is how the integrator handles the slight noisiness of the function at the low end of the range, when sinh(x) and cosh(x) become large but with opposite signs, so numerical cancellation becomes an issue.

Indeed. Evaluating the integrand numerically, the HP 50g returns 2. for positive values of x but 1.99999690584 for x = -10.

[robve]

(10-09-2021 02:27 AM)Paul Dale Wrote:  

(10-09-2021 02:22 AM)Valentin Albillo Wrote:  Are you saying that a complex state-of-the-art Tanh-Sinh method requires 27 evaluations of the integrand (the constant 2) to integrate the constant 2 from -10 to 10 ? Really !?

It's neither complex nor state of the art.

How many function evaluations should be performed?

Clearly one evaluation is sufficient to get the correct answer here. Two might be prudent to verify that it is a constant function. Three perhaps?

There is never a guarantee that a given function is constant by probing a few points. That should be obvious, no?

Naive implementations of Gauss-Legendre quadrature methods may stop at two or three points to fit a constant (polynomial). However, it would be detrimental to give up that early. Evidently, calculators run quite a bit longer on this function for a good reason.

My point about Tanh-Sinh was simply that this method doesn't require a huge number of points to converge accurately and it is used in the WP 34S... so make your own conclusions. Tanh-Sinh simply works well for this type of function based on its syntactical structure (duh). More importantly, the function is NOT numerically constant in its present non-simplified form, which is to be expected: noise increases quickly towards -10 and beyond. This noise is already in the order of 10^-2 at x=-18. Tabulating the error in Excel from -10 to -8.9 gives:

-5.31248E-08
-6.01538E-08
-5.39791E-08
-1.77564E-08
-1.86315E-08
-1.64802E-08
-5.18465E-10
-8.54069E-09
-1.47008E-08
1.40072E-08
-6.26628E-09
-5.77673E-09

This amount of noise is sufficient to trip up quadrature methods. Of course, the noise may differ with non-IEEE 754 double precision and different implementations of the constituent functions EXP, SINH and COSH. However, a general consequence of the noise and floating point limitations, you can either get lucky to get to the exact result 40 with a few points or unlucky, which can cost you a great deal of time wasted to evaluate points in the presence of noise. The WP 34S and qthsh points reported show exactly what I mean. Also, attempting to push the accuracy beyond 10^8 (or about) is pointless. There is too much noise to make definitive conclusions.

- Rob

[robve]

(10-09-2021 02:22 AM)Valentin Albillo Wrote:  Someone here is not paying attention.

Are you saying that a complex state-of-the-art Tanh-Sinh method requires 27 evaluations of the integrand (the constant 2) to integrate the constant 2 from -10 to 10 ? Really !?

What a lousy implementation !

You do realize that it is a fallacy to assume exp(x) and cosh(x)+sinh(x) equate numerically, in addition to algebraically?

For one, sinh(x) typically uses expm1(x)=exp(x)-1 for higher accuracy sinh(x) = (z+z/(z+1))/2 with z=expm1(x). This asymmetrically affects cosh(x)+sinh(x) for x<0. The difference grows with smaller x<0, which is easy to see when we plug in expm1(x).

- Rob

[Albert Chan]

(10-09-2021 08:37 PM)robve Wrote:  There is never a guarantee that a given function is constant by probing a few points. That should be obvious, no?

Tanh-Sinh (and other DE variants) does not see a constant function ... it see a bell-shaped curve.

CAS> 'int(1,x,-1,1)' (x = tanh(sinh(t)))

\(\int _{-\infty }^{+\infty }\mathrm{cosh}\left(t\right)\cdot (1-\left(\mathrm{tanh}\left(\mathrm{sinh}\left(t\right)\right)\right)^{2})\, dt\)

(10-09-2021 09:50 PM)robve Wrote:  You do realize that it is a fallacy to assume exp(x) and cosh(x)+sinh(x) equate numerically, in addition to algebraically?

It is interesting that my MAPM implementation of raw_exp() actually use above "identity"
Of course, |x| is limited to tiny value, below 1E-4, so catastrophic cancellation is not an issue.

sinh(x) converge with half as many terms, compared with exp(x)

exp(x) = x + x^2/2! + x^3/3! + ...
sinh(x) = x + x^3/3! + x^5/5! + ...