Groups | Search | Server Info | Keyboard shortcuts | Login | Register [http] [https] [nntp] [nntps]

Groups > comp.graphics.apps.gnuplot > #1313

Re: Airy function wrong

Followups directed to: comp.graphics.apps.gnuplot

Show all headers | View raw

Alex van der Spek wrote:

> Just finished compiling the NetLib FORTRAN routines from Donald Amos for
> the Airy function. This triggered me to see if gnuplot had the Airy
> function. It does! Unfortunately:
> 
> gnuplot> print airy(0.1)
> 0.32848812610163
> 
> The correct value according to Abramowitz Stegun is
> 0.32920313
> 
> The Amos routines return more digits
> 0.329203129943538
> 
> There must be a glaring error in the gnuplot implementation for such a
> large discrepancy to occur. If I can be of help to sort out where this
> might be I will be glad to do so.
> 
> Alex van der Spek

The gnuplot implementation is documented as follows:
/* ------------------------------------------------------------
   Airy Function Ai(x)

   After:
   "Two-Point Quasi-Fractional Approximations to the Airy Function Ai(x)"
   by Pablo Martin, Ricardo Perez, Antonio L. Guerrero
   Journal of Computational Physics 99, 337-340 (1992)

   Beware of a misprint in equation (5) in this paper: The second term in
   parentheses must be multiplied by "x", as is clear from equation (3)
   and by comparison with equation (6). The implementation in this file
   uses the CORRECT formula (with the "x").

   This is not a very high accuracy approximation, but sufficient for
   plotting and similar applications. Higher accuracy formulas are
   available, but are much more complicated (typically requiring iteration).

   Added: janert (PKJ) 2009-09-05
   ------------------------------------------------------------ */

Figure 1a of the 1992 paper shows a maximum error between 0.003 and 0.012
in the region of x=0.1   So your result (delta = 0.007) is consistent 
with the accuracy claimed for the algorithm by the original authors .  

If you need a more accurate implementation and are willing to code it,
I'm sure it would be welcomed.  I suspect that might require first 
implementing a more general routine for Bessel functions; so far
gnuplot only handles J0 and J1.

	Ethan

Back to comp.graphics.apps.gnuplot | Previous | Next — Previous in thread | Next in thread | Find similar

Thread

Airy function wrong "Alex van der Spek" <zdoor@xs4all.nl> - 2012-08-09 16:25 +0200
  Re: Airy function wrong sfeam <sfeam@users.sourceforge.net> - 2012-08-09 15:31 -0700
    Re: Airy function wrong "Alex van der Spek" <zdoor@xs4all.nl> - 2012-08-11 20:41 +0200
      Re: Airy function wrong sfeam <sfeam@users.sourceforge.net> - 2012-08-14 11:29 -0700

csiph-web