Date: Thu, 24 Oct 2002 22:54:35 +1000 (EST) From: Bruce Evans <bde@zeta.org.au> To: Loren James Rittle <rittle@latour.rsch.comm.mot.com> Cc: current@FreeBSD.ORG Subject: Re: Lack of real long double support (was Re: libstdc++ does not contain fabsl symbol) Message-ID: <20021024205944.R1451-100000@gamplex.bde.org> In-Reply-To: <200210240833.g9O8XB1W033884@latour.rsch.comm.mot.com>
next in thread | previous in thread | raw e-mail | index | archive | help
On Thu, 24 Oct 2002, Loren James Rittle wrote:
> ... Anyways, that work exposed some issues.
>
> We have this in the system header:
>
> #define LDBL_MANT_DIG DBL_MANT_DIG
> #define LDBL_MIN_EXP DBL_MIN_EXP
> #define LDBL_MAX_EXP DBL_MAX_EXP
> [...]
This seems to be correct. Long double precision is not really supported
either at complie time or runtime. The default precision (on i386's)
is 53 bits, so (normalized) long doubles have a hard time getting any
larger than DBL_MAX or any smaller than DBL_MIN (you can create them
as bits or using a hardware transcendental function, but any arithmetic
on them will convert them to double precision).
> Yet, the FP hardware is actually configured by default to provide
> `long double' as:
>
> #define LDBL_MANT_DIG 53
I think you mean 64 (the hardware default).
> #define LDBL_MIN_EXP (-16381)
> #define LDBL_MAX_EXP 16384
>
> Indeed, FP code using long double generated by gcc 2.95.X, installed
> as the FreeBSD 4 system compiler, uses the full exponent range of
> -16381 to 16384. However, printf(), etc loses on that exponent range
> and reports Inf. I suspect this is why the system header misreports
> the FP hardware, thus the header describes the largest printable long
> double value, not the largest that could be held in a calculation.
gcc can't actually support the full range, since it doesn't control
the runtime environement (it could issue a fninit before main() to
change the default, but it shouldn't and doesn't). The exponent
range is lost long before printf() is reached. E.g.,
long double x= DBL_MAX;
long double y = 2 * x;
gives +Inf for y since the result is doesn't fit in 53-bit precision.
The system header correctly reports this default precision. Any header
genrated by the gcc build should be no different, since the build should
run in the target environment.
> Anyways, two questions for FreeBSD maintainers. How should gcc, as
> provided from the FSF, describe the long double FP format for
> FreeBSD/i386 4.x? Shall we assume that no changes for FreeBSD 5.x
> will occur?
It should use whatever is the default format for the host environment,
It still has enquire.c for doing this automatically. I don't know
when enquire.c actually gets used. The FreeBSD build just ignores it
and uses src/sys/${MACHINE_ARCH}/float.h. My first encounter with this
bug suite was near enquire.c 12-14 years ago. enquire gave different
results for -msoft-float and hardware FP because my soft-float libraries
only had 53-bit precision but hardware FP used 64-bit precision.
> ...
> However, that also brings to mind: Are we happy with the current
> default FP hardware setup especially if support for long double is to
> be improved? At the moment we support long double in name only.
> There is no extra dynamic range over a double; as might be implied
> (and allowed by the hardware).
Not really. Assembly is still required to get full control over precision.
I'm still waiting for (C) language support to catch up. Been waiting for
about 14 years so far. C99 clarified the semantics of floating point
precision but is not supported by gcc (yet?).
> Also, as an aside, having just learned everything there is to know
> about IEEE FP from a witty 80-slide presentation on the WWW by one of
> the original designers of the spec@Berkeley, I'd say that he would be
> sad that we default to use 53-bit precision mode. I'd say that it is
> dumb to claim we even support long double unless we change that to
> 64-bit precision mode...
I think he would be even unhappier with gcc's support for it. We default
to 53-bit precision partly because Kahan's (his?) paranioa(1) is unhappy
with 64-bit precision.
> Regarding the comment in sys/i386/include/npx.h:
>
> * 64-bit precision often gives bad results with high level languages
> * because it makes the results of calculations depend on whether
> * intermediate values are stored in memory or in FPU registers.
>
> This is pure bunk when considered in broad context. Using 53-bit
> precision mode with high level languages while making actual
> calculations yields more outright undetectable precision-related
> errors due to algorithm design by non-numerical analysts. I know
> which error I'd rather *not* find in my naively-written and compiled
> FP code.
I think it makes no provably significant difference for naively-written
code. I'd rather get the same (perhaps wrong) results on all machines
from naively-written code.
> People that write their FP code to correctly use epsilon (the one
> appropriate for the type of the float they are using) should never see
> this problem, no?
No. They would see strange behaviour cause by compiler bugs.
> I don't see why FreeBSD should cater to people that
> don't know how to write FP code properly (i.e. anyone that expected
> exact results across compilation switches, etc) at the expense of
> being able to write code that fully utilizes the best mode available
> from the CPU.
IEEE754 provides many guarantees about exact results. FreeBSD caters
mainly to people that _do_ know how to write low level FP code (i.e.,
anyone that knows when to expect exact results) but who don't know about
compiler bugs.
The (i386 gcc) compiler bugs are mainly that in:
double x, y, z, t;
...
t = (double)(x + y) + z;
`t' may have long double precision despite the explicit assignment to it,
and the intermediate value (double)(x + y) may have long double precision
despite its explict cast to double.
Bruce
To Unsubscribe: send mail to majordomo@FreeBSD.org
with "unsubscribe freebsd-current" in the body of the message
Want to link to this message? Use this URL: <https://mail-archive.FreeBSD.org/cgi/mid.cgi?20021024205944.R1451-100000>