FAQFAQ   SearchSearch   MemberlistMemberlist   UsergroupsUsergroups 
 ProfileProfile   PreferencesPreferences   Log in to check your private messagesLog in to check your private messages   Log inLog in 
Forum index » Science and Technology » Math » num-analysis
Correction factor in computing exp()?
Post new topic   Reply to topic Page 1 of 1 [3 Posts] View previous topic :: View next topic
Author Message
lcw1964
science forum beginner


Joined: 16 Jul 2006
Posts: 2

PostPosted: Sun Jul 16, 2006 7:22 am    Post subject: Correction factor in computing exp()? Reply with quote

I have recently developed an interest in Cody's approximations of
erf() and erfc(). Here is a snippet of the Fortran code from
http://www.netlib.org/specfun/erf:

YSQ = AINT(Y*SIXTEN)/SIXTEN
DEL = (Y-YSQ)*(Y+YSQ)
RESULT = EXP(-YSQ*YSQ) * EXP(-DEL) * RESULT

The point of this is to simply compute exp(-y*y) and multiply to the
result of a rational approximation calculation. Instead of doing it
straight ahead, the code multiplies y by 16, takes the integer part,
multiplies the 16 back, calls it YSQ, then computes some curious factor
DEL. Then the routine takes exp() of these separated pieces and
multiplies the whole thing together.

Does anyone know the computational point of doing this? Is or was there
some bug in the Fortran implementation of exp() that required some sort
of correction? Someone familiar with the GSL in C advises me that this
turns up in that code too. I have Cody's original 1969 paper on which
this code is based (indeed, Cody is the author of the code too), and
there is no mention of this little adjustment, or correction, or
whatever it is.

Many thanks in advance. This is driving me nuts with head
scratching....

Les
Back to top
lcw1964
science forum beginner


Joined: 16 Jul 2006
Posts: 2

PostPosted: Sun Jul 16, 2006 7:37 am    Post subject: Re: Correction factor in computing exp()? Reply with quote

lcw1964 wrote:

Quote:

Many thanks in advance. This is driving me nuts with head
scratching....


I have even just written a little code in Maple the carries out the
computation to whatever precision I like, then compares it to a direct
calculation of exp(-y*y). Gives the exact same answer for a wide range
of arguments--something I could''ve readily concluded if I had just
worked thru the basic algebra, I gather, but I am a bit slothful.

But still, I don't know why the programmer would go to this trouble.
Evidently, it is intended to improve the accuracy or precision of the
exp() calculation.

Still eager for wisdom. Fill me, please!

Les
Back to top
dave_and_darla@Juno.com
science forum beginner


Joined: 17 Oct 2005
Posts: 49

PostPosted: Sun Jul 16, 2006 3:42 pm    Post subject: Re: Correction factor in computing exp()? Reply with quote

lcw1964 wrote:
Here is a snippet of the Fortran code from
Quote:
http://www.netlib.org/specfun/erf:

YSQ = AINT(Y*SIXTEN)/SIXTEN
DEL = (Y-YSQ)*(Y+YSQ)
RESULT = EXP(-YSQ*YSQ) * EXP(-DEL) * RESULT

Does anyone know the computational point of doing this?

It enhances accuracy on computers that don't have guard digits. On such
a machine, typically, Y*Y will be accurate to about 1/2 ULP (unit in
the last place). But YSQ*YSQ+DEL should have more accuracy because
YSQ*YSQ and Y-YSQ can be computed exacly, and Y-YSQ will be smaller
than Y if abs(Y) > 1/16. Thus, the 1/2 ULP introduced in Y+YSQ will be
multiplied by a small number.

Example: Suppose the machine uses IEEE single precision and Y =
4.00xxx. Then YSQ will be 4.0 and Y-YSQ = 0.00xxx, which is less than
Y/400. Then DEL < (Y/400)*(2*Y) = Y*Y/200, and YSQ*YSQ+DEL will have 7
to 8 bits more accuracy than Y*Y. You can check this by comparing Y*Y -
YSQ*YSQ to DEL.

Dave
Back to top
Google

Back to top
Display posts from previous:   
Post new topic   Reply to topic Page 1 of 1 [3 Posts] View previous topic :: View next topic
The time now is Thu Dec 14, 2017 5:01 pm | All times are GMT
Forum index » Science and Technology » Math » num-analysis
Jump to:  

Similar Topics
Topic Author Forum Replies Last Post
No new posts Markov transition probabilities with Error correction meris300@gmail.com Probability 0 Thu Jul 20, 2006 8:22 pm
No new posts Cognition Factor....Party soon over? _Schwann_ New Theories 6 Tue Jul 18, 2006 10:53 pm
No new posts Computing pdf/cdf for X^2 when X is normally distributed Konrad Viltersten Math 8 Wed Jul 12, 2006 2:58 pm
No new posts Computing the Gerstenhaber Bracket Aaron Bergman Research 2 Fri Jul 07, 2006 2:50 am
No new posts Student Factor and degree of freedom dh Prediction 2 Thu Jun 29, 2006 3:17 pm

Copyright © 2004-2005 DeniX Solutions SRL
Other DeniX Solutions sites: Electronics forum |  Medicine forum |  Unix/Linux blog |  Unix/Linux documentation |  Unix/Linux forums  |  send newsletters
 


Powered by phpBB © 2001, 2005 phpBB Group
[ Time: 0.0207s ][ Queries: 16 (0.0036s) ][ GZIP on - Debug on ]