(PC-12xx~14xx) qthsh Tanh-Sinh quadrature

[robve]

(04-19-2021 10:16 PM)Albert Chan Wrote:  

(04-18-2021 12:56 PM)robve Wrote:  There is one bad case exp(-3*(x*x)-4/(x*x)) with 228 more points. The predictor pushes too hard on this one where the diff is close to zero already at the begin of the search, meaning there is no need to search for d in this case.

Note: this quote was removed, but it is instructive to see why it failed ...

Sharp eye to notice Rather than posting every small change in a new post which is harder and more time consuming to go through, this small change in this post avoids this problematic case for the reason I had already mentioned earlier in that post, the case when the initial diffs are too close to zero.

This small update also has another advantage. It allows the search for an optimized d from any initially-specified d instead of d=1. Perhaps this means that the method may have a fixed point that is the optimal d by invoking the predictor repeatedly with improved d and by using smaller steps. Interesting to investigate.

(04-19-2021 10:16 PM)Albert Chan Wrote:  And, there are many more bad cases, say, g(x) = exp(-m*x), m ≥ 1

Nope. The diff sign does not change for this integrand. And for larger m the search quits early.

(04-19-2021 10:16 PM)Albert Chan Wrote:  Perhaps estimate of d based on peak of f(v)*v is safer ?

I have not found a single counter example yet that does not work with the predictor. I also spent many hours to try various improvements and test them with >208 integrals. I have about a dozen tweaks of the method with mixed results. Not interesting to share.

A small change was needed to avoid the "bad" case and others similar to it, which are not difficult to recognize from the valid cases by inspecting the diff h sequence. Note that v (now more appropriately renamed r) should not be too small initially and is r=2 for the \( 2^i \) sequence. I also tried regula falsi to improve d, but unfortunately this does not always help. I also found that larger estimates of d (in the 1000s) tend to overshoot the optimal d. In this case starting with a larger d and then repeating the predictor recursively should improve the estimate closer to optimal. However, keeping things simple pays off, because moving d at least in the right direction helps convergence, even when we don't hit the optimal d, so we don't evaluate more points than necessary (or more than Exp-Sinh needs!). In any case, the geometric sequence of r can be \( 2^i \) or \( e^i \) or another to take larger steps to reduce the number of evaluations. \( 2^i \) is just cheap to compute.

(04-19-2021 10:16 PM)Albert Chan Wrote:  f(x) = exp(-3*(x*x)-4/(x*x))
Searching for peak of f(v)*v (we can also do numerically)

I already tried searching for the "cliff" (or point after the peak). I had mentioned the cliff before as relevant. But it is easy to see that the search should continue until the sign changes. Only at the point when the sign changes the following comparison can be safely made to determine the direction that d should move in:

if (fabs(fl) < fabs(fr))
  d /= r;
else
  d *= r;

Doing this at the peak, or right after, fails the method, causing d to move in the wrong direction. However, the peak or cliff can be closer to where the optimal d is located for large d (in the 100s or 1000s), before the sign changes, but this adds more complexity. A simple condition is to check that the diff h is not zero and in that case move r back one step:

if (h != 0)
  r /= 2;

If the diff is zero, we are at a reasonably good or optimal r to adjust d. This can happen very early for example 1/(1+(x*x)) and log((1+(x*x))/x)/(1+(x*x)) or late with for example x/sinh(x). For the last integrand we get d=4 (131 points, est.rel.err=2e-11) which improves the accuracy with the same number of points. Applying a regula falsi correction to r (instead of r/=2 see the code) we get d=6.7 (69 points, est.rel.err=6e-10). That's incredible

- Rob