optimizer-proof minus zero generation?

[David G. Hough on validgh]

/*
The following program provides 
ways of generating IEEE minus zero, of which the last two 
should resist most optimizers.
One uses Fortran equivalence, complicated slightly by the uncertainty of
size of long vs. double, and another uses a run time computation to
produce an infinity, which should not be optimizable but could run into 
trouble at run time 
on systems that incorrectly enable trapping by default on IEEE
overflow exceptions.    Since what follows is a main program, it may rely
on default IEEE trapping modes; if converted to a library
subroutine, that factors would have to be considered.    This program 
seems to work as I intend with SPARCompilers C 2.0.1 on SunOS 4.1.3, 
unless compiled with "-fnonstd" either directly or via the macro "-fast";
those options elicit non-IEEE-standard behavior in order to be less shocking
to programs coming from VAX/VMS.

The only purpose of this program is to show that minus zeros can be
generated in a way that is fairly portable among current IEEE implementations,
even among compilers that optimize agressively or have odd ideas about
compile-time expression evaluation.
This program doesn't
do anything useful on non-IEEE systems, and won't serve any purpose
in the not-too-distant future when compile-time expression evaluation is
better regulated for NCEG-compliant compilers.
*/

#include <stdio.h>
#include <math.h>

main()
{
	double          bigger, big, diff;
	typedef union {
		double          d;
		long            i[2];
	}               equivalence;
	equivalence     p1, m1, m0;

	m0.d = -0.0;		/* how simple can you get */
	(void) printf(" minus zero %g %X %X \n", m0.d, m0.i[0], m0.i[1]);

	m0.d = -fabs(0.0);	/* should work if optimizer doesn't know
				 * about fabs */
	(void) printf(" minus zero %g %X %X \n", m0.d, m0.i[0], m0.i[1]);

	p1.d = 1;
	m1.d = -1;

	/* Get sign bit - big or little endian. */
	m1.i[0] ^= p1.i[0];
	m1.i[1] ^= p1.i[1];

	/* Apply sign bit to +0. */
	m0.d = 0;
	m0.i[0] ^= m1.i[0];
	m0.i[1] ^= m1.i[1];

	(void) printf(" minus zero %g %X %X \n", m0.d, m0.i[0], m0.i[1]);

	/*
	 * Following won't work if floating-point overflow or invalid trap
	 * enabled.
	 */
	bigger = 1e300;
	do {
		big = bigger;
		bigger = big * big;
		diff = bigger - big;	/* diff will be NaN normally, +-0 in
					 * RZ or RN mode. */
	}
	while (bigger != big);

	m0.d = -1 / bigger;	/* -0 from -1/+inf usually... */
	if (m0.d != 0)
		m0.d = -diff * diff;
	/* -0 in a way no optimizer can guess for RN/RZ modes. */
	(void) printf(" minus zero %g %X %X \n", m0.d, m0.i[0], m0.i[1]);

	return 0;
}