Another C99 challenge

[jacob navia]

This program tests the setround/getround
primitives.

The program divides 1 by 3 with two different settings
of fesetround: upward and downard. Obviously
the result should be different in the last bit.
-----------------------------------------cut here
#include <stdio.h>
#include <math.h>
#include <fenv.h>
#include <assert.h>
#include <stdint.h>

#pragma STDC FP_CONTRACT OFF
#pragma STDC FENV_ACCESS ON
double foo(int round_dir,double d)
{
double result;
int save_round;
int setround_ok;
save_round = fegetround();
setround_ok = fesetround(round_dir);
assert(setround_ok == 0);
result = 1.0/d;
fesetround(save_round);
return result;
}

static void print_double(double d)
{
uint32_t *pint;
pint = (uint32_t *)&d;
// For little endian machines (x86)
// use 1 then zero
// For Big endian machines (power pc)
// use 0 then 1
printf("0x%x%x\n",pint[1],pint[0]);
}

int main(void)
{
double d1 = foo(FE_UPWARD,3.0);
double d2 = foo(FE_DOWNWARD,3.0);
print_double(d1);
print_double(d2);
}
------------------------------------------cut here

Output:
0x3fd5555555555556
0x3fd5555555555555

Passes under windows lcc-win, gcc, and AIX xlc

Note that in principle you should be able to put
the standard pragmas within a compound statement
and they should be local to that compound statement.

This is quite difficult to implement in real
world compilers. Neither gcc nor xlc (AIX)
have completely rewritten their code generation
phase to satisfy this crazy requirement.

lcc-win ignores those pragmas for the moment.

Besides that detail, C99 compilers support
the setting/unsetting of the floating point environment,
a bonus for complicated calculations!

--
jacob navia
jacob at jacob point remcomp point fr
logiciels/informatique
http://www.cs.virginia.edu/~lcc-win32

[Keith Thompson]

jacob navia <ja***@nospam.comwrites:
[...]

static void print_double(double d)
{
uint32_t *pint;
pint = (uint32_t *)&d;
// For little endian machines (x86)
// use 1 then zero
// For Big endian machines (power pc)
// use 0 then 1
printf("0x%x%x\n",pint[1],pint[0]);
}

[...]

This is, of course, non-portable. It assumes that double is exactly
64 bits, and it assumes that accessing half of a double object as a
uint32_t won't cause alignment problems. (It's not likely that
uint32_t will have stricter alignment requirements than a 64-bit
double, but it's permitted.)

Why not just print the floating-point value in hexadecimal using "%A"?

--
Keith Thompson (The_Other_Keith) ks***@mib.org <http://www.ghoti.net/~kst>
Nokia
"We must do something. This is something. Therefore, we must do this."
-- Antony Jay and Jonathan Lynn, "Yes Minister"

[Flash Gordon]

jacob navia wrote, On 27/08/08 23:09:

This program tests the setround/getround
primitives.

The program divides 1 by 3 with two different settings
of fesetround: upward and downard. Obviously
the result should be different in the last bit.
-----------------------------------------cut here
#include <stdio.h>
#include <math.h>
#include <fenv.h>
#include <assert.h>
#include <stdint.h>

#pragma STDC FP_CONTRACT OFF
#pragma STDC FENV_ACCESS ON
double foo(int round_dir,double d)
{
double result;
int save_round;
int setround_ok;
save_round = fegetround();
setround_ok = fesetround(round_dir);
assert(setround_ok == 0);
result = 1.0/d;

Now replace the above by:
result = 1.0/3.0;

fesetround(save_round);
return result;
}

static void print_double(double d)
{
uint32_t *pint;
pint = (uint32_t *)&d;
// For little endian machines (x86)
// use 1 then zero
// For Big endian machines (power pc)
// use 0 then 1
printf("0x%x%x\n",pint[1],pint[0]);

As Keith points out, replace the above by
prints("%A\n",d);

}

int main(void)
{
double d1 = foo(FE_UPWARD,3.0);
double d2 = foo(FE_DOWNWARD,3.0);
print_double(d1);
print_double(d2);
}
------------------------------------------cut here

Output:
0x3fd5555555555556
0x3fd5555555555555

Passes under windows lcc-win, gcc, and AIX xlc

Now explain my original question or tell us what output you get with my
modifications and explain why this output is valid for my system:

0X1.5555555555555P-2
0X1.5555555555555P-2

Note that in principle you should be able to put
the standard pragmas within a compound statement
and they should be local to that compound statement.

This is quite difficult to implement in real
world compilers.

It is also potentially quite useful in real applications. It means you
can easily get the behaviour you need for one part of the code without
affecting the rest of it.

Neither gcc nor xlc (AIX)
have completely rewritten their code generation
phase to satisfy this crazy requirement.

Crazy or not it is something where all three fail.

lcc-win ignores those pragmas for the moment.

Well, one of them can be ignored with no impact on conformance.
Specifically there is no requirement for you to use contracted
expressions so you can comply by ignoring the directive as long as you
never contract the expressions.

Note that I believe that forbidding contracted expressions means that
the compiler is required to behave even with a contracted floating point
expression as if it performed it at runtime using the selected rounding
mode. Of course, I could be wrong.

Besides that detail, C99 compilers support
the setting/unsetting of the floating point environment,
a bonus for complicated calculations!

However, I believe that they do not implement the requirements of for
not contracting expressions when it is forbidden.
--
Flash Gordon