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:
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.