Reproducible floating-point results

1. Abstract

The C++ standard, hereafter referred to as [C++], has very little to say about floating-point specification, which means that it is extremely hard to write floating-point expressions that, given the same inputs and rounding mode, will produce the same results on all implementations. This paper proposes the introduction of a new type which specifies sufficient conformance with the provisions specified in ISO/IEC 60559:2020, hereafter referred to as [559], to guarantee this reproducibility.

2. Motivation

From [559] Clause 11:

Programs that can be reliably translated into an explicit or implicit sequence of reproducible operations on reproducible formats produce reproducible results. That is, the same numerical and reproducible status flag results are produced.

Reproducible results require cooperation from language standards, language processors and users. A language standard should support reproducible programming.

C++ does not support reproducible programming. In fact, the standard does not specify floating point behaviour at all. In [C++], [basic.fundamental] paragraph 12:

Except as specified in 6.8.3, the object and value representations and accuracy of operations of floating-point types are implementation-defined.

Even in 6.8.3 [basic.extended.fp], we read:

If the implementation supports an extended floating-point type [basic.fundamental] whose properties are specified by the ISO/IEC 60559 floating-point interchange format binary16, then the typedef-name std​::​float16_t is defined in the header <stdfloat> and names such a type, the macro __STDCPP_FLOAT16_T__ is defined [cpp.predefined], and the floating-point literal suffixes f16 and F16 are supported [lex.fcon].

The subsequent paragraphs continue with the other binary interchange formats specified in Clause 3 of [559]. All that is specified is that the interchange format is supported. Behaviour is not specified and thus remains implementation defined.

C++ programmers might expect that any code they write shall behave identically on any platform, in much the same way that 2 + 2 will always equal 4. However, the very broad latitude granted here means that this is not necessarily the case, and that floating-point calculations are not necessarily reproducible on all platforms.

Worse than that, even the same function within the same executable compiled with the same compiler flags but called in different contexts can produce subtly different results due to different optimisations after inlining and due to link time optimisation. Even moving a function out of a header file to a source file is not guaranteed to help matters.

Several domains would benefit from floating-point reproducibility.

Game development

The video game industry is worth over $340bn worldwide [VGI]. Games typically provide a simulation which players interact with, sometimes from different machines over the internet. As each player interacts with the simulation, those interactions must be reflected on all machines so that the players continue to interact with the same simulation.

Games can be deployed to multiple platforms, typically to mobile phones, consoles, desktop and laptop computers. These platforms all enjoy different tool chains. Since simulations are specified using floating-point arithmetic, this prevents players on different platforms playing together unless the simulation is extremely basic, and the entire simulation can be transmitted trivially across the internet to each player. The divergence in floating-point results across implementations means that the simulation cannot remain synchronised across different platforms.

Reproducibility eliminates this divergence, opening up the possibility of all players being able to compete with one another regardless of their hardware platform.

Scientific computing

Since behaviour is not specified, implementations are free to fuse and reorder instructions to benefit execution speed. Sometimes this can have catastrophic results. From [GOL] 1.4:

Catastrophic cancellation occurs when the operands are subject to rounding errors. For example, in the quadratic formula, the expression b² – 4ac occurs. The quantities b² and 4ac are subject to rounding errors since they are the results of floating-point multiplications. Suppose they are rounded to the nearest floating-point number and so are accurate to within 1/2 ulp. When they are subtracted, cancellation can cause many of the accurate digits to disappear, leaving behind mainly digits contaminated by rounding error. Hence the difference might have an error of many ulps.

The paper then goes on to describe some strategies for mitigating this problem, involving rearrangement of operations. However, since there is no guarantee in the standard that the order of evaluation will be honoured, such rearrangement may be futile and catastrophic cancellation cannot be avoided.

Reproducibility requires the order of evaluation to be specified unambiguously. As a result, catastrophic cancellation can be avoided.

3. Status quo

The obvious question is “why are floating-point operations not reproducible on different implementations?”

Consider the following example available on [CEX]:

int main() {
	float line_0x(0.f);
	float line_0y(7.f);

	float normal_x(0.57f);
	float normal_y(0.8f);

	float px(10.f);
	float py(0.f);

	float v2p_x = px - line_0x;
	float v2p_y = py - line_0y;

	float distance = v2p_x * normal_x + v2p_y * normal_y;

	float direction_x = normal_x * distance;
	float direction_y = normal_y * distance;

	float proj_vector_x = px - direction_x;
	float proj_vector_y = py - direction_y;

	std::cout << "distance:      " << distance << std::endl;
	std::cout << "proj_vector_y: " << proj_vector_y << std::endl;
}

When built with version 24.5 of the Nvidia C++ compiler, without specifying any compiler options, the output is:

distance:  	   0.0999997
proj_vector_y: -0.0799998

When built with version 18.1.0 of the clang C++ compiler, without specifying any compiler options, the output is:

distance:  	   0.0999999
proj_vector_y: -0.0799999

Worse, if -march=skylake is passed to the clang C++ compiler, the output is:

distance:  	   0.1
proj_vector_y: -0.08

Analysing the memory of each running program reveals that, indeed, the bit patterns vary. Application of optimisations leads to variations between implementations in the nature, number and order of floating-point operations being executed during the program. This leads to different numbers of rounding operations, in different orders, and thus different values on different platforms.

4. Prior art

Before considering an approach, we should consider the approaches taken by the other languages standardised through ISO/IEC.

ISO/IEC 8652:2023 Ada

The Ada standard does not normatively reference any revision of ISO/IEC 60559. It defines its own floating-point format and operations. However, there is an expectation that some future version of Ada will indeed support a revision of ISO/IEC 60559.

ISO/IEC 13211-1:1995/Cor 3:2017 Prolog

The Prolog standard makes no reference to any revision of ISO/IEC 60559 whatsoever.

ISO/IEC 1989:2023 COBOL

The COBOL standard specifies language support for the facilities defined by ISO/IEC 60559:2011 in several paragraphs throughout sections such as 8.8.1 Arithmetic expressions, 8.8.4 Conditional expressions, 14.6.8 Alignment and transfer of data into data items and in Annex D.17 Forms of arithmetic. Typically, the specification is similar to the form:

All additions, subtractions, multiplications and divisions performed in the development of the result shall be performed in accordance with the corresponding rules in ISO/IEC/IEEE 60559.

Particularly, the specification directs implementers to follow ISO/IEC 60559:2011. There is no option to avoid supporting these facilities; all arithmetic support is specified with reference to ISO/IEC 60559:2011.

ISO/IEC 1539-1:2023 Fortran

The FORTRAN standard specifies, in Clause 17, language support for the facilities defined by [559]. From the first paragraph of that clause:

The intrinsic modules IEEE_EXCEPTIONS, IEEE_ARITHMETIC, and IEEE_FEATURES provide support for the facilities defined by ISO/IEC 60559:2020. Whether the modules are provided is processor dependent. If the module IEEE_FEATURES is provided, which of the named constants defined in this document are included is processor dependent.

The decision to support [559] facilities is not up to the implementer. Their only constraint is whether the processor supports [559] facilities. Again, all operations are specified with reference to [559].

ISO/IEC 9899:2018 C

In its annexes F though H, the C standard specifies language support for the facilities defined by [559]. From Annex F, Introduction, para 1:

This annex specifies C language support for the IEC 60559 floating-point standard. The IEC 60559 floating-point standard is specifically Floating-point arithmetic (ISO/IEC 60559:2020), also designated as IEEE Standard for Floating-Point Arithmetic (IEEE 754–2019). IEC 60559 generally refers to the floating-point standard, as in IEC 60559 operation, IEC 60559 format, etc.

It does not mandate support for floating-point arithmetic as defined by [559], and implementations are free to ignore that standard when implementing floating-point arithmetic. The annexes explicitly describe how an implementation SHOULD support [559] floating point arithmetic, not how an implementation SHALL support it.

ISO/IEC DIS 14882:2023 C++

[C++] barely references [559]. In fact, only when defining the optional extended floating point types in section 6.8.3 [basic.extended.fp] is there any meaningful reference. The type trait std::is_iec559 is true for those types, and signals that those types offer some representations from that specification.

5. Possible approaches

There are several ways in which reproducible floating-point operations could be achieved.

Compiler options

One approach is to look at what options compilers already offer for reproducibility. Surveying clang, GCC and MSVC, there is no explicit support for reproducibility. However, they all offer opportunities for controlling the way that floating-point maths arithmetic is implemented.

For example, MSVC offers a number of command line flags prefixed /fp. Of note are two options, /fp:fast and /fp:strict. /fp:fast allows the compiler to reorder, combine or simplify floating-point operations to optimise floating-point code for speed and space. /fp:strict preserves the source ordering and rounding properties of floating-point code when it generates object code, and observes [559] specifications when handling special values. Similarly, GCC offers -fno-fast-math as a command line option, which offers similar behaviour to MSVC’s /fp:strict option, although it is not sufficient on its own to get strict floating-point behaviour.

However, this relies on details of the implementation and is not future-proof. Additionally, optimisation flags can interfere with these floating-point arithmetic choices. This is a deeply unergonomic solution that is vulnerable to future decisions of the implementers.

Fences

Another approach might be to signal to the compiler that ordering and rounding properties be untouched in selected places by using a new kind of fence. This would be specified to prevent the compiler from reordering, combining or simplifying floating-point operations. It would ensure that the same number of rounding operations took place regardless of implementation, the same rounding modes were used regardless of implementation, and that the operations took place in the same order, regardless of implementation. Writing a larger floating-point expression as a series of single expression statements, interleaved with this new kind of fence, would detail to the compiler exactly what was wanted.

Again, this is a rather unergonomic solution. Additionally, given that the accuracy of operations of floating-point types are implementation-defined, this would require that latitude to be tightened and for floating-point behaviour to be fully specified. If that specification is to take place, then specifying reproducibility as defined in [559] would be the least astonishing thing to do.

Specification overhaul

The next approach to consider is following the example of prior art. Looking at COBOL and Fortran, their standards dictate that all floating-point arithmetic conforms to a revision of ISO/IEC 60559. Following this pattern, the C++ standard would need to be respecified to ensure that the object and value representations and accuracy of operations of the standard and extended optional floating-point types conform fully to a revision of ISO/IEC 60559.

While this solution offers much better ergonomics for the programmer, it would degrade performance of existing code, requiring programmers to interact more carefully with the command line options of their compilers.

#pragma options

[C++] does not require compilers to support any pragmas, but many implementations support the ISO C standard pragmas used to control floating-point behaviour:

#pragma STDC FENV_ACCESS
#pragma STDC FP_CONTRACT
#pragma STDC CX_LIMITED_RANGE

These could be introduced into the language, along with an additional pragma, for example:

While this seems like a simple solution, it would still require re-specification of existing types, which takes us very close to the specification overhaul described above. Also, it would require that any function called from within the block covered by the pragma would also be modified with such a pragma.

Function qualification

Just as class member functions can be qualified as const, a new function qualifier could be introduced which indicates that all floating-point operations within the scope of the function should be reproducible. This has the advantage of being a transitive property, like const, and users would be able to call other functions that are similarly qualified and retain reproducibility.

Unfortunately, we are in the same position with pragmas: this approach would still require re-specification of existing types.

A new core type

We have std::int_fast32_t, std::int_least32_t and the related types of differing width specifications. We also have the optional extended floating-point types std::float32_t and so on. We could take the same approach with floating-point types.

A new type could be introduced that supports all the specifications of Clause 11 of [559], particularly reproducible operations, format, status flags, and an implicit reproducible attribute, with a name like std::float_strict32_t. The specification of this type would prohibit reordering, fusion and variation in rounding. All rounding should be roundTiesToEven, as specified in Clause 4 of [559].

The requirement for reproducibility would mean that it would not be feasible to specify implementations of the contents of the <cmath> header. Users would be expected to provide all operations beyond addition, subtraction, multiplication and division.

This solution contains all the constraints of the problem to a single type, rather than requiring annotation distributed through code. However, introducing a new type to the core language is a considerable undertaking.

A new library type

There is a proposal in flight which partially solves this problem, [P2746R5] Deprecate and Replace Fenv Rounding Modes. From its abstract:

We argue that floating point rounding modes as specified by fesetround() are largely unusable, at least in C++. Furthermore, this implicit argument to floating point operations, which is almost never used, and largely opaque to the compiler, causes both language design and compilation problems. We thus make the somewhat drastic proposal to deprecate it, in spite of its long history, and to replace it with a much better-behaved facility.

In summary, that proposal proposes the elimination of functions fesetround and fegetround and replacing them with a struct which offers a very small set of the [559]-specified operations along with a rounding mode provided at construction. The value is not stored in the struct; rather, the member functions which are binary operations take two operands and apply the operation with the stored rounding mode. The struct effectively defines the rounding attribute as an object with a lifetime rather than as a #pragma qualifying a section of source code.

By additionally specifying that the operations cannot be fused or reordered, rounding mode would be respected, the correct number of roundings would be applied in the same order on all implementations, and all operations would be reproducible. This carries the penalty that a function call may be incurred for every operation, rather than some function calls being optimised away via FMA, for example. However, just as there is a trade-off between accuracy and speed, so is there a trade-off between reproducibility and speed. As described earlier, compiler flags specify optimisations that improve floating-point performance at the cost of accuracy. None of these optimisations would necessarily be available to this type.

Alternatively, a type can be proposed which is a companion type to the type proposed in [P2746R5]. At Six Impossible Things Before Breakfast, we have implemented several classes that we use to emulate reproducible floating-point operations. In the next section we shall look at the design of such a type.

6. API design

[559] conformance

From [559] Clause 11:

Reproducible floating-point numerical and status flag results are possible for reproducible operations, with reproducible attributes, operating on reproducible formats, defined for each language as follows:

• A reproducible operation is one of the operations described in Clause 5 or is a supported operation from 9.2, 9.3, 9.5, or 9.6.
• A reproducible attribute is an attribute that is required by a language standard in all implementations (see Clause 4).
• A reproducible status flag is one raised by the invalid operation, division by zero, or overflow exceptions (see 7.2, 7.3, and 7.4).
• A reproducible format is an arithmetic format that is also an interchange format (see Clause 3).

Programs that can be reliably translated into an explicit or implicit sequence of reproducible operations on reproducible formats produce reproducible results. That is, the same numerical and reproducible status flag results are produced.

Reproducible results require cooperation from language standards, language processors, and users. A language standard should support reproducible programming. Any conforming language standard supporting reproducible programming shall:

• Support the reproducible-results attribute.
• Support a reproducible format by providing all the reproducible operations for that format.
• Provide means to explicitly or implicitly specify any sequence of reproducible operations on reproducible formats supported by that language.

and shall explicitly define:

• Which language element corresponds to which supported reproducible format.
• How to specify in the language each reproducible operation on each supported reproducible format.
• One or more unambiguous expression evaluation rules that shall be available for user selection on all conforming implementations of that language standard, without deferring any aspect to implementations. If a language standard permits more than one interpretation of a sequence of operations from this standard it shall provide a means of specifying an unambiguous evaluation of that sequence (such as by prescriptive parentheses).
• A reproducible-results attribute, as described in 4.1, with values to indicate when reproducible results are required or reproducible results are not required. Language standards define the default value. When the user selects reproducible results required:
1. Execution behavior shall preserve the literal meaning (see 10.4) of the source code.
2. Conversions to and from external decimal character sequence shall not limit the maximum supported precision H (see 5.12.2).
3. Language processors shall indicate where reproducibility of operations that can affect the results of floating-point operations can not be guaranteed.
4. Only default exception handling (see Clause 7) shall be used.
If a language supports separately compiled routines (e.g., library routines for common functions) there must be some mechanism to ensure reproducible behavior.

To encapsulate enough of these requirements into a single type to achieve reproducibility, we need to decide which set of operations to support, which attributes to support and how, how to support the status flags and which interchange formats to support.

Operations

[559] Clause 5 lists a large number of operations over several sections. Many of these operations have analogues in [C++] within the <cmath> header, for example rounding and succession (5.3.1), some of the arithmetic operations (5.4.1), trigonometric functions (9.2) and so on.

[C++] does not specify the algorithms or rounding used in the calculation of these functions. It will not be feasible to expect reproducibility across implementations when the functions in <cmath> are invoked.

Therefore, we propose the set of operations to support should be limited to the fundamental operations of add, subtract, multiply and divide. All other functions should be implemented by the user using these operations.

We should also consider how to specify the addition of different reproducible types, since this may require conversion, which must be achieved reproducibly.

Attributes

The relevant attributes to support are rounding-direction and reproducibility. The rounding-direction attribute is discussed by [P2746R5]; it is bound up as a member of the type specified therein rather than explicitly specified for a block of code. Similarly, this proposal specifies the reproducibility attribute implicitly as part of this type.

However, the rounding attribute must be specified in some way, either dynamically as a data member of the type or statically as part of the type via a template parameter. There are four rounding attributes specified in [559] Clause 4 for a binary format implementation.

One attractive feature of the type would be for it to be a drop-in replacement for regular floating-point types. To do this, it would have to be the same size as those types. Specifying attributes dynamically would increase the size of the type beyond the size of the representation of the value being contained, so these should either be enumerated and provided as a template parameter, or, like the reproducibility parameter, be implicitly part of the type, thus requiring four separate types, one for each rounding mode.

If operations are requested on reproducible types with different rounding modes, which one should be applied? Clearly, the operand being modified should have its rounding mode respected.

Status flags

Status flags will be raised by the operations, so these need to be specified with them. At the moment, there is no specification in [C++] for responding to floating-point operations that behave unexpectedly; [basic.fundamental] paragraph 12 disposes of this by leaving everything implementation defined.

However, in <cfenv> there is provision for setting exception flags. The macros FE_DIVBYZERO, FE_INVALID and FE_OVERFLOW cover the required status flags. Although not specified by Clause 11 of [559], FE_INEXACT could also be specified for signalling.

Formats

Clause 3 of [559] defines four binary interchange formats binary16, binary32, binary64 and binary128. These are analogous to the optional extended floating-point types described in [C++] at [basic.extended.fp] named std::float16_t, std::float32_t, std::float64_t, and std::float128_t.

Class template or multiple similar classes?

With four formats and four rounding modes to choose from, we must decide if we want to specify sixteen similar types (one for each combination), or specify a class template.

This hints at a problem described earlier, optimisation. If the specification can forbid participation in optimisations which may reorder or fuse operations, or otherwise potentially change the result of an operation, then we have a solution. Otherwise, the implementation will need to be provided as a separate library out of reach of such optimisations, including link time optimisation.

7. Proposed wording

Feature test macro

8. Open questions

1. What other domains have uses for reproducibility?
2. [P2746R5] also specifies square root as an operation. Should this also be supported in this API? How about FMA? Any other operations?
3. Should the rounding attribute be specified statically or dynamically? If statically, should there be a separate type for each attribute, or should the attribute be a class template parameter?
4. Are all the rounding modes for binary arithmetic required? Are they used?
5. Can specification be provided which will prevent code-motion optimisations?
6. Should a single class template, over format and rounding, be provided, or sixteen individual classes?
7. Can the type proposed in [P2746R5] fully meet all the requirements of reproducibility, or can the proposal be easily modified so to do?
8. Clause 11 of [559] does NOT specify the inexact status flag as one of the reproducible flags. Is the overspecification in the wording problematic?
9. Can core language types be specified “as-if” they were a full specialisation of std::strict_float? E.g. std::float32_strict_rte_t is specified “as-if” it were std::strict_float<std::float32_t, std::round_ties_to_even>?
10. Should the non-member binary operators be hidden friends?
11. Should there be member functions which perform the operations and offer the user the opportunity to override the rounding mode?
12. Should <cmath> be aware of this type? How might those functions be consistently specified?

9. References

[559] ISO/IEC JTC 1/SC 25 (May 2020). ISO/IEC 60559:2020 — Information technology — Microprocessor Systems — Floating-Point arithmetic. ISO. pp. 1–74.

[C++] C++23 Working Draft https://www.open-std.org/jtc1/sc22/wg21/docs/papers/2023/n4950.pdf

[CEX] Compiler Explorer https://godbolt.org/z/W8shs6dao

[GOL] What Every Computer Scientist Should Know About Floating-Point Arithmetic, David Goldberg, https://dl.acm.org/doi/10.1145/103162.103163

[P2746R5] Deprecate and Replace Fenv Rounding Modes, Hans Boehm, https://www.open-std.org/jtc1/sc22/wg21/docs/papers/2024/p2746r5.pdf

[VGI] Video game industry - Statistics & Facts https://www.statista.com/topics/868/video-games/