Logo Questions Linux Laravel Mysql Ubuntu Git Menu
 

Numerical Precision in Fortran 95:

I have the following Fortran code:

Program Strange
   Real(Kind=8)::Pi1=3.1415926535897932384626433832795028841971693993751058209;
   Real(Kind=8)::Pi2=3.1415926535897932384626433832795028841971693993751058209_8;

   Print*, "Pi1=", Pi1;
   Print*, "Pi2=", Pi2;

End Program Strange

I compile with gfortran, and the output is:

 Pi1=   3.1415927410125732     
 Pi2=   3.1415926535897931

Of course the second is correct, but should this be the case? It seems like Pi1 is being input to memory as a single precision number, and then put into a double precision memory slot. But this seems like an error to me. Am I correct?

like image 895
user14717 Avatar asked May 04 '13 03:05

user14717


People also ask

How do you define precision in Fortran?

The decimal precision is the number of significant digits, and the decimal exponent range specifies the smallest and largest representable number. The range is thus from 10-r to 10+r.

How do you write double precision in Fortran 90?

If you do not specify a format, GNU Fortran will print real numbers using about 9 digits, even if you do calculations in double precision. If you want to print in double precision you must use write and format statements. When double precision is used the maximum number of digits possible is 17.

How many digits is double precision in Fortran?

Double precision provides greater range (approximately 10**(-308) to 10**308) and precision (about 15 decimal digits) than single precision (approximate range 10**(-38) to 10**38, with about 7 decimal digits of precision).

What is Signal comment in Fortran 95?

A comment indicator or a statement number may precede the tab. If a tab is the first nonblank character, then: If the character after the tab is anything other than a nonzero digit, then the text following the tab is an initial line. If there is a nonzero digit after the first tab, the line is a continuation line.


2 Answers

I do know a bit of Fortran ! @Dougal's answer is correct though the snippet he quotes from is not, embedding the letter d into a real literal constant is not required (since Fortran 90), indeed many Fortran programmers now regard that approach as archaic. The snippet is also misleading in advising the use of 3.1415926535d+0 to initialise a 64-bit floating-point value for pi, it doesn't set enough of the digits to their correct values.

The statement:

Real(Kind=8)::Pi1=3.1415926535897932384626433832795028841971693993751058209

defines Pi1 to be a real variable of kind 8. The literal real value 3.1415926535897932384626433832795028841971693993751058209 is, however, a real value of default kind, most likely to be a 4-byte real on most current compilers. That seems to explain your output but do check your documentation.

On the other hand, the literal real value Pi2=3.1415926535897932384626433832795028841971693993751058209_8 is, by the suffixing of the kind specification, declared to be of kind=8 which is the same as the kind of the variable it is assigned to.

Three more points:

1) Don't fall into the trap of thinking that kind=8 means the same thing as 64-bit floating-point number or double. For many compilers it does, for some it doesn't. Kind numbers are not portable between Fortran implementations. They are, according to the standard, arbitrary positive integers. Better, with a modern compiler, would be to use the predefined constants from the intrinsic module iso_fortran_env, e.g.

use, intrinsic :: iso_fortran_env
...
real(real64) :: pi = 3.14159265358979323846264338_real64

There are other portable approaches to setting variable kinds using functions such as selected_real_kind.

2) Since the value of pi is unlikely to change during the execution of your program you might care to make it a parameter thus:

real(real64), parameter :: pi = 3.14159265358979323846264338_real64

3) It isn't necessary (or usual) to end Fortran statements with a ';' unless you want to have more than one statement on the same line in the source file.

like image 194
High Performance Mark Avatar answered Sep 28 '22 08:09

High Performance Mark


I don't really know Fortran, but this page says:

The letter "d" must be embedded in the literal, otherwise, the compiler's pre-processor would round it off to be a Single Precision literal. For example, 3.1415926535 would be read as 3.141593 while 3.1415926535d+0 would be stored with all the digits intact. The letter "d" for double precision numbers has the same meaning as "e" for single precision numbers.

So it seems like your guess is correct.

like image 25
Danica Avatar answered Sep 28 '22 08:09

Danica