From: John Reid Subject: Procedures technical report Date: 1 November 1995 0. NOTES This is a draft Technical Report using procedures in a module for exception handling. It is based on Keith Bierman's message of 7 September. For IEEE hardware, it supports exceptions and other aspects of IEEE arithmetic. Non-IEEE processors are required to support only OVERFLOW and DIVIDE_BY_ZERO, and may support more. We could reduce this proposal to one that supports exceptions only. This would involve removal of the data type IEEE_RND_VALUE and the reduction of the set of procedures to IEEE_GET_FLAG(FLAG), IEEE_CLEAR_FLAGS, IEEE_SET_FLAGS, IEEE_GET_FSR, IEEE_HALT, and IEEE_SET_FSR. More work is needed yet on the interaction of these features with PURE and FORALL. 1. RATIONALE Exception handling is required for the development of robust and efficient numerical software. In particular, it is necessary in order to be able to write portable scientific libraries. In numerical Fortran programming, current practice is to employ whatever exception handling mechanisms are provided by the system/vendor. This clearly inhibits the production of fully portable numerical libraries and programs. It is particularly frustrating now that IEEE arithmetic is so widely used, since built into it are the five conditions: overflow, invalid, divide_by_zero, underflow, and inexact. Our aim to to provide support for these conditions. This proposal involves a standard module, IEEE_ARITHMETIC, containing three derived types, some named constants of these types, and a collection of simple procedures. They allow the flags to be tested, cleared, set, saved, or restored. To facilitate maximum performance, each of the proposed functions does very little processing of arguments. In many cases, a processor may generate only a few inline machine code instructions rather than library calls. In order to allow for the maximum number of processors to provide the maximum value to users, we do NOT require complete IEEE conformance. What a vendor must do is to provide the module, and support the overflow and divide-by-zero flags and the basic inquiry functions. A user must utilize the inquiry functions to determine if they can count on specific features of the IEEE standard, before using other inquiry or value setting procedures. Note that a processor ought not implement these as "macros", as IEEE conformance is often controlled by compiler switches. A processor which offers a switch to turn off a facility should adjust the values returned for these inquiries. For example, a processor which allows gradual underflow to be turned off (replaced with flush to zero) should return .FALSE. for IEEE_SUPPORT_DENORMAL(X) when a source file is processed with that option on. Naturally it should return .TRUE. when that option is not in effect. The most important use of a floating-point exception handling facility is to make possible the development of much more efficient software than is otherwise possible. The following "hypotenuse" function, sqrt(x**2 + y**2), illustrates the use of the facility in developing efficient software. REAL FUNCTION HYPOT(X, Y) ! In rare circumstances this may lead to the setting of the OVERFLOW flag USE IEEE_ARITHMETIC REAL X, Y REAL SCALED_X, SCALED_Y, SCALED_RESULT LOGICAL OLD_OVERFLOW, OLD_UNDERFLOW INTRINSIC SQRT, ABS, EXPONENT, MAX, DIGITS, SCALE ! Store the old flags and clear them OLD_OVERFLOW = IEEE_GET_FLAG(OVERFLOW) OLD_UNDERFLOW = IEEE_GET_FLAG(UNDERFLOW) CALL IEEE_CLEAR_FLAGS(UNDERFLOW, OVERFLOW) ! Try a fast algorithm first HYPOT = SQRT( X**2 + Y**2 ) IF (IEEE_GET_FLAG(UNDERFLOW) .OR. IEEE_GET_FLAG(OVERFLOW) ) THEN CALL IEEE_CLEAR_FLAGS(UNDERFLOW, OVERFLOW) IF ( X==0.0 .OR. Y==0.0 ) THEN HYPOT = ABS(X) + ABS(Y) ELSE IF ( 2*ABS(EXPONENT(X)-EXPONENT(Y)) > DIGITS(X)+1 ) THEN HYPOT = MAX( ABS(X), ABS(Y) )! one of X and Y can be ignored ELSE ! scale so that ABS(X) is near 1 SCALED_X = SCALE( X, -EXPONENT(X) ) SCALED_Y = SCALE( Y, -EXPONENT(X) ) SCALED_RESULT = SQRT( SCALED_X**2 + SCALED_Y**2 ) HYPOT = SCALE( SCALED_RESULT, EXPONENT(X) ) ! may cause overflow END IF END IF IF(OLD_OVERFLOW) CALL IEEE_SET_FLAGS(OVERFLOW) IF(OLD_UNDERFLOW) CALL IEEE_SET_FLAGS(UNDERFLOW) END FUNCTION HYPOT An attempt is made to evaluate this function directly in the fastest possible way. (Note that with hardware support, exception checking is very efficient; without language facilities, reliable code would require programming checks that slow the computation significantly.) The fast algorithm will work almost every time, but if an exception occurs during this fast computation, a safe but slower way evaluates the function. This slower evaluation may involve scaling and unscaling, and in (very rare) extreme cases this unscaling can cause overflow (after all, the true result might overflow if X and Y are both near the overflow limit). If the overflow or underflow flag is set on entry, it is reset on return, so that earlier exceptions are not lost. Can all this be accomplished without the help of an exception handling facility? Yes, it can - in fact, the alternative code can do the job, but of course it is much less efficient. That's the point. The HYPOT function is special, in this respect, in that the normal and alternative codes try to accomplish the same task. This is not always the case. In fact, it very often happens that the alternative code concentrates on handling the exceptional cases and is not able to handle all of the non-exceptional cases. When this happens, a program which cannot take advantage of hardware flags could have a structure like the following: if ( in the first exceptional region ) then handle this case else if ( in the second exceptional region ) then handle this case : else execute the normal code end But this is not only inefficient, it also "inverts" the logic of the computation. For other examples, see Hull, Fairgrieve and Tang (1994) and Demmel and Li (1994). The HYPOT algorithm is used in example 4 of section 15.10. The code for the HYPOT function can be generalized in an obvious way to compute the Euclidean Norm, sqrt( x(1)**2 + x(2)**2 + ... + x(n)**2 ), of an n-vector; the generalization of the alternative code is not so obvious (though straightforward) and will be much slower relative to the normal code than is the case with the HYPOT function. Of the five IEEE conditions, underflow and inexact are obviously milder than the other three (inexact is particularly mild and signals almost all the time in typical codes). We therefore group the other three together as USUAL. There is a need for further intrinsic conditions in connection with reliable computation. Examples are a. INSUFFICIENT_STORAGE for when the processor is unable to find sufficient storage to continue execution. b. INTEGER_OVERFLOW and INTEGER_DIVIDE_BY_ZERO for when an intrinsic integer operation has a very large result or has a zero denominator. c. INTRINSIC for when an intrinsic procedure has been unsuccessful. d. SYSTEM_ERROR for when a system error occurs. This proposal has been designed to allow such enhancements in the future. References Demmel, J.W. and Li, X. (1994). Faster Numerical Algorithms via Exception Handling. IEEE Transactions on Computers, 43, no. 8, 983-992. Hull, T.E., Fairgrieve, T.F., and Tang, T.P.T. (1994). Implementing complex elementary functions using exception handling. ACM Trans. Math. Software 20, 215-244. 2. TECHNICAL SPECIFICATION This proposal involves a standard module, IEEE_ARITHMETIC, containing three derived types: 1. IEEE_FLAG, with the constants INVALID, OVERFLOW, DIVIDE_BY_ZERO, USUAL, UNDERFLOW, INEXACT, and ALL. 2. IEEE_RND_VALUE, with the constants NEAREST, TO_ZERO, UP, DOWN, and OTHER 3. IEEE_FSR_VALUE, for saving the current floating point status. and a collection of simple procedures. There must be flags that are initially clear and that are set when an exception occurs. The value of a flag is determined by the function IEEE_GET_FLAG(flag) Flags may be cleared by the subroutine IEEE_CLEAR_FLAGS(flag-list) and may be set by the subroutine IEEE_SET_FLAGS(flag-list) The minimum requirement is for the support of OVERFLOW and DIVIDE_BY_ZERO. Whether the other exceptions are supported may be determined by the inquiry function IEEE_SUPPORT_FLAG(FLAG) The extent of further support for the IEEE standard may be determined by the inquiry functions: IEEE_DATATYPE (X) Inquire if the processor supports IEEE arithmetic for reals of the same kind type parameter as the argument. IEEE_SUPPORT_ALL() Inquire if processor supports all the IEEE facilities defined in this standard for all kinds of reals. IEEE_SUPPORT_DENORMAL(X) Inquire if the processor supports IEEE gradual underflow for reals of the same kind type parameter as the argument. IEEE_SUPPORT_FLAG(FLAG) Inquire if the processor supports an exception or set of exceptions. IEEE_SUPPORT_HALTING() Inquire if the processor supports the ability to control during program execution whether to abort or continue execution after an exception. IEEE_SUPPORT_INF(X) Inquire if processor supports the IEEE infinity facility for reals of the same kind type parameter as the argument. IEEE_SUPPORT_NAN(X) Inquire if processor supports the IEEE Not-A-Number facility for reals of the same kind type parameter as the argument. IEEE_SUPPORT_ROUNDING (RND_VALUE) Inquire if processor supports changing to a particular rounding mode for IEEE kinds of reals. IEEE_SUPPORT_SQRT(X) Inquire if the processor supports IEEE square root for reals of the same kind type parameter as the argument. Implementations with the appropriate extent of IEEE support provide: A. Elemental functions: IEEE_COPYSIGN(X,Y) IEEE copysign function. IEEE_INFINITY(X) Generate IEEE infinity with the same sign as X. IEEE_IS_DENORMAL(X) Determine if value is IEEE denormalized. IEEE_IS_INF(X) Determine if value is IEEE Infinity. IEEE_IS_NAN(X) Determine if value is IEEE Not-a-Number. IEEE_IS_NORMAL(X) Whether a value is normal, that is, neither an Infinity, a NaN, nor denormalized. IEEE_LOGB(X) Unbiased exponent in the IEEE floating point format. IEEE_NAN(X) Generate IEEE Not-a-Number. IEEE_NEXTAFTER(X,Y) Returns the next representable neighbor of X in the direction toward Y. IEEE_SCALB (X,I) Returns X * 2**I. IEEE_SQRT(X) Square root as defined in the IEEE standard. IEEE_UNORDERED(X,Y) IEEE unordered function. True if X or Y is a NaN and false otherwise. B. Scalar functions: IEEE_GET_FSR() Save the current values of the set of flags that define the current state of the floating point environment (usually the floating point status register). The result is of type IEEE_FSR_VALUE. IEEE_GET_ROUNDING_MODE() Save the current IEEE rounding mode. The result is of type IEEE_RND_VALUE. C. Subroutines: IEEE_CLEAR_FLAGS (FLAG1 [, FLAG2, ...]) Clear one or more exception flags. IEEE_HALT(FLAG,HALTING) Controls continuation or halting on exceptions. FLAG is of type IEEE_FLAG and value INVALID, OVERFLOW, DIVIDE_BY_ZERO, USUAL, UNDERFLOW, INEXACT, or ALL. IEEE_SET_FLAGS(FLAG1 [, FLAG2, ...]) Set one or more exception flags. IEEE_SET_FSR(FSR_VALUE) Restore the values of the set of flags that define the current state of the floating point environment (usually the floating point status register). FSR_VALUE is of type IEEE_FSR_VALUE and has been set by a call of IEEE_GET_FSR. IEEE_SET_ROUNDING_MODE(RND_VALUE) Set the current IEEE rounding mode. RND_VALUE is of type IEEE_RND_VALUE, with value NEAREST, TO_ZERO, UP, or DOWN. Note: An IEEE processor is required to support single precision, and an extended precision. This is usually REAL and DOUBLE PRECISION but it might not always be. In addition, some vendors support a "REAL*16" which does not follow the "logical extension" to the IEEE 754 standard (e.g. IBM's RS6K "doubled double" implementation of REAL*16). Users should be able to determine if a particular datatype is going to behave as "one would expect". But implementors should not be constrained to supply only datatypes which conform to IEEE 754 (for performance, or compatibility with older processors or eventual improvements in state of the art in floating point implementation). If a processor supports a particular IEEE exception, the associated flag is always clear. 3. EDITS TO THE DRAFT STANDARD (N1122) 130/11. Add: 'If any exception (15) is set, the processor shall issue a warning on the unit identified by * in a WRITE statement, indicating which exceptions are set.'. 288+. Add <<15. Intrinsic module for support of exceptions and IEEE arithmetic>> The intrinsic module IEEE_ARITHMETIC provides support for the floating point exceptions associated with overflow and divide by zero for real or complex data with any kind parameter. If the type and kind parameter implements IEEE arithmetic (IEEE 754-1985), the module provides more extensive support, including the exceptions invalid, underflow, and inexact. There are inquiry functions that permit the extent of support of IEEE 754-1985 to be determined. <<15.1 Derived data types defined in the module>> The module contains three data types, whose components are private. No operation is defined for any of them and only intrinsic assignment is available for them. They are: (i) TYPE (IEEE_FLAG) is used to identify exception flags or sets of flags. The constants INVALID, OVERFLOW, DIVIDE_BY_ZERO, USUAL, UNDERFLOW, INEXACT, and ALL are defined for it. (ii) TYPE (IEEE_RND_VALUE) is used to identify rounding modes. The constants NEAREST, TO_ZERO, UP, DOWN, and OTHER are defined for it. (iii) TYPE (IEEE_FSR_VALUE) is used for saving the current floating point status. No constants are defined for it. <<15.2 The exceptions>> The exceptions are: OVERFLOW This exception occurs when the result for an intrinsic real operation has a very large processor-dependent absolute value, or the real or imaginary part of the result for an intrinsic complex operation has a very large processor-dependent absolute value. DIVIDE_BY_ZERO This exception occurs when a real or complex division has a nonzero numerator and a zero denominator. INVALID This exception occurs when a real or complex operation is invalid. UNDERFLOW This exception occurs when the result for an intrinsic real operation has a very small processor-dependent absolute value, or the real or imaginary part of the result for an intrinsic complex operation has a very small processor-dependent absolute value. INEXACT This exception occurs when the result of a real or complex operation is not exact. Each exception has a flag whose value is either clear or set. The value may be determined by the function IEEE_GET_FLAG. Its initial value is clear and it becomes set when the associated exception occurs. It may also be set by the subroutine IEEE_SET_FLAGS or the subroutine IEEE_SET_FSR. Once set, it remains set unless cleared by an invocation of the subroutine IEEE_CLEAR_FLAGS or the subroutine IEEE_SET_FSR. If any exception is set when the program terminates, the processor shall issue a warning on the unit identified by * in a WRITE statement, indicating which conditions are set. If a flag is set on entry to a pure procedure, the flag shall not be cleared unless it is reset before return. In a sequence of statements that contains no exception handling statements, if the execution of a process would cause an exception to be set but after execution of the sequence no value of a variable depends on the process, whether the exception is set is processor dependent. For example, when Y has the value zero, whether the code X = 1.0/Y X = 3.0 sets DIVIDE_BY_ZERO is processor dependent. An exception must not set if this could arise only during execution of a process further to those required or permitted by the standard. For example, the statement IF (F(X)>0.) Y = 1.0/Z must not set DIVIDE_BY_ZERO when both F(X) and Z are zero and the statement WHERE(A>0.) A = 1.0/A must not set DIVIDE_BY_ZERO. On the other hand, when X has the value 1.0 and Y has the value 0.0, the expression X>0.00001 .OR. X/Y>0.00001 is permitted to cause the setting of DIVIDE_BY_ZERO. Note: It is intended that implementations be free to use code motion techniques. The processor need not support INVALID, UNDERFLOW, and INEXACT. If an exception is not supported, its flag is always clear. The function IEEE_SUPPORT_FLAG may be used to inquire whether a particular flag is supported. The constant USUAL is associated with the list OVERFLOW, DIVIDE_BY_ZERO, INVALID. The constant ALL is associated with the list of all the flags. <<15.3 The rounding modes>> IEEE 754-1985 specifies four rounding modes: NEAREST rounds the exact result to the nearest representable value. TO_ZERO rounds the exact result towards zero to the next representable value. UP rounds the exact result towards +infinity to the next representable value. DOWN rounds the exact result towards -infinity to the next representable value. The function IEEE_GET_ROUNDING_MODE may be used to inquire which rounding mode is in operation. Its value is one of the above four or OTHER if the rounding mode does not conform to IEEE 754-1985. If the processor supports the alteration of the rounding mode during execution, the subroutine IEEE_SET_ROUNDING_MODE may be used to alter it. The function IEEE_SUPPORT_ROUNDING may be used to inquire whether this facility is available for a particular mode. <<15.4 Halting>> If the processor supports the ability to control during program execution whether to abort or continue execution after an exception, such control may be exercised by invocation of the subroutine IEEE_HALT. Halting is not precise and may occur any time after the exception has occurred. The function IEEE_SUPPORT_HALTING may be used to inquire whether this facility is available. <<15.5 The floating point status>> The values of all the supported flags for exceptions, rounding mode, and halting may be saved in a scalar variable of type TYPE(IEEE_FSR_VALUE) with the function IEEE_GET_FSR and restored with the subroutine IEEE_SET_FSR. There are no facilities for finding the values of particular flags held within such a variable. Note. Some processors hold all these flags in a floating point status register that can be saved and restored as a whole much faster than all individual flags can be saved and restored. These procedures are provided to exploit this feature. <<15.6 Exceptional values>> IEEE 754-1985 specifies the following exceptional floating point values: Denormalized values have very small absolute values and lowered precision. Infinite values are created by overflow or division by zero. Not-a-Number (NaN) values are undefined values or values created by an invalid operation. The functions IEEE_IS_DENORMAL, IEEE_IS_INF, IEEE_IS_NAN are provided to test whether a value is denormalized, infinite, or Nan. The functions IEEE_INFINITY and IEEE_NAN are provided to generate an infinity or a Nan. The function IEEE_IS_NORMAL is provided to test whether a value is none of these. The functions IEEE_SUPPORT_DENORMAL, IEEE_SUPPORT_INF, and IEEE_SUPPORT_NAN may be used to inquire whether this facility is available for a particular kind of real. <<15.7 IEEE arithmetic>> The function IEEE_DATATYPE may be used to inquire whether IEEE arithmetic is available for a particular kind of real. Complete conformance with the IEEE standard is not required, but the normalized numbers must be exactly those of IEEE single, IEEE double, or the corresponding quad with an exponent width of 15 bits and a mantissa width of 112 bits; the arithmetic operators must be implemented with at least one of the IEEE rounding modes; and the functions copysign, scalb, logb, nextafter, and unordered must be provided by the functions IEEE_COPYSIGN, IEEE_SCALB, IEEE_LOGB, IEEE_NEXTAFTER, and IEEE_UNORDERED. IEEE 754-1985 specifies a square root function that returns -0 for the square root of -0. The function IEEE_SQRT is provided for this and the function IEEE_SUPPORT_SQRT may be used to inquire whether this facility is available for a particular kind of real. The inquiry function SUPPORT_ALL is provided to inquire whether the processor supports all the IEEE facilities defined in this standard for all kinds of reals. <<15.8 Tables of the procedures>> In this section, the procedures are tabulated with the names of their arguments and a short description. <<15.8.1 Inquiry functions>> IEEE_DATATYPE (X) Inquire if the processor supports IEEE arithmetic for reals of the same kind type parameter as the argument. IEEE_SUPPORT_ALL() Inquire if processor supports all the IEEE facilities defined in this standard for all kinds of reals. IEEE_SUPPORT_DENORMAL(X) Inquire if the processor supports IEEE gradual underflow for reals of the same kind type parameter as the argument. IEEE_SUPPORT_FLAG(FLAG) Inquire if the processor supports an exception or set of exceptions. IEEE_SUPPORT_HALTING() Inquire if the processor supports the ability to control during program execution whether to abort or continue execution after an exception. IEEE_SUPPORT_INF(X) Inquire if processor supports the IEEE infinity facility for reals of the same kind type parameter as the argument. IEEE_SUPPORT_NAN(X) Inquire if processor supports the IEEE Not-A-Number facility for reals of the same kind type parameter as the argument. IEEE_SUPPORT_ROUNDING (RND_VALUE) Inquire if processor supports changing to a particular rounding mode for IEEE kinds of reals. IEEE_SUPPORT_SQRT(X) Inquire if the processor supports IEEE square root for reals of the same kind type parameter as the argument. <<15.8.2 Scalar functions>> IEEE_GET_FLAG(FLAG) Get an exception flag. IEEE_GET_FSR() Save the current values of the set of flags that define the current floating point status, including all the exception flags. IEEE_GET_ROUNDING_MODE() Save the current IEEE rounding mode. <<15.8.3 Elemental functions>> IEEE_COPYSIGN(X,Y) IEEE copysign function. IEEE_INFINITY(X) Generate IEEE infinity with the same sign as X. IEEE_IS_DENORMAL(X) Determine if value is IEEE denormalized. IEEE_IS_INF(X) Determine if value is IEEE Infinity. IEEE_IS_NAN(X) Determine if value is IEEE Not-a-Number. IEEE_IS_NORMAL(X) Whether a value is normal, that is, neither an Infinity, a NaN, nor denormalized. IEEE_LOGB(X) Unbiased exponent in the IEEE floating point format. IEEE_NAN(X) Generate IEEE Not-a-Number. IEEE_NEXTAFTER(X,Y) Returns the next representable neighbor of X in the direction toward Y. IEEE_SCALB (X,I) Returns X * 2**I. IEEE_SQRT(X) Square root as defined in the IEEE standard. IEEE_UNORDERED(X,Y) IEEE unordered function. True if X or Y is a NaN and false otherwise. <<15.8.4 Subroutines>> IEEE_CLEAR_FLAGS (FLAG1 [, FLAG2, ...]) Clear one or more exception flags. IEEE_HALT(FLAG,HALTING) Controls continuation or halting on exceptions. IEEE_SET_FLAGS(FLAG1 [, FLAG2, ...]) Set one or more exception flags. IEEE_SET_FSR(FSR_VALUE) Restore the values of the set of flags that define the the floating point status. IEEE_SET_ROUNDING_MODE(RND_VALUE) Set the current IEEE rounding mode. <<15.9 Specifications of the procedures>> IEEE_CLEAR_FLAGS (FLAG1 [, FLAG2, ...]) Description. Clear one or more exception flags. Class. Subroutine. Arguments. Each argument shall be of type TYPE(IEEE_FLAG) and with value OVERFLOW, DIVIDE_BY_ZERO, INVALID, USUAL, UNDERFLOW, INEXACT, or ALL, which are named PARAMETERS defined in the module. They are INTENT(IN) arguments and specify the IEEE flags to be cleared. Specifying USUAL is equivalent to specifying OVERFLOW, DIVIDE_BY_ZERO, INVALID. Specifying ALL is equivalent to specifying all the flags. Example. CALL IEEE_CLEAR_FLAGS(OVERFLOW,UNDERFLOW) clears the overflow and underflow flags. IEEE_COPYSIGN(X,Y) Description. IEEE copysign function. Class. Elemental function. Arguments. The arguments shall be of type real and such that both IEEE_DATATYPE(X) and IEEE_DATATYPE(Y) have the value true. Result Characteristics. Same as X. Result Value. The result has the value of X with the sign of Y. This is true even for IEEE special values, such as Nan and Inf (on processors supporting such values). Examples. The value of COPYSIGN(X,1.0) is ABS(X) even when X is NaN. The value of COPYSIGN(-X,1.0) is X copied with its sign reversed, not 0-x; the distinction is germane when X is +0, -0, or NaN. IEEE_DATATYPE(X) Description. Inquire if the processor supports IEEE arithmetic for reals of the same kind type parameter as the argument. Class. Inquiry function. Argument. X shall be scalar and of type real. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports IEEE arithmetic for real and complex variables of the same kind type parameter as X; otherwise, it has the value false. If the result is true, the normalized numbers shall be exactly those of IEEE single, IEEE double, or the corresponding quad with an exponent width of 15 bits and a mantissa width of 112 bits; the arithmetic operators shall be implemented with at least one of the IEEE rounding modes; and the functions IEEE_IS_NORMAL, IEEE_COPYSIGN, IEEE_SCALB, IEEE_LOGB, IEEE_NEXTAFTER, and IEEE_UNORDERED shall be provided. Example. IEEE_DATATYPE(1.0) has the value .TRUE. if default reals are implemented as in the IEEE standard except that underflowed values flush to 0.0 instead of being denormal. IEEE_GET_FLAG(FLAG) Description. Get an exception flag. Class. Transformational function. Argument. FLAG shall be scalar and of type TYPE(IEEE_FLAG). Its value shall be OVERFLOW, DIVIDE_BY_ZERO, INVALID, USUAL, UNDERFLOW, INEXACT, or ALL, which are named PARAMETERS defined in the module. It specifies the IEEE flag to be inspected. Result Characteristics. Default logical scalar. Result Value. Case (i): If the value of FLAG is not USUAL or ALL, the result value is true if the exception flag is set and is false otherwise. Case (ii): If the value of FLAG is USUAL, the result value is true if INVALID, OVERFLOW, or DIVIDE_BY_ZERO is set and is false otherwise. Case (iii): If the value of FLAG is ALL, the result value is true if any flag is set and is false otherwise. Example. The value of IEEE_GET_FLAG(OVERFLOW) is true if the overflow flag is set and is false if it is clear. IEEE_GET_FSR() Description. Save the current values of the set of flags that define the current floating point status, including all the exception flags. Class. Transformational function. Arguments. None. Result Characteristics. Scalar of type TYPE(IEEE_FSR_VALUE). Result Value. The result value is that of the floating point status. Example. To store all the exceptions flags, do a calculation involving exception handling, and restore them later: USE IEEE_ARITHMETIC TYPE(IEEE_FSR_VALUE) FSR_VALUE : FSR_VALUE = IEEE_GET_FSR() ! Store the flags CALL IEEE_CLEAR_FLAGS (ALL) : ! calculation involving exception handling CALL IEEE_SET_FSR(FSR_VALUE) ! Restore the flags Note. The result can be used only in an IEEE_SET_FSR invocation. IEEE_GET_ROUNDING_MODE() Description. Save the current IEEE rounding mode. Class. Transformational function. Arguments. None. Result Characteristics. Scalar of type TYPE(IEEE_RND_VALUE). Result Value. The result value is that of the floating point rounding mode, and may be NEAREST, TO_ZERO, UP, DOWN, or OTHER. Example. To store the rounding mode, do a calculation with round to nearest, and restore the rounding mode later: USE IEEE_ARITHMETIC TYPE(IEEE_RND_VALUE) RND_VALUE : RND_VALUE = IEEE_GET_ROUNDING_MODE() ! Store the rounding mode CALL IEEE_SET_ROUNDING_MODE (NEAREST) : ! calculation with round to nearest CALL IEEE_SET_ROUNDING_MODE(RND_VALUE) ! Restore the rounding mode Note. The result can legally be used only in an IEEE_SET_ROUNDING_MODE invocation. IEEE_HALT(FLAG,HALTING) Description. Controls continuation or halting on exceptions. Class. Subroutine. IEEE_SUPPORT_HALTING() shall be true. Arguments. FLAG shall be scalar and of type TYPE(IEEE_FLAGS). It is of INTENT(IN) and shall have one of the values: OVERFLOW, DIVIDE_BY_ZERO, INVALID, USUAL, UNDERFLOW, INEXACT, or ALL, which are named PARAMETERS defined in the module. USUAL specifies the exceptions INEXACT, OVERFLOW, and DIVIDE_BY_ZERO. ALL specifies all exceptions. HALTING shall be scalar and of type logical. It is of INTENT(IN). If the value is true, the exception or exceptions specified by FLAG will cause halting. Otherwise, execution will continue after these exceptions. Example. CALL IEEE_HALT(DIVIDE_BY_ZERO,.TRUE.) causes halting after a divide_by_zero exception. Note: Halting is not precise and may occur some time after the exception has occurred. IEEE_INFINITY(X) Description. Generate IEEE infinity. Class. Elemental function. Argument. X shall be of type real and such that IEEE_SUPPORT_INF(X) has the value true. Result Characteristics. Same as X. Result Value. The result value is the IEEE infinity with the sign of X. Example. IEEE_INFINITY(-1.0) has the value -Infinity. IEEE_IS_DENORMAL(X) Description. Whether a value is IEEE denormalized. Class. Elemental function. Argument. X shall be of type real and such that IEEE_SUPPORT_DENORMAL(X) has the value true. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the value of X is IEEE denormalized; otherwise, it has the value false. Example. IEEE_IS_DENORMAL(TINY(X)/2) has the value true. IEEE_IS_INF(X) Description. Whether a value is IEEE Infinity. Class. Elemental function. Argument. X shall be of type real and such that IEEE_SUPPORT_INF(X) has the value true. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the value of X is an IEEE Infinity (either positive or negative); otherwise, it has the value false. Example. IEEE_IS_INF(IEEE_INFINITY(-1.0)) has the value true. Note: the sign of X can be determined directly via comparison (e.g. (X > 0) is true for +Inf) or indirectly via SIGN or COPYSIGN. IEEE_IS_NAN(X) Description. Whether a value is IEEE Not-a-Number. Class. Elemental function. Argument. X shall be of type real and such that IEEE_SUPPORT_NAN(X) has the value true. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the value of X is an IEEE NaN; otherwise, it has the value false. Example. IEEE_IS_NAN(IEEE_SQRT(-1.0)) has the value true. IEEE_IS_NORMAL(X) Description. Whether a value is normal, that is, neither an Infinity, a NaN, nor denormalized. Class. Elemental function. Argument. X shall be of type real and such that IEEE_DATATYPE(X) has the value true. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the value of X is neither an IEEE Infinity (either positive or negative), an IEEE NaN, nor denormalized; otherwise, it has the value false. Example. IEEE_IS_NORMAL(IEEE_SQRT(-1.0)) has the value false. Note: While one could code this via a combination of the primitives, it is often significantly more efficient to just test it once. It also makes for more readable code. IEEE_LOGB(X) Description. Unbiased exponent in the IEEE floating point format. Class. Elemental function. Argument. X shall be of type real and such that IEEE_DATATYPE (X) has the value true. Result Characteristics. Default integer. Result Value. Case (i): If the value of X is neither zero, infinity, nor NaN, the result has the value of the unbiased exponent of X. Case(ii): If X==0, the result is -HUGE(X). Case(iii): If X has an infinite value, the result is +Inf. Case(iv): If X has a NaN value, the result is X. Note: Cases (iii) and (iv) cannot occur on processors that lack Inf and NaN values. Example. IEEE_LOGB(-1.1) has the value 0. Comment to Keith: I have made the result integer. IEEE_NAN(X) Description. Generate IEEE Not-a-Number. Class. Elemental function. Argument. X shall be of type real and such that IEEE_SUPPORT_NAN(X) has the value true. Result Characteristics. Same as X. Result Value. The result value is an IEEE signaling NaN. Example. IEEE_NAN(X) is an IEEE signaling NaN. Note: IEEE defines an entire family of values as NaN. A processor may choose to allow users to specify the bit pattern. If so, the value of X is appropriately transformed to be a signaling NaN with all the low order bits taken from the value of X. Comment: one can argue that IEEE_NAN and IEEE_INF are not needed, because users can code expressions to generate them. However, experience shows that optimizers often make it difficult to do this easily and portably. It is better to have these functions in the Standard to promote portability, without impacting optimization. IEEE_NEXTAFTER(X,Y) Description. Returns the next representable neighbor of X in the direction toward Y. Class. Elemental function. Arguments. The arguments shall be of type real and such that IEEE_DATATYPE(X) and IEEE_DATATYPE(Y) have the value true. Result Characteristics. Same as X. Result Value. Case (i): If X == Y, the result is X without any exception ever being signaled. Case (ii): If X /= Y, the result has the value of the next representable neighbor of X in the direction of Y. If either X or Y is a NaN, the result is one or the other of the input NaNs. Overflow is signaled when X is finite but NEXTAFTER(X,Y) is infinite; underflow is signaled when NEXTAFTER(X,Y) is denormalized; in both cases, INEXACT is signaled. Example. The value of IEEE_NEXTAFTER(1.0,2.0) is 1.0+EPSILON(X). IEEE_SCALB(X,I) Description. Returns X * 2**I. Class. Elemental function. Arguments. X shall be of type real and such that IEEE_DATATYPE(X) has the value true. I shall be of type integer. Result Characteristics. Same as X. Result Value. Case (i): If X * 2**I is within range, the result has this value. Case (ii): If X * 2**I is too large and IEEE infinities are supported, the result value is infinity with the sign of X. Case (iii): If X * 2**I is too large and IEEE infinities are not supported, the result value is SIGN(HUGE,X). Example. The value of IEEE_SCALB(1.0,2) is 4.0. IEEE_SET_FLAGS(FLAG1 [, FLAG2, ...]) Description. Set one or more supported exception flags. Class. Subroutine. Arguments. The arguments shall be of type TYPE(IEEE_FLAG) and with value OVERFLOW, DIVIDE_BY_ZERO, INVALID, USUAL, UNDERFLOW, INEXACT, or ALL, which are named constants defined in the module. They are INTENT(IN) arguments and specify the exception flags to be set. USUAL is equivalent to the list OVERFLOW, DIVIDE_BY_ZERO, INVALID. ALL is equivalent to the list of all flags. If a flag is not supported on the processor, it remains clear. Example. CALL IEEE_SET_FLAGS(OVERFLOW) sets the overflow flag. IEEE_SET_FSR(FSR_VALUE) Description. Restore the values of the set of flags that define the the floating point status. Class. Subroutine. Argument. FSR_VALUE shall be scalar and of type TYPE(IEEE_FSR_VALUE). Its value shall have been set in a previous invocation of IEEE_GET_FSR. Example. To store all the exceptions flags, do a calculation involving exception handling, and restore them later: USE IEEE_ARITHMETIC TYPE(IEEE_FSR_VALUE) FSR_VALUE : FSR_VALUE = IEEE_GET_FSR() ! Store the flags CALL IEEE_CLEAR_FLAGS (ALL) : ! calculation involving exception handling CALL IEEE_SET_FSR(FSR_VALUE) ! Restore the flags Note: getting and setting may be expensive operations. It is the programmer's responsibility to do it when necessary to assure correct results. IEEE_SET_ROUNDING_MODE (RND_VALUE) Description. Set the current IEEE rounding mode. Class. Subroutine. IEEE_SUPPORT_ROUNDING(RND_VALUE) shall be true. Argument. RND_VALUE shall be scalar and of type TYPE(IEEE_RND_VALUE). Its value shall be one of NEAREST, TO_ZERO, UP, and DOWN, which are named constants of the module. Example. To store the rounding mode, do a calculation with round to nearest, and restore the rounding mode later: USE IEEE_ARITHMETIC TYPE(IEEE_RND_VALUE) RND_VALUE : RND_VALUE = IEEE_GET_ROUNDING_MODE() ! Store the rounding mode CALL IEEE_SET_ROUNDING_MODE (NEAREST) : ! calculation with round to nearest CALL IEEE_SET_ROUNDING_MODE(RND_VALUE) ! Restore the rounding mode IEEE_SQRT(X) Description. Square root as defined in the IEEE standard. Class. Elemental function. Argument. X shall be of type real and such that IEEE_SUPPORT_SQRT(X) has the value true. Result characteristics. Same as X. Result value. The result has the value of the square-root of X as defined in the IEEE Standard. Example. IEEE_SQRT(-0.) has the value -0. Note: On most IEEE processors, SQRT is probably the same as IEEE_SQRT. However, processors without adequate hardware support may choose to return different values in hard cases, in which case a careful programmer may wish to use IEEE_SQRT for safety. IEEE_SUPPORT_ALL() Description. Inquire if processor supports all the IEEE facilities defined in this standard for all kinds of reals. Class. Inquiry function. Arguments. None. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the results of all the functions IEEE_DATATYPE, IEEE_SUPPORT_DENORMAL, IEEE_SUPPORT_FLAG, IEEE_SUPPORT_HALTING, IEEE_SUPPORT_INF, IEEE_SUPPORT_NAN, IEEE_SUPPORT_ROUNDING, and IEEE_SUPPORT_SQRT are always true; otherwise, it has the value false. Example. IEEE_SUPPORT_ALL() has the value false if the processor supports both IEEE and non-IEEE kinds of reals. IEEE_SUPPORT_DENORMAL(X) Description. Inquire if the processor supports IEEE gradual underflow for reals of the same kind type parameter as the argument. Class. Inquiry function. Argument. X shall of type real and such that IEEE_DATATYPE(X) has the value true. It may be scalar or array valued. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports IEEE gradual underflow for real variables of the same kind type parameter as X; otherwise, it has the value false. Example. IEEE_SUPPORT_DENORMAL(X) has the value true if the processor supports gradual underflow for X. IEEE_SUPPORT_FLAG(FLAG) Description. Inquire if the processor supports an exception or set of exceptions. Class. Inquiry function. Argument. FLAG shall be scalar and of type TYPE(IEEE_FLAGS). It shall have one of the values: INEXACT, UNDERFLOW, OVERFLOW, INVALID, DIVIDE_BY_ZERO, USUAL and ALL. USUAL specifies the exceptions INEXACT, OVERFLOW, and DIVIDE_BY_ZERO. ALL specifies the exceptions INEXACT, UNDERFLOW, OVERFLOW, INVALID, and DIVIDE_BY_ZERO. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports detection of the specified exception or exceptions; otherwise, it has the value false. Example. IEEE_SUPPORT_FLAG(ALL) has the value true if the processor supports all the exceptions. IEEE_SUPPORT_HALTING() Description. Inquire if the processor supports the ability to control during program execution whether to abort or continue execution after an exception. Class. Inquiry function. Arguments. None Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports the ability to control during program execution whether to abort or continue execution an exception; otherwise, it has the value false. Example. IEEE_SUPPORT_HALTING() has the value true if the processor supports halting after an exception. IEEE_SUPPORT_INF(X) Description. Inquire if processor supports the IEEE infinity facility for reals of the same kind type parameter as the argument. Class. Inquiry function. Argument. X shall of type real and such that IEEE_DATATYPE(X) has the value true. It may be scalar or array valued. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports IEEE infinities (positive and negative) for real variables of the same kind type parameter as X; otherwise, it has the value false. Example. IEEE_SUPPORT_INF(X) has the value true if the processor supports IEEE infinities for X. IEEE_SUPPORT_NAN(X) Description. Inquire if processor supports the IEEE Not-A-Number facility for reals of the same kind type parameter as the argument. Class. Inquiry function. Argument. X shall of type real and such that IEEE_DATATYPE(X) has the value true. It may be scalar or array valued. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports IEEE NaNs for real variables of the same kind type parameter as X; otherwise, it has the value false. Example. IEEE_SUPPORT_NAN(X) has the value true if the processor supports IEEE NaNs for X. IEEE_SUPPORT_ROUNDING(RND_VALUE) Description. Inquire if processor supports changing to a particular rounding mode for IEEE kinds of reals. Class. Inquiry function. Argument. RND_VALUE shall be of type TYPE(IEEE_RND_VALUE) and value one of NEAREST, TO_ZERO, UP, and DOWN, which are named constants of the module. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports the IEEE rounding mode defined by RND_VALUE for IEEE kinds of reals; otherwise, it has the value false. Example. IEEE_SUPPORT_ROUNDING(TO_ZERO) has the value true if the processor supports rounding to zero for IEEE kinds of reals. IEEE_SUPPORT_SQRT(X) Description. Inquire if the processor supports IEEE square root for reals of the same kind type parameter as the argument. Class. Inquiry function. Argument. X shall of type real and such that IEEE_DATATYPE(X) has the value true. It may be scalar or array valued. Result Characteristics. Default logical scalar. Result Value. The result has the value true if the processor supports IEEE square root for real variables of the same kind type parameter as X; otherwise, it has the value false. Example. IEEE_SUPPORT_SQRT(X) has the value true if the processor supports IEEE square root for X. IEEE_UNORDERED(X,Y) Description. IEEE unordered function. True if X or Y is a NaN. and false otherwise. Class. Elemental function. Arguments. The arguments shall be of type real and such that IEEE_DATATYPE(X) and IEEE_DATATYPE(Y) have the value true. Result Characteristics. Same as X. Result Value. The result has the value true if X or Y is a NaN or both are NaNs; otherwise, it has the value false. Example. IEEE_UNORDERED(0.,IEEE_SQRT(-1.0)) has the value true. <<15.10 Examples>> Example 1: MODULE DOT ! Module for dot product of two real arrays of rank 1. USE IEEE_ARITHMETIC LOGICAL MATRIX_ERROR INTERFACE OPERATOR(.dot.) MODULE PROCEDURE MULT END INTERFACE CONTAINS REAL FUNCTION MULT(A,B) REAL, INTENT(IN) :: A(:),B(:) INTEGER I LOGICAL OLD_OVERFLOW IF (SIZE(A)/=SIZE(B)) THEN MATRIX_ERROR = .TRUE. RETURN END IF OLD_OVERFLOW = IEEE_GET_FLAG(OVERFLOW) CALL IEEE_CLEAR_FLAGS(OVERFLOW) MULT = 0. DO I = 1, SIZE(A) MULT = MULT + A(I)*B(I) END DO IF (IEEE_GET_FLAG(OVERFLOW) ) THEN MATRIX_ERROR = .TRUE. ELSE IF(OLD_OVERFLOW) CALL IEEE_SET_FLAGS(OVERFLOW) END IF END FUNCTION MULT END MODULE DOT This module provides the dot product of two real arrays of rank 1. If the sizes of the arrays are different, an immediate return occurs with MATRIX_ERROR true. If OVERFLOW occurs during the actual calculation, the overflow flag will be set and MATRIX_ERROR will be true. Example 2: USE IEEE_ARITHMETIC TYPE(IEEE_FSR_VALUE) FSR_VALUE FSR_VALUE = IEEE_GET_FSR() CALL IEEE_CLEAR_FLAGS (OVERFLOW, INVALID, DIVIDE_BY_ZERO) ! First try the "fast" algorithm for inverting a matrix: MATRIX1 = FAST_INV (MATRIX) ! This must not alter MATRIX. IF (IEEE_GET_FLAG(OVERFLOW) .OR. IEEE_GET_FLAG(INVALID) .OR. & IEEE_GET_FLAG(DIVIDE_BY_ZERO) ) THEN ! "Fast" algorithm failed; try "slow" one: CALL IEEE_CLEAR_FLAGS (OVERFLOW, INVALID, DIVIDE_BY_ZERO) MATRIX1 = SLOW_INV (MATRIX) IF (IEEE_GET_FLAG(OVERFLOW) .OR. IEEE_GET_FLAG(INVALID) .OR. & IEEE_GET_FLAG(DIVIDE_BY_ZERO) ) THEN WRITE (*, *) 'Cannot invert matrix' STOP END IF END IF CALL IEEE_SET_FSR(FSR_VALUE) In this example, the function FAST_INV may cause a condition to signal. If it does, another try is made with SLOW_INV. If this still fails, a message is printed and the program stops. Note, also, that it is important to clear the flags before the second try. Example 3: USE IEEE_ARITHMETIC : ! First try a fast algorithm for inverting a matrix. CALL IEEE_CLEAR_FLAGS (OVERFLOW) DO K = 1, N : IF (IEEE_GET_FLAG(OVERFLOW)) EXIT END DO IF (IEEE_GET_FLAG(OVERFLOW)) THEN ! Alternative code which knows that K-1 steps have executed normally. : END IF Here the code for matrix inversion is in line and the transfer is made more precise by adding extra tests of the flag. Example 4: REAL FUNCTION HYPOT(X, Y) ! In rare circumstances this may lead to the setting of the OVERFLOW flag USE IEEE_ARITHMETIC REAL X, Y REAL SCALED_X, SCALED_Y, SCALED_RESULT LOGICAL OLD_OVERFLOW, OLD_UNDERFLOW INTRINSIC SQRT, ABS, EXPONENT, MAX, DIGITS, SCALE ! Store the old flags and clear them OLD_OVERFLOW = IEEE_GET_FLAG(OVERFLOW) OLD_UNDERFLOW = IEEE_GET_FLAG(UNDERFLOW) CALL IEEE_CLEAR_FLAGS(UNDERFLOW, OVERFLOW) ! Try a fast algorithm first HYPOT = SQRT( X**2 + Y**2 ) IF (IEEE_GET_FLAG(UNDERFLOW) .OR. IEEE_GET_FLAG(OVERFLOW) ) THEN CALL IEEE_CLEAR_FLAGS(UNDERFLOW, OVERFLOW) IF ( X==0.0 .OR. Y==0.0 ) THEN HYPOT = ABS(X) + ABS(Y) ELSE IF ( 2*ABS(EXPONENT(X)-EXPONENT(Y)) > DIGITS(X)+1 ) THEN HYPOT = MAX( ABS(X), ABS(Y) )! one of X and Y can be ignored ELSE ! scale so that ABS(X) is near 1 SCALED_X = SCALE( X, -EXPONENT(X) ) SCALED_Y = SCALE( Y, -EXPONENT(X) ) SCALED_RESULT = SQRT( SCALED_X**2 + SCALED_Y**2 ) HYPOT = SCALE( SCALED_RESULT, EXPONENT(X) ) ! may cause overflow END IF END IF IF(OLD_OVERFLOW) CALL IEEE_SET_FLAGS(OVERFLOW) IF(OLD_UNDERFLOW) CALL IEEE_SET_FLAGS(UNDERFLOW) END FUNCTION HYPOT An attempt is made to evaluate this function directly in the fastest possible way. (Note that with hardware support, exception checking is very efficient; without language facilities, reliable code would require programming checks that slow the computation significantly.) The fast algorithm will work almost every time, but if an exception occurs during this fast computation, a safe but slower way evaluates the function. This slower evaluation may involve scaling and unscaling, and in (very rare) extreme cases this unscaling can cause overflow (after all, the true result might overflow if X and Y are both near the overflow limit). If the overflow or underflow flag is set on entry, it is reset on return, so that earlier exceptions are not lost.