Getting REAL with Fortran¶
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:
DOUBLE PRECISION Y
What we knew most generally when comparing these two types was that
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
would occupy 32 bits and be IEEE single precision.
Literal constants exist in corresponding ways:
As extensions beyond FORTRAN 77, some compilers would accept declaration styles with byte counts:
C A variable taking 4 bytes
C A variable taking 8 bytes
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
We still have at least two approximation methods for real numbers and these match those of FORTRAN 77:
double precision y1
x2 are both of the same ("default") kind and
are of the same kind.
Literal constants have corresponding designations of kind in addition to
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 = ...
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)
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.
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:
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¶
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
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
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)
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
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
RANGE tell us about the
particular representation for a given entity.
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
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
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
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
compiler flags such as
-freal-8-real-16 may cause
to correspond to the value
REAL128 rather than
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.
This post was inspired by discussion on the UK RSE Slack workspace.