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

double (*volatile explicit_ldexp)(double, int) = ldexp;

void
test(unsigned n, double expected, uint64_t e, uint64_t u)
{
	int e0, e1;
	double s0, s0a, s1a, s2a, s0b, s1b, s2b, s;
	const char *msg;

	e0 = e - 64;
	e1 = 64;
	assert(e == e0 + e1);

	u |= 0x8000000000000001ULL;
	s0 = (double)u;

	s0a = s0;
	s1a = s0a * explicit_ldexp(1, -64);
	s2a = s1a * explicit_ldexp(1, -e);

	s0b = s0;
	s1b = explicit_ldexp(s0b, -64);
	s2b = explicit_ldexp(s1b, -e);

	s = s2a;
	msg = (s == expected ? "PASS" : "FAIL");

	printf("%u %s %a = %a (a=%a b=%a)\n", n, msg, expected, s, s2a, s2b);
}

int
main(void)
{

	test(0, 4.9406564584124654e-324, 1074, 0x0000000000000400ULL);
	test(1, 4.9406564584124654e-324, 1074, 0x8000000000000400ULL);
	test(2, 9.8813129168249309e-324, 1073, 0x3fffffffffffffffULL);
	test(3, 9.8813129168249309e-324, 1073, 0xbfffffffffffffffULL);
	test(4, 2.2250738585072014e-308, 1022, 0x7ffffffffffff400ULL);

	return 0;
}