rse,

Getting REAL with Fortran

Iain Barrass Iain Barrass Follow Jun 12, 2019 · 11 mins read
Getting REAL with Fortran
Share this

Fortran provides a variety of intrinsic representations of real numbers. In this post we look at what these representations are and how we choose a particular representation for our work.

When writing a Fortran program to work with floating point numbers we have a choice between a number of representations of the real type. A compiler will offer at least two such kinds. Here we’ll look at how one chooses a kind parameter for a variable.

Representing real numbers in FORTRAN 77

Back in the days of FORTRAN 77 we had exactly two types representing real numbers:

      REAL X
      DOUBLE PRECISION Y

What we knew most generally when comparing these two types was that y would take up twice the storage space as x and would have a greater decimal precision in the approximation of the real number.

If we had a specific compiler and computer system in mind we would know additional information about these variables (which may be modified by flags or options when compiling). For example, we may know that x would occupy 32 bits and be IEEE single precision.

Literal constants exist in corresponding ways:

      X=0.
      X=0E0
      Y=0D0

As extensions beyond FORTRAN 77, some compilers would accept declaration styles with byte counts:

C A variable taking 4 bytes
      REAL*4 X
C A variable taking 8 bytes
      REAL*8 Y

Representing real numbers in modern Fortran

Fortran 90 introduced the concept of kind type parameters which continues in to modern Fortran. In a type declaration for a real entity we can specify the value of the parameter:

  integer, parameter :: rk=1
  real(kind=27) x
  real(rk) y

We still have at least two approximation methods for real numbers and these match those of FORTRAN 77:

  real x1
  real(KIND(0e0)) x2
  double precision y1
  real(KIND(0d0)) y2

Here, x1 and x2 are both of the same (“default”) kind and y1 and y2 are of the same kind.

Literal constants have corresponding designations of kind in addition to the forms 0., 0e0 and 0d0 as above:

  x = 0._27
  y = 1.35e7_rk

Selecting a real’s kind parameter

We ideally use a named constant to specify the kind parameter when declaring a real variable:

  integer, parameter :: rk = ...
  real(kind=rk) x

How do we determine the value of the named constant rk we want to use?

We could say rk=8, because we know that our compiler uses that value to denote the one we want. We shouldn’t do that, though: not only are we limiting the portability of the code to other compilers or systems, we don’t even know that our compiler will always use this value for the real we want.

Instead, we can consider:

  • the value corresponding to default real or double precision;
  • interoperability with the companion processor;
  • the characteristics of the approximation method;
  • the storage size of the real.

Kinds of other entities

As seen above, our named constant could be defined by asking about the kind parameter (using the intrinsic function KIND) of another entity, such as a default real or double precision literal constant:

integer, parameter :: rk_default = KIND(0.)
integer, parameter :: rk_double = KIND(0d0)

Interoperable kinds

The intrinsic module iso_c_binding has named constants for specifying kinds for interoperability with corresponding C types:

use, intrinsic :: iso_c_binding
integer, parameter :: rk_cfloat = C_FLOAT
integer, parameter :: rk_cdouble = C_DOUBLE
integer, parameter :: rk_clongdouble = C_LONG_DOUBLE

Any of these named constants may have negative value, signifying that there is no matching Fortran real for the corresponding C type.

Characteristics of the approximation

The intrinsic function SELECTED_REAL_KIND returns, if one is available, the kind of a real having certain properties. The function takes up to three arguments, corresponding to: decimal precision, decimal exponent range; radix. The function IEEE_SELECTED_REAL_KIND of the intrinsic module ieee_arithmetic is similar, but the resulting kind is that of a real satisfying the requirements of “IEEE arithmetic”.

use, intrinsic :: ieee_arithmetic, only : IEEE_SELECTED_REAL_KIND
integer, parameter :: rk_precise = SELECTED_REAL_KIND(20,300)
integer, parameter :: rk_ieee = IEEE_SELECTED_REAL_KIND(15,307,2)

Each of these functions may return, instead of a valid kind number, a negative value which cannot be used as a kind type parameter. These negative results indicate that there exists no approximation method with the desired specification.

Storage size

The intrinsic module iso_fortran_env provides named constants to select a kind number based on the storage size of the corresponding real: defined in Fortran 2008 and 2018 are: real32, real64 and real128. If such a constant has a positive value then that value is a kind number for a real with corresponding storage size. real(real32) x, for example, declares a real variable occupying 32 bits.

use, intrinsic :: iso_fortran_env, only : rk_real64 => REAL64

If there is no real with a given storage size then the constants will have negative value (and cannot be used in declaring a real entity).

Uniqueness of kind numbers from function or named constant

The functions SELECTED_REAL_KIND and IEEE_SELECTED_REAL_KIND have a single result. Several approximation methods may match the request and in such cases the approximation with the smallest matching decimal precision is denoted.

Equally, the compiler may have more than one real kind with a given storage size. Here one of those kind values is given by the corresponding constant but other than the storage size there is no requirement on which is the matching approximation method.

Supported kind numbers

Fortran 2008 specified the named constant REAL_KINDS in the intrinsic module iso_fortran_env. This constant is an array containing the values of the kind numbers of the supported approximation methods:

  use, intrinsic :: iso_fortran_env, only : REAL_KINDS
  integer, parameter :: rk_guess = MAXVAL(REAL_KINDS)
  real(rk_guess) x
end

Properties of the approximation method

With our Fortran programs we often care about the properties of a real variable which relate to numerical aspects such as precision and exponent range. The function SELECTED_REAL_KIND guarantees that the chosen kind meets these minimum criteria (or even has a specific radix). The function IEEE_SELECTED_REAL kind behaves similarly but goes further in giving us an IEEE-compliant kind. However, the resulting kinds may far exceed our requested minimum support which may cause problems in other ways.

For reals selected by the storage size constants, or the non-standard forms like real*8, we know the storage size but that tells us little about the numeric properties. In particular, a real with storage size 128-bits may not be using all of those bits for the numeric model, and if it is it may be using either IEEE quadruple precision or IEEE double-double precision.

The intrinsic functions PRECISION, RADIX, RANGE tell us about the particular representation for a given entity.

Portability

There are two aspects of portability to consider when choosing a real’s kind:

  • compiling and working in a variety of places;
  • getting the “same or equivalent” results in a variety of places.

Portability between compilers/systems

For portability between compilers it is best to avoid using literal constants as kind type parameters. Although values such as 8 are commonly used by compilers to denote an 8-byte real, this is not guaranteed. The NAG Fortran compiler, for example, does not use this numbering scheme by default:

Error: KIND value (8) does not specify a valid representation method

The non-standard style real*4 should also be avoided.

The named constants REAL32, REAL64 and REAL128 were added in the Fortran 2008 revision, so care should be taken when considering portability with “old” compilers. Requesting specific storage sizes may limit portability, in particular with REAL128.

The function IEEE_SELECTED_REAL_KIND is part of Fortran 2003’s IEEE arithmetic facility and naturally requires support for IEEE types. The radix argument to IEEE_SELECTED_REAL_KIND was added in Fortran 2008 and so use may be limiting.

Reproducibility of results

The Fortran language standard is written in such a way as not to impose a precise numerical model on the real intrinsic type. This means that it is not (currently) possibly to write a Fortran program which would have a real entity having exactly the same representation when interpreted by any Fortran processor.

In practice, however, most compilers on commonly available hardware will offer binary IEEE floating point representations binary32 and binary64 and these would be the default real and double precision kinds. SELECTED_REAL_KIND(15,307) is a reasonably portable way of selecting a binary64 representation and this may be preferred over simply KIND(0d0): compiler flags such as gfortran’s -freal-8-real-16 may cause KIND(0d0) to correspond to the value REAL128 rather than REAL64.

Conclusion

In this post we’ve looked at various ways of finding a valid and desired kind parameter for a real entity. We’ve seen that we have a lot of choice and that using a literal constant such as in real(8) x is not necessary (and is harmful).

How we choose a kind type parameter is partly determined by how we intend to use the entity: it may be for particular numerical reasons, for interoperability with a C library, to be a particular size, or just to be most portable between compilers with minimal changes.

When we have selected a kind type parameter, we’ve seen that we can query the numerical model to find its characteristics, such as precision, range and radix. Although we can’t guarantee exact matches between all systems these inquiry functions will at least tell us how the entity behaves. We may even use our build system to choose a particular type for the closest match on a new build target.

In another post we’ll look at how, even with exactly matching representations for real approximations, we will need to be aware of differences in numerical results.

Acknowledgements

This post was inspired by discussion on the UK RSE Slack workspace.

Iain Barrass
Written by Iain Barrass
RSE Team Leader. He prefers Fortran, with a background in infectious diseases