Home > Archive > Fortran > June 2005 > Fortran functions
You are viewing an archived Text-only version of the thread.
To view this thread in it's original format and/or if you want to reply to
this thread please [click here]
|
|
| Dieter Britz 2005-06-06, 3:58 pm |
| I am wondering how the inbuilt functions in Fortran are
evaluated, e.g. SQRT, SIN, EXP, LOG. Can anyone here tell me?
--
Dieter Britz, Kemisk Institut, Aarhus Universitet, Danmark.
| |
| Arjen Markus 2005-06-06, 3:58 pm |
| Dieter Britz wrote:
>
> I am wondering how the inbuilt functions in Fortran are
> evaluated, e.g. SQRT, SIN, EXP, LOG. Can anyone here tell me?
> --
> Dieter Britz, Kemisk Institut, Aarhus Universitet, Danmark.
That will depend on the compiler. I remember once seeing an
IBM manual on the way their VM FORTRAN compiler did those
computations. I think I have seen a similar document from
SUN. With g77 and gfortran and ... you can actually
inspect the source code.
But perhaps you want to see guarantees about accuracies
and such?
Regards,
Arjen
| |
| David Frank 2005-06-06, 3:58 pm |
|
"Dieter Britz" <britz@chem.au.dk> wrote in message
news:d81cst$df7$1@news.net.uni-c.dk...
>I am wondering how the inbuilt functions in Fortran are
> evaluated, e.g. SQRT, SIN, EXP, LOG. Can anyone here tell me?
> --
> Dieter Britz, Kemisk Institut, Aarhus Universitet, Danmark.
>
The Intel CPU has a superior 80-bit FPU with on-board instructions that can
speed up those functions
IF the Fortran compiler has option to use them, CVF compiler (perhaps
others) has /fast switch that significantly speeds up the evaluations. My
own testing years ago satisfied me that there is NO loss of accuracy for
64-bit results using
the on-board FPU instructions.
| |
| Dieter Britz 2005-06-06, 3:58 pm |
| Arjen Markus wrote:
> Dieter Britz wrote:
>
>
>
> That will depend on the compiler. I remember once seeing an
> IBM manual on the way their VM FORTRAN compiler did those
> computations. I think I have seen a similar document from
> SUN. With g77 and gfortran and ... you can actually
> inspect the source code.
>
> But perhaps you want to see guarantees about accuracies
> and such?
No, although that might be interesting to know. I want to know
rather generally how they do it. There must be preferred methods
for each function. For example, I could do SQRT using Newton's
method, but I doubt that that is the way most compiler writers
would do it. Are there, despite the individual compiler preferences,
any general methods that tend to be used for these functions by most?
--
Dieter Britz, Kemisk Institut, Aarhus Universitet, Danmark.
| |
| Tim Prince 2005-06-06, 3:58 pm |
|
"Dieter Britz" <britz@chem.au.dk> wrote in message
news:d81fpb$edn$1@news.net.uni-c.dk...
> Arjen Markus wrote:
>
> No, although that might be interesting to know. I want to know
> rather generally how they do it. There must be preferred methods
> for each function. For example, I could do SQRT using Newton's
> method, but I doubt that that is the way most compiler writers
> would do it. Are there, despite the individual compiler preferences,
> any general methods that tend to be used for these functions by most?
Variations on Newton's method for sqrt() are generated in line only for
machines which have specific instruction level support for that. On certain
compilers there are flags to switch among methods for implementations of
sqrt(), which tends IMHO to show lack of generality of any method.
Compilers I can think of at the moment provide at least an option to make
sqrt() comply with IEEE standard http://grouper.ieee.org/groups/754/.
Many compilers share library implementation with C and other compilers.
Thus, you might be interested in texbooks such as
The Standard C. Library
Author: P J Plauger
ISBN: 0131315099 / Paperback
Publisher: Prentice Hall PTR / 1991-06
List Price: £31.95
| |
| Charles Russell 2005-06-06, 3:58 pm |
| Dieter Britz wrote:
I want to know[color=darkred]
> rather generally how they do it. There must be preferred methods
> for each function.
A standard reference has been: William J. Cody, Jr. and William Waite,
"Software Manual for the Elementary Fuctions," Prentice-Hall 1980. But
perhaps these old methods have been superseded.
| |
| James Van Buskirk 2005-06-06, 3:58 pm |
| "David Frank" <dave_frank@hotmail.com> wrote in message
news:FkXoe.26$pa3.19@newsread2.news.atl.earthlink.net...
> "Dieter Britz" <britz@chem.au.dk> wrote in message
> news:d81cst$df7$1@news.net.uni-c.dk...
[color=darkred]
> The Intel CPU has a superior 80-bit FPU with on-board instructions that
can
> speed up those functions
> IF the Fortran compiler has option to use them, CVF compiler (perhaps
> others) has /fast switch that significantly speeds up the evaluations. My
> own testing years ago satisfied me that there is NO loss of accuracy for
> 64-bit results using
> the on-board FPU instructions.
If you are curious about how the Intel CPU implements these
functions, I think the most relevant references may be:
P. T. P. Tang, "Table-driven implementation of the exponential
function in IEEE floating-point arithmetic," ACM Transactions
on Mathematical Software, vol. 16, no. 2, June 1989, pp. 144-157.
P. T. P. Tang, "Table-driven implementation of the logarithm
function in IEEE floating-point arithmetic," ACM Transactions
on Mathematical Software, vol. 16, no. 2, December 1990, pp.
378-400.
Hmmm... I lifted these references out of the Pentium Processor
User's Manual, Volume 3: Architecture and Programming Manual,
Appendix G: Report on Transcendental Functions, but there seems
to be a mistake there because the volume and issue numbers
are the same, but the months and years are different. I would
tend to believe the months and years in the above.
--
write(*,*) transfer((/17.392111325966148d0,6.5794487871554595D-85, &
6.0134700243160014d-154/),(/'x'/)); end
| |
| glen herrmannsfeldt 2005-06-06, 3:58 pm |
| Dieter Britz wrote:
(snip)
(snip)
[color=darkred]
> No, although that might be interesting to know. I want to know
> rather generally how they do it. There must be preferred methods
> for each function. For example, I could do SQRT using Newton's
> method, but I doubt that that is the way most compiler writers
> would do it. Are there, despite the individual compiler preferences,
> any general methods that tend to be used for these functions by most?
Software sqrt on processors with floating point does usually use
Newton's method. The trick is doing a good first approximation.
First thing is to divide the exponent by two, and remember if it
was odd. Then usually a simple approximation, sometimes a linear
approximation, is enough. Two cycles of Newton's method for single
precision, four for double.
For log, the exponent is separated and used directly, with a
series expansion on the fraction. The reverse for exp.
For sin and cos, first argument reduction is done, then some
identities to further reduce the argument, followed by a
polynomial expansion.
Traditionally they were done in assembly, with easy access to
the individual bits. More recently they are usually done in
a high level language, using functions to access the parts of
the floating point value. Such functions are best written in assembly,
but more portable in a high level language.
-- glen
| |
| Bart Vandewoestyne 2005-06-06, 3:58 pm |
| In article <d81cst$df7$1@news.net.uni-c.dk>, Dieter Britz wrote:
> I am wondering how the inbuilt functions in Fortran are
> evaluated, e.g. SQRT, SIN, EXP, LOG. Can anyone here tell me?
Absolutely coincidence, but today I have been looking at the erf(x)
description in the following book:
J.F. Hart, E.W. Cheney et al, Computer Approximation, John Wiley, N.Y.
(1968)
I looked at it because somewhere in some code from Allan Miller for erf(x),
this book was referred to.
I don't know in how far this book from 1968 is superseded by newer
and better methods though...
Regards,
Bart
--
"Share what you know. Learn what you don't."
| |
| robert.corbett@sun.com 2005-06-07, 3:57 am |
| > Variations on Newton's method for sqrt() are generated in line only for
> machines which have specific instruction level support for that.
No. I used to use a compiler that inlined all of the simple elementary
functions: SQRT, SIN, COS, TAN, EXP and ALOG. I know that the
SQRT used a Newton-Raphson algorithm, not a machine instruction.
The other functions used tables that came from a library, but the code
was all inlined.
Bob Corbett
| |
| robert.corbett@sun.com 2005-06-07, 3:57 am |
| Peter Tang' papers are excellent, but they do not describe
how Intel CPUs implement log and exp. The correct
volume and issue numbers are TOMS 15, 2 and TOMS 16,4.
The dates given are correct.
BTW, don't lose your Perntium manual. It has the best description
of Intel's pre-SSE floating-point arithmetic of any of Intel's manuals.
Bob Corbett
| |
| Tim Prince 2005-06-07, 4:02 pm |
|
"James Van Buskirk" <not_valid@comcast.net> wrote in message
news:ft2dnbG3tLC06TnfRVn-1g@comcast.com...
> "David Frank" <dave_frank@hotmail.com> wrote in message
> news:FkXoe.26$pa3.19@newsread2.news.atl.earthlink.net...
>
>
>
> can
>
I am wondering why no one pointed out that CVF initialized precision mode to
53 bits, so that sqrt() (but not other built-in transcendentals) would be
rounded the same as the current SSE2 sqrt(). I never saw an application
which employed the CVF library calls to set precision back up to 64 bits. I
can't tell whether Dave means 53-bit precision when he says "64-bit
results," but 53-bit precision mode generally eliminates the advantage of
80-bit registers.
Beginning with Pentium III, the internal math functions generally carried
64-bit precision, but implementation of exponentiation by the simplest
method could easily expose inaccuracy in range reduction, as the Intel
documents point out.
| |
| glen herrmannsfeldt 2005-06-07, 4:02 pm |
| Tim Prince wrote:
(snip)
> I am wondering why no one pointed out that CVF initialized precision mode to
> 53 bits, so that sqrt() (but not other built-in transcendentals) would be
> rounded the same as the current SSE2 sqrt(). I never saw an application
> which employed the CVF library calls to set precision back up to 64 bits. I
> can't tell whether Dave means 53-bit precision when he says "64-bit
> results," but 53-bit precision mode generally eliminates the advantage of
> 80-bit registers.
The 80 bit registers with 64 bit fractions sometimes help and sometimes
hurt. There have been a number of posts here where the wrong answers
resulted from some calculations being done with 64 bit fractions being
compared to or subtracted from ones rounded (through intermediate
stores) to 53 bits. On average, I am not sure which is better.
> Beginning with Pentium III, the internal math functions generally carried
> 64-bit precision, but implementation of exponentiation by the simplest
> method could easily expose inaccuracy in range reduction, as the Intel
> documents point out.
Also, you say the CVF sets 53 bit mode, but does it, or the OS, restore
the previous mode when the program is finished? If not, other programs
might unexpectedly be run in 53 bit mode.
-- glen
| |
| Tim Prince 2005-06-08, 4:00 am |
|
"glen herrmannsfeldt" <gah@ugcs.caltech.edu> wrote in message
news:e9Cdne3RgvwiTjjfRVn-hQ@comcast.com...
>
> Also, you say the CVF sets 53 bit mode, but does it, or the OS, restore
> the previous mode when the program is finished? If not, other programs
> might unexpectedly be run in 53 bit mode.
If you follow a CVF program by a gfortran program, or even run them at the
same time, causing them to swap, you will see that the precision mode for
gfortran has returned to the default 64-bit (full 80-bit register). 53-bit
mode has to be set for each process. Programs generated by compilers which
use the 53-bit mode as a default have the code to switch the mode buried
either in the beginning of the main program or in a library initialization.
It is possible, of course, for the OS to set the floating point modes when
launching each program. I believe some versions of BSD set 53-bit mode
before starting a program. Recent versions of Windows are defined to set
SSE abrupt underflow mode prior to starting the program, so you would expect
a Windows gfortran program to inherit abrupt underflow regardless of whether
it is compiled in.
| |
| Arjen Markus 2005-06-08, 4:00 pm |
| Dieter Britz wrote:
>
> No, although that might be interesting to know. I want to know
> rather generally how they do it. There must be preferred methods
> for each function. For example, I could do SQRT using Newton's
> method, but I doubt that that is the way most compiler writers
> would do it. Are there, despite the individual compiler preferences,
> any general methods that tend to be used for these functions by most?
>
>
What about the publications on the evaluation of these functions
with arbitrary precision - Bailey's work comes to mind. These
publications generally have one or two things to say on efficient
and accurate evaluations. And they emphasize performance too.
Regards,
Arjen
|
|
|
|
|