> I don’t think the standard specifies how double is converted to single,
Correct in the default case, though IEEE_REAL will be provided in Fortran 2015 to enable explicit specification of an ISO/IEC/IEEE 60559 conformant conversion.
> and I think that the typical implementation just disregards the low-order bits of the mantissa.
No, the typical implementation uses IEEE arithmetic with unbiased round to nearest. Unless you've changed the rounding mode to "round towards zero".
Note that unlike VAX D and F format, you cannot convert IEEE double to IEEE single by just dropping the second word, as the exponent size in IEEE double is larger than that of IEEE single. Doing bit-fiddling to change the formats would be a lot more work (and is highly unlikely to be any faster) than simply using the machine instructions for conversion, e.g. on "modern" x86 the SSE2 instruction CVTSD2SS, which rounds the result according to the current rounding mode.
In other words, it is highly likely the machine is already giving you what you're asking for. You could use IEEE_GET_ROUNDING_MODE to check that you really are in round-to-nearest mode.
> My other thought is to do the summation in single precision but using Kahan. Won’t improve accuracy, but it might eliminate the bias?
Kahan single precision summation is unlikely to be significantly faster than double precision summation (promoting the single precision operands before adding), and could actually be slower (e.g. if single precision addition is less than twice as fast as double precision). The double precision version would also be simpler and more accurate. And single precision Kahan sums take the same amount of storage as plain double precision sums...
Cheers,
--
..............Malcolm Cohen, NAG Oxford/Tokyo.
|