[GiNaC-list] Re: Airy functions

Alexei Sheplyakov varg at theor.jinr.ru
Mon Dec 24 11:13:33 CET 2007


On Sun, Dec 23, 2007 at 02:53:55PM +0100, Jan Bos wrote:

> To answer your questions a fundamental property of the algorithm is
> needed: 
> Depending on the argument a so called "steepest descent" contour is
> determined along which the integrand is monotonically decreasing. 

> This ensures that a) roundoff does not occur,

Roundoff errors occur (almost) always, unless you've got a magic hardware
which operates on real numbers (as opposed to floating point ones).
However, in a carefully designed code these errors are smaller (or of the
same order) than the errors of the algorithm itself (i.e. due to truncation
of a series, replacing an integral with a sum, etc). I don't quite understand
how your code solves this problem.

> b) the range over which the integrand has to be evaluated to obtain a certain
> accuracy is minimized and 

Anyway, the higher accuracy, the bigger this interval is. With IEEE floating
point numbers (or any other numbers with the fixed size of the mantissa)
arbitrary high accuracy does not make sense, so this is (probably) not a problem.
However, for arbitrary precision calculation the accuracy can be arbitrary high,
so the interval could be arbitrary big.

> c) the dominant factor (possibly very larger or small) can be placed outside
> the integral.

> Furthermore all derivatives at +/- infinity are zero. For functions with
> this property the trapezoidal rule has exponential convergence,
> that is, halving the stepsize more than doubles the number of significant
> digits of the result.

That would be correct in real number arithmetics. With floating point numbers
halving the step size can result in bigger roundoff error (becouse of doubled
number of operations). Thus, one need to find an optimal step size, so both
roundoff error and the error of algorithm itself are small enough. If
the target accuracy is very high[1], such a step size might not exist at all.
Hence my question -- is numerical integration suitable for *arbitary*
precision calculation?

[1] and/or the integrand is weird, i.e. sin(1/x) when x -> 0

> For an integrator with non-constant step size this property would be lost. 

I don't quite understand. Actually, if the function is monotoneous along
the integration contour, I'd expect adaptive integration algorithm to
converge faster (and give smaller roundoff errors). So, could you elaborate,

Best regards,

All science is either physics or stamp collecting.

-------------- next part --------------
A non-text attachment was scrubbed...
Name: not available
Type: application/pgp-signature
Size: 827 bytes
Desc: Digital signature
Url : http://www.cebix.net/pipermail/ginac-list/attachments/20071224/68d58023/attachment.pgp

More information about the GiNaC-list mailing list