Hi,
there is several ways to declare a double precision real in Fortran:
real*8 :: w
real(8) :: x
real(kind(0d0)) :: y
double precision :: z
My habits are to use double precision. Is it a good practice ?
The following is more or less, my go-to module that I use it everywhere:
module mod_types
use iso_fortran_env, only: int8, int16, int32, int64, real32, real64, &
input_unit, output_unit, error_unit
implicit none
integer, parameter :: sp = real32
integer, parameter :: dp = real64
integer, parameter :: stdin = input_unit
integer, parameter :: stdout = output_unit
integer, parameter :: stderr = error_unit
end module
An example use is:
module foo
use mod_types, only: wp=>dp
implicit none
real(wp)::bar
end module
wp stands for working precision and if for some reason I want to change precision, I just change this wp=>sp so I do not need to change manually all real variables. I believe something in that line is considered the “standard way” now.
7 Likes
I also use the constants defined in iso_fortran_env, though I have recently learned (from this discussion in stdlib) that there are still concerns over the portability of these constants since they only guarantee a byte width but nothing about precision.
In his post “It takes all kinds”, Steve Lionel recommends use of the intrinsic function selected_real_kind above all other methods.
It’s also worth mentioning the constants provided by iso_c_binding which are required for c interoperability:
use iso_c_binding, only: sp=>c_float, dp=>c_double
8 Likes
Note that the real*8 syntax is an extension and not standard. Also with real(8) you are assuming that “real kind 8” is an 8-byte real (which conceivably isn’t even double precision), but there’s nothing in the standard that says kind values correspond to the byte size. In fact with the NAG compiler this isn’t the case. They number their kinds consecutively from 1 (though there is a compiler flag you can use to change that).
1 Like
If one truly knows what decimal precision and exponent range they need for their problem, they can use selected_real_kind and let the system decide the kind (and hope it has one). But if you’re like me, all I really know is that the hardware supports 32-bit and 64-bit floating point and that I want to use 64-bit. So in the past I’d figure out a precision and range that would ensure that selected_real_kind gave me the 64-bit kind. Thankfully that silliness isn’t needed anymore with the constants from iso_fortran_env. I hate to disagree with Steve, but if you know you want a specific fp size, just use the constant and don’t play games with selected_real_kind.
4 Likes
You can disagree with me, that’s fine. All I wanted to say in that post was that use of real32, real64 etc. was little better than real(4) and real(8). If you’re aware that all you’re doing is asking for a physical storage size (and no other properties), and document that, then OK. But I would be happier to see a module that defines constants using selected_real_kind() - you need do this only once.
Of course, if you’re interfacing with C, then use the kinds from ISO_C_BINDING.
P.S. Thank you for reading my blog posts!
10 Likes
Thank you all for your valuable answers (and “It Takes All KINDs” was so interesting and pleasant to read). This is a very interesting discussion, I am learning or re-learning many things.
Perhaps using iso_fortran_env or selected_real_kind depends on the kind of computer scientist you are:
- If you work only on a PC, like myself,
iso_fortran_envcan probably do the job. That kind of machine is relatively stable in time (concerning the reals in the CPU). Although, if GPU acceleration is democratized things can become different. And in 20 years, I don’t know what kind of CPU we will use. - If you compute on supercomputers, I guess it is more complicated. Your code can run on very different architectures and
selected_real_kindcan surely be very useful. Moreover, you probably work in a bigger team and your code will perhaps still used by other people in 20 years, so keeping the memory of the needed precision will be helpful.
Now, what I am interested in is getting the best compromise between precision and speed. On my PC, I can use 32 bits real, 64 bits, 80 bits, 128 bits:
- but 32 bits is not interesting because I have half the precision of 64 bits with exactly the same computing time. Just a relic from the good old days.
- REAL(10) in gfortran is 80 bits, the computing time is 50% higher but I gain around 5 digits. It can be interesting if I really need a greater precision.
- 128 bits reals doubles the number of digits, but it is not “hardware” (hard wired) reals. It is based on libquadmath and the CPU time is 10 times longer ! Not an option for me…
My question would now be: is it possible with Fortran to know the biggest “hardware” reals offered by the CPU ?
Another question: I thought 80 bits were hard wired, no ? But in that case why the computing time is 50 % higher than with 64 bits reals ?
2 Likes
80-bit floating is only possible with the old x87 instructions. You can’t take advantage of any of the vector SSE or AVX instructions that are generally much faster.
3 Likes
Since you mention speed, it should be pointed out that there are performance advantages to using 32bit, even though the fundamental operations (+ / *) are computed in the same time as 64bit:
-
Functions implemented with polynomial series will be quicker since fewer terms are needed for the reduced precision;
-
If your program can be vectorized, then throughput will be higher with 32bit than 64bit since more values can be placed in the vector registers;
-
Many intensive programs are memory-bound these days, so if your program operates on a large-enough dataset (i.e. one that doesn’t fit in the cache), and it involves a lot of memory accesses per flop, then you will benefit from ~2x speedup when using 32bit simply by being able to load more values into the processor.
Obviously, these aren’t options if you need double precision (though there are so-called mixed-precision solvers that attempt to exploit this).
6 Likes
@sblionel
Thank you for you explanation. I suspected it was that, but I was surprised that those 80 bits reals works also with gfortran -mfpmath=sse. I guess that even with this flag the compiler choose to use x87 when needed.
@lkedward
Thank you very much for you answer, you are right, I was convinced there was no difference, based on arithmetical measurements.
When computing 2e9 cos(), the 32 bits version is 2.66x times faster than the 64 bits.
Finally, it’s another good reason to adopt a practice like writing real(wp)::bar. By just changing the wp value I could study the difference in computing time.
2 Likes
Hi everyone,
About the nag compiler, you’re right, kind=8 is not working with the default option. However, when the “-kind=byte” option is used, old codes (and new ones) can be compile with real*8 or real (kind=8).
Why does this have to be so complicated? (it’s a rhetorical question…I know how we got here).
But the fact is, this is all very confusing for beginners and people coming from other languages. This is why every day on stack overflow we see questions about why isn’t this working right:
real*8 pi = 3.14159265358979323
Making people have to use iso_fortran_env or selected_real_kind just to declare a real number is just user-hostile. Nobody wants to type all that (hence real*8).
Is there nothing that can be done about it?
In most other languages that people might come to Fortran from, double precision is the default. Can we just make double precision the default in Fortran? I feel like that would fix 80% of the problems people have.
It was said many times. If one hates iso-s and kinds, one can still write
DOUBLE PRECISION x=3.14159265358979323d0
On comp.lang.fortran I asked if we can Standardize real*8 as real(kind=real64)? There was not much support for doing so. I don’t think making double precision the default is feasible, because it would break the codes calling procedures with single precision dummy arguments where the actual argument is a literal such as 3.14159265358979323.
Note, also, that the C world is not that-better. Indeed, the literal FP constants in C default to double but it is only because the pre-ANSI C did all FP arithmetic in double precision. But then if you want to define a long double value, you have to add L suffix, otherwise it will get truncated to double precision before being assigned to long double variable.
PS. for some reason, this is not the case with integer C literal constants, for which the compiler assigns the first (“shortest”, starting from int) type which contains the value.
I am pretty sure that because there are multiple ways of defining types in Fortran and people are unsure as to which on is “the correct one”, my research group eventually gave up and are using the way more satanic flags -fdefault-real-8 and -fdefault-double-8
1 Like
That is great comment 
Intel Fortran can add flags to define the size of default integer, real, double precision,
gfortran can easily do such things too, as @gnikit pointed out.
I like this style.
real(8) :: pi = 3.14159265358979323d0
I do the following,
integer, parameter :: r8=selected_real_kind(15,9)
real(kind=r8), parameter :: pi=4.0_r8*atan(1.0_r8)
1 Like