2

Different platforms have varying FP capabilities with varying parameters and behaviors, as a result there is a degree of variance between the calculation results they produce, which cascade and amplify on each intermediate step.

I am in a situation where it is critical for (+-*/ only) calculations to produce identical results on each and every different target platform, using different compiler vendors, so I wonder if there is a standard way to do that. I am not asking about arbitrary high precision floating point numbers but standard 64 bit IEEE double, and a performance hit is expected and tolerable.

dtech
  • 47,916
  • 17
  • 112
  • 190
  • You could do the floating point calculation in software using integers. – eerorika Mar 07 '19 at 14:24
  • Have you tried a fixed point library? I think that should give the same results on all implementations. – NathanOliver Mar 07 '19 at 14:25
  • I was hoping there might be some way to force software computation for a compilation unit... – dtech Mar 07 '19 at 14:27
  • Floating point arithmetic is mostly inaccurate anyway, and there is not overall standard of compliance by C. You could lower your sights to a specific accuracy and have a reasonable expectation of identical results. Before calculators were available that was what was done: the accuracy was specified, for example an engineer would work to the best accuracy that he could and present the result to, say, 4 significant figures. – Weather Vane Mar 07 '19 at 14:30
  • @WeatherVane reproducibility takes precedence to accuracy in this usage scenario – dtech Mar 07 '19 at 14:31
  • Even with fixed point arithmetic, you will need to specify the rounding method, for intermediate values too, but if you supply the library with the software it should produce identical results on different platforms. – Weather Vane Mar 07 '19 at 14:33
  • 1
    I'm starting to think this is too broad. What set of operations are you performing on the floating point types? – Bathsheba Mar 07 '19 at 14:36
  • @Bathsheba +-*/ only – dtech Mar 07 '19 at 14:38
  • Well the performance hit won't be as bad as with math functions. – Weather Vane Mar 07 '19 at 14:39
  • 2
    gcc uses the mpfr library for compile time floating point computations. Might look into it as an option. – Shawn Mar 07 '19 at 14:41
  • @dtech It will be possible then, you need to check your compiler documentation to make sure it follows IEEE754 faithfully. E.g. with MSVC, play around here: https://learn.microsoft.com/en-us/cpp/build/reference/fp-specify-floating-point-behavior?view=vs-2017 – Bathsheba Mar 07 '19 at 14:45
  • 1
    All major compilers have options which affects floating point code generation. If you use a strict model, then the basic operations (+, -, *, /, sqrt) should give the same result. The only wide-spread exception is x86 fpu, which uses a 80-bit precision by default. But, this can be changed to double precision as well (plus, gcc has the option `-ffloat-store`. Or just use sse instead of fpu.). – geza Mar 07 '19 at 14:54

1 Answers1

0

Even if you have a 64 bit IEEE754 double, there are a few extra things you need to check.

  1. Make sure you have strict floating point. Don't allow your compiler to use, for example, 80 bits for intermediate calculations.

  2. Various operations (all the arithmetic operations such as the ones you mention, std::sqrt, &c.) are required by IEEE754 to return the best number possible. (Should you need others then make sure that all your operations are mentioned in the IEEE754 standard and your platform obeys that faithfully - it might not).

  3. Shy away from other functions (such as trigonometric functions), for which there is no guarantee of precision, even under IEEE754.

In your specific case it appears that (1) is sufficient, along with perhaps (for C++)

static_assert(std::numeric_limits<double>::is_iec559, "IEEE 754 floating point required");
Bathsheba
  • 231,907
  • 34
  • 361
  • 483
  • C compilers are allowed to default to `#pragma STDC FP_CONTRACT on`, allowing them to contract `x*y + z` into an FMA with no rounding of the `x*y` intermediate result. AFAIK, the same thing applies in C++. GCC does this even when the mul and add aren't part of the same expression, but in separate statements, in violation of ISO C. [Is floating point expression contraction allowed in C++?](//stackoverflow.com/q/49278125) indicates that this applies to C++. [How to use Fused Multiply-Add (FMA) instructions with SSE/AVX](//stackoverflow.com/a/34461738) has more – Peter Cordes Mar 07 '19 at 23:18
  • But once you take care of that, and avoid 80-bit x87 by checking [`FLT_EVAL_METHOD`](http://en.cppreference.com/w/c/types/limits/FLT_EVAL_METHOD) to rule out high-precision temporaries, you might be ok for the basic ops that are required to produce correctly-rounded results. See also [Does any floating point-intensive code produce bit-exact results in any x86-based architecture?](//stackoverflow.com/q/27149894) for the stuff I dug up answering a similar question. – Peter Cordes Mar 07 '19 at 23:22