Maths constants
in [Fortran]
|
From: glen herrmannsfeldt on 23 Jan 2010 01:21 Tim Prince <TimothyPrince(a)sbcglobal.net> wrote: (snip, I wrote) >> x87 has the FLDPI instruction to load the internal value of PI >> in extended real (64 bit significand) into a register. It would >> be nice to have a way to do that from Fortran. (snip) > Among the problems here is that x87 hasn't been a usual mode of > operation for Fortran since SSE2 was introduced 8 years ago. If you use > an inverse trig function to get Pi, you're at the mercy of the writer of > the math library. While I understand that the compiler generates SSE2 instructions for floating point operations, is it possible that library routines, written in assembler, use x87 instructions? A quick look through glibc2 makes it look that way, but I am not completely sure yet. In addition, it seems that x87 does not have an arccosine instruction, but uses FPATAN instead. That is, one computes atan2(sqrt(1-x**2),x) -- glen
From: FX on 23 Jan 2010 09:22 > But I found it interesting that the same compiler (actually runtime > library) crashed when confronted with a nearby value: It doesn't really crash :) Minor point: GCC doesn't come with its own math library, so gfortran uses the system math library (with some low-quality fall-backs in case the system has a poor math library). -- FX
From: Tim Prince on 23 Jan 2010 09:44 glen herrmannsfeldt wrote: > > While I understand that the compiler generates SSE2 instructions > for floating point operations, is it possible that library routines, > written in assembler, use x87 instructions? A quick look through > glibc2 makes it look that way, but I am not completely sure yet. It's true, for years, if your compiler didn't provide its own math library for SSE2, you were likely calling on an x87 library. Part of the point of this entire thread is that you can't guess whether this kind of thing will happen. With glibc, you must actually trace your execution if you want to find out what is happening in math libraries. There are still games being played with math libraries for glibc. > > In addition, it seems that x87 does not have an arccosine instruction, > but uses FPATAN instead. That is, one computes > > atan2(sqrt(1-x**2),x) > That was a common way of doing it with x87 instructions. It didn't work so well for double precision when 53-bit precision mode was set, as Microsoft practice requires. atan2(sqrt((1+x)*(1-x)),x) is a little better.
From: David Thompson on 28 Jan 2010 11:41 On Wed, 20 Jan 2010 20:10:38 +0100, Giorgio Pastore <pastgio(a)units.it> wrote: > Les Neilson wrote: > ... > > The first difficulty would be specifying a *finite* list of constants > > (there are rather a lot across the sciences that I can think of and > > presumably plenty more of which I am unaware or have forgotten - > > someone's favourite constant being left out might cause a diplomatic > > incident!) > > But, why such a difficulty should be more serious for Fortran community > than for C-speaking people ? > I would guess that at least the issue of compatibility with C should > prompt a more careful consideration of the possibility of including > mathematical constants in the language. > Standard C does not include math/science constants, although Unix (POSIX/SUS) does. See 14.8 in comp.lang.c FAQ at the usual places and http://c-faq.com . C does have 'machine' constants like INT_MAX, and 'programming' ones like BUFSIZ and EDOM.
From: glen herrmannsfeldt on 28 Jan 2010 15:10
David Thompson <dave.thompson2(a)verizon.net> wrote: > On Wed, 20 Jan 2010 20:10:38 +0100, Giorgio Pastore wrote: (snip) >> But, why such a difficulty should be more serious for Fortran community >> than for C-speaking people ? (snip) > Standard C does not include math/science constants, although Unix > (POSIX/SUS) does. See 14.8 in comp.lang.c FAQ at the usual places and > http://c-faq.com . C does have 'machine' constants like INT_MAX, and > 'programming' ones like BUFSIZ and EDOM. I suppose that"C-speaking people" is general enough to include POSIX. (Presuming POSIX was designed by C-speaking people and not Fortran speaking people.) Java has a standard system class with some math constants including pi. Not strictly part of the language, but close enough for this discussion. It seems, though, that the most common use of pi is with trig. functions which then (usually) divide by pi as part of the argument reduction process. A sind() function in degrees, or some other integer multiple of the period, (maybe octant or quadrant) would help in that case. One of the features of PL/I that still hasn't been adopted by Fortran... -- glen |